This notebook works through the logistic regression lecture in R, in the same order as the slides. Every model, test and number shown in the lecture is computed here from the data, so you can run each step yourself and see where the numbers come from.

We use the same data as in Topic 1: clinical_sample.csv, with 280 participants measured on self-esteem, negative affect, depression, gender (0 = man, 1 = woman) and treatment. The data file should be in the same folder as this notebook.

1 The data and the binary outcome

In Topic 1, depression was a continuous outcome. Here we turn it into a yes/no outcome: high_depression = 1 if a person’s depression score is above the sample mean, and 0 otherwise.

d <- read.csv("clinical_sample.csv")
str(d)
## 'data.frame':    280 obs. of  6 variables:
##  $ id         : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ self_esteem: num  19 28.2 37.6 34.3 28.7 23 30.5 21.9 23.2 29.4 ...
##  $ neg_affect : num  33.1 19.5 24.8 20.1 30.8 23.2 19.7 24.1 18.9 21.3 ...
##  $ treatment  : int  0 0 1 1 0 0 1 0 1 1 ...
##  $ gender     : int  1 1 0 0 0 0 0 1 1 0 ...
##  $ depression : num  16.8 14.2 0.2 6.1 7.1 13 10.3 9.8 6.4 1.9 ...
d <- d %>% mutate(high_depression = as.integer(depression > mean(depression)))

mean(d$depression)
## [1] 9.110714
table(d$high_depression)
## 
##   0   1 
## 145 135

The mean depression score is 9.11, and 135 of the 280 participants (48.2%) score above it.

2 Why not just fit a straight line?

We first try ordinary regression (OLS) on the 0/1 outcome, exactly as in Topic 1.

ols <- lm(high_depression ~ self_esteem, data = d)
coef(ols)
## (Intercept) self_esteem 
##  1.16278993 -0.02373394
# predictions inside the observed range of self-esteem (10 to 40) ...
predict(ols, newdata = data.frame(self_esteem = c(10, 40)))
##         1         2 
## 0.9254505 0.2134323
# ... and a little outside it
predict(ols, newdata = data.frame(self_esteem = c(5, 55)))
##          1          2 
##  1.0441202 -0.1425769
# where the line crosses 1 and 0
(1 - coef(ols)[1]) / coef(ols)[2]
## (Intercept) 
##     6.85895
(0 - coef(ols)[1]) / coef(ols)[2]
## (Intercept) 
##     48.9927

Inside the observed range the line stays between 0 and 1, but only just (from 0.93 down to 0.21). A little outside it, the line predicts “probabilities” above 1 and below 0: it crosses 1 at a self-esteem of about 6.9 and 0 at about 49. A straight line has no floor and no ceiling, which is why we need a different model.

3 From probability to log-odds

3.1 Odds and log-odds

The odds are \(p/(1-p)\), and the log-odds (the logit) are \(\log\big(p/(1-p)\big)\). In R, log() is the natural logarithm, the same \(\log\) used throughout the lecture. R also has ready-made functions: qlogis() turns a probability into log-odds, and plogis() turns log-odds back into a probability (the inverse logit).

p <- 0.8
p / (1 - p)          # odds
## [1] 4
log(p / (1 - p))     # log-odds, by hand
## [1] 1.386294
qlogis(p)            # log-odds, with R's logit function
## [1] 1.386294
plogis(1.386294)     # back from log-odds to probability
## [1] 0.7999999

3.2 The stretched y-axis, tick by tick

These are the seven probabilities from the step-by-step slides. Equal steps on the log-odds ruler (the last column) are not equal steps on the probability ruler (the first column): the gaps shrink towards 0 and 1.

ticks <- tibble(p = plogis(-3:3)) %>%
  mutate(odds = p / (1 - p),
         log_odds = log(odds))
ticks %>% kable(digits = 3)
p odds log_odds
0.047 0.050 -3
0.119 0.135 -2
0.269 0.368 -1
0.500 1.000 0
0.731 2.718 1
0.881 7.389 2
0.953 20.086 3

A probability of 0.5 is odds of 1 and log-odds of 0. At the ends, the log-odds of a probability of exactly 1 or 0 are \(+\infty\) and \(-\infty\):

qlogis(c(0, 1))
## [1] -Inf  Inf

This is why every observed person (who is exactly 0 or 1) sits at \(\pm\infty\) on the log-odds scale, and why the straight line cannot be fitted by least squares on that scale.

4 Maximum likelihood: how the straight line is fitted

4.1 The likelihood, person by person

A candidate straight line gives each person a predicted probability \(p_i\). A person with high depression (\(Y=1\)) contributes \(p_i\); a person without it (\(Y=0\)) contributes \(1-p_i\). The log-likelihood adds up the logs of all 280 contributions.

loglik <- function(b0, b1) {
  p_i <- plogis(b0 + b1 * d$self_esteem)
  contribution <- ifelse(d$high_depression == 1, p_i, 1 - p_i)
  sum(log(contribution))
}

For example, under the line \(\text{log-odds} = 2.818 - 0.101 \times \text{self-esteem}\), participants 1 and 3 contribute:

d[c(1, 3), c("id", "self_esteem", "high_depression")] %>%
  mutate(predicted_p = plogis(2.818324 - 0.100899 * self_esteem),
         contribution = ifelse(high_depression == 1, predicted_p, 1 - predicted_p)) %>%
  kable(digits = 3)
id self_esteem high_depression predicted_p contribution
1 1 19.0 1 0.711 0.711
3 3 37.6 0 0.274 0.726

4.2 Trying candidate lines

Three candidate lines, each scored by its log-likelihood (closer to 0 is better). The flat line ignores self-esteem: its intercept is simply the log-odds of the overall proportion with high depression.

flat_b0 <- qlogis(mean(d$high_depression))

candidates <- tibble(
  line = c("flat line (ignores self-esteem)", "too steep", "maximum likelihood"),
  B0   = c(flat_b0, 0.300 * 27.932, 2.818324),
  B1   = c(0, -0.300, -0.100899)
) %>%
  mutate(LL = map2_dbl(B0, B1, loglik))

candidates %>% kable(digits = c(0, 3, 3, 2))
line B0 B1 LL
flat line (ignores self-esteem) -0.071 0.000 -193.90
too steep 8.380 -0.300 -210.10
maximum likelihood 2.818 -0.101 -185.21

The maximum-likelihood line has the highest log-likelihood. glm() does this search for us and finds exactly this line:

mod1 <- glm(high_depression ~ self_esteem, data = d, family = binomial(link = "logit"))
coef(mod1)
## (Intercept) self_esteem 
##    2.818325   -0.100899
logLik(mod1)
## 'log Lik.' -185.2054 (df=2)

4.3 From likelihood to deviance

R reports deviance \(= -2 \times\) log-likelihood. The deviance of our model is the “Residual deviance”; the deviance of the flat line (the intercept-only “null” model) is the “Null deviance”.

mod0 <- glm(high_depression ~ 1, data = d, family = binomial(link = "logit"))

-2 * as.numeric(logLik(mod1))   # residual deviance of Model 1
## [1] 370.4109
-2 * as.numeric(logLik(mod0))   # null deviance
## [1] 387.8052
c(deviance(mod1), deviance(mod0))
## [1] 370.4109 387.8052

5 Model 1: does self-esteem predict high depression?

summary(mod1)
## 
## Call:
## glm(formula = high_depression ~ self_esteem, family = binomial(link = "logit"), 
##     data = d)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  2.81832    0.73311   3.844 0.000121 ***
## self_esteem -0.10090    0.02526  -3.994 6.49e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 387.81  on 279  degrees of freedom
## Residual deviance: 370.41  on 278  degrees of freedom
## AIC: 374.41
## 
## Number of Fisher Scoring iterations: 4

The fitted equation is \(\log\big(p/(1-p)\big) = 2.818 -0.101 \times \text{self-esteem}\).

5.1 The intercept

The intercept is the predicted log-odds at self-esteem = 0. To get the odds, exponentiate the log-odds; to get the probability, use odds / (1 + odds). Remember that self-esteem = 0 is outside the observed range (10 to 40), so this is an extrapolation.

b0 <- coef(mod1)[1]
b0                        # log-odds
## (Intercept) 
##    2.818325
exp(b0)                   # odds
## (Intercept) 
##    16.74877
exp(b0) / (1 + exp(b0))   # probability
## (Intercept) 
##   0.9436581

5.2 The slope, on all three scales

slope_table <- tibble(self_esteem = c(10, 20, 30)) %>%
  mutate(log_odds    = predict(mod1, newdata = ., type = "link"),
         odds        = exp(log_odds),
         probability = odds / (1 + odds))
slope_table %>% kable(digits = 2)
self_esteem log_odds odds probability
10 1.81 6.11 0.86
20 0.80 2.23 0.69
30 -0.21 0.81 0.45
diff(slope_table$log_odds)      # the same change for each 10-point step
##        2        3 
## -1.00899 -1.00899
diff(slope_table$probability)   # a different change for each step
##          2          3 
## -0.1692331 -0.2420210
exp(10 * coef(mod1)[2])         # each 10-point step multiplies the odds by this
## self_esteem 
##   0.3645872

The log-odds drop by the same amount for every 10 points of self-esteem; the probability does not. The odds are multiplied by the same number, 0.365, at every step – the idea behind the odds ratio.

5.3 Where Model 1 reaches each tick of the log-odds ruler

The self-esteem value at which Model 1 predicts a given log-odds value is (log-odds \(- B_0\)) / \(B_1\):

ticks %>%
  mutate(self_esteem = (log_odds - coef(mod1)[1]) / coef(mod1)[2]) %>%
  kable(digits = 3)
p odds log_odds self_esteem
0.047 0.050 -3 57.665
0.119 0.135 -2 47.754
0.269 0.368 -1 37.843
0.500 1.000 0 27.932
0.731 2.718 1 18.021
0.881 7.389 2 8.110
0.953 20.086 3 -1.801

5.4 Using the equation: self-esteem 15 and 35

tibble(self_esteem = c(15, 35)) %>%
  mutate(log_odds    = predict(mod1, newdata = ., type = "link"),
         odds        = exp(log_odds),
         probability = predict(mod1, newdata = ., type = "response")) %>%
  kable(digits = 2)
self_esteem log_odds odds probability
15 1.30 3.69 0.79
35 -0.71 0.49 0.33

6 Model 2: does gender add to self-esteem?

mod2 <- glm(high_depression ~ self_esteem + gender, data = d, family = binomial(link = "logit"))
summary(mod2)
## 
## Call:
## glm(formula = high_depression ~ self_esteem + gender, family = binomial(link = "logit"), 
##     data = d)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  2.62169    0.74943   3.498 0.000468 ***
## self_esteem -0.10397    0.02583  -4.025  5.7e-05 ***
## gender       0.60978    0.25086   2.431 0.015066 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 387.81  on 279  degrees of freedom
## Residual deviance: 364.42  on 277  degrees of freedom
## AIC: 370.42
## 
## Number of Fisher Scoring iterations: 4

The fitted equation is \(\log\big(p/(1-p)\big) = 2.622 -0.104 \times \text{self-esteem} + 0.61 \times \text{gender}\). At any given self-esteem, a woman’s log-odds are 0.61 higher than a man’s.

6.1 Predictions at self-esteem 35

tibble(who = c("woman", "man"), self_esteem = 35, gender = c(1, 0)) %>%
  mutate(log_odds    = predict(mod2, newdata = ., type = "link"),
         odds        = exp(log_odds),
         probability = predict(mod2, newdata = ., type = "response")) %>%
  kable(digits = 2)
who self_esteem gender log_odds odds probability
woman 35 1 -0.41 0.67 0.40
man 35 0 -1.02 0.36 0.27

7 Model 3: does the self-esteem effect differ by gender?

We center self-esteem first, so that 0 means “average self-esteem”:

d <- d %>% mutate(c_self_esteem = self_esteem - mean(self_esteem))
mean(d$self_esteem)
## [1] 28.67821
mod3 <- glm(high_depression ~ c_self_esteem * gender, data = d, family = binomial(link = "logit"))
summary(mod3)
## 
## Call:
## glm(formula = high_depression ~ c_self_esteem * gender, family = binomial(link = "logit"), 
##     data = d)
## 
## Coefficients:
##                      Estimate Std. Error z value Pr(>|z|)   
## (Intercept)          -0.36927    0.17398  -2.123  0.03379 * 
## c_self_esteem        -0.12363    0.03812  -3.243  0.00118 **
## gender                0.61259    0.25091   2.441  0.01463 * 
## c_self_esteem:gender  0.03715    0.05169   0.719  0.47227   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 387.81  on 279  degrees of freedom
## Residual deviance: 363.90  on 276  degrees of freedom
## AIC: 371.9
## 
## Number of Fisher Scoring iterations: 4

The interaction (c_self_esteem:gender) is small and not significant: no evidence that the self-esteem effect differs between men and women. The self-esteem slope for women is the men’s slope plus the interaction:

coef(mod3)["c_self_esteem"] + coef(mod3)["c_self_esteem:gender"]
## c_self_esteem 
##   -0.08647125

8 Odds ratios

An odds ratio is \(\exp(B)\): the number the odds are multiplied by when the predictor increases by one unit.

exp(coef(mod3))
##          (Intercept)        c_self_esteem               gender c_self_esteem:gender 
##            0.6912402            0.8837109            1.8452026            1.0378529
exp(confint(mod3))
##                          2.5 %    97.5 %
## (Intercept)          0.4887982 0.9687986
## c_self_esteem        0.8174421 0.9499398
## gender               1.1317452 3.0308171
## c_self_esteem:gender 0.9377623 1.1493832

A confidence interval for an odds ratio that does not include 1 means a significant effect. Self-esteem and gender are significant; the interaction is not. Note that exp() of the intercept is not an odds ratio: it is the odds for the reference person (a man with average self-esteem).

8.1 Model 3 on the log-odds, odds and probability scales

tibble(who = c("average man", "average woman"),
       c_self_esteem = 0, gender = c(0, 1)) %>%
  mutate(log_odds    = predict(mod3, newdata = ., type = "link"),
         odds        = exp(log_odds),
         probability = odds / (1 + odds)) %>%
  kable(digits = 3)
who c_self_esteem gender log_odds odds probability
average man 0 0 -0.369 0.691 0.409
average woman 0 1 0.243 1.275 0.561
# the woman's odds are the man's odds times the odds ratio for gender
exp(coef(mod3)[1]) * exp(coef(mod3)["gender"])
## (Intercept) 
##    1.275478

8.2 Odds ratios multiply

For Model 1, one point of self-esteem multiplies the odds by the odds ratio; 5 or 10 points multiply them by the odds ratio 5 or 10 times.

OR <- exp(coef(mod1)["self_esteem"])
OR
## self_esteem 
##   0.9040244
OR^5
## self_esteem 
##   0.6038106
OR^10
## self_esteem 
##   0.3645872
(OR - 1) * 100   # percentage change in the odds per point
## self_esteem 
##   -9.597562

9 Model selection

9.1 Likelihood-ratio tests

Each test compares a simpler model with a more complex one. \(\chi^2 = -2(LL_{\text{simpler}} - LL_{\text{complex}})\), with degrees of freedom equal to the number of extra parameters.

lrtest(mod0, mod1)   # does self-esteem improve on the null model?
## Likelihood ratio test
## 
## Model 1: high_depression ~ 1
## Model 2: high_depression ~ self_esteem
##   #Df LogLik Df  Chisq Pr(>Chisq)    
## 1   1 -193.9                         
## 2   2 -185.2  1 17.394  3.037e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lrtest(mod1, mod2)   # does gender improve on Model 1?
## Likelihood ratio test
## 
## Model 1: high_depression ~ self_esteem
## Model 2: high_depression ~ self_esteem + gender
##   #Df  LogLik Df Chisq Pr(>Chisq)  
## 1   2 -185.21                      
## 2   3 -182.21  1 5.991    0.01438 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lrtest(mod2, mod3)   # does the interaction improve on Model 2?
## Likelihood ratio test
## 
## Model 1: high_depression ~ self_esteem + gender
## Model 2: high_depression ~ c_self_esteem * gender
##   #Df  LogLik Df  Chisq Pr(>Chisq)
## 1   3 -182.21                     
## 2   4 -181.95  1 0.5173      0.472

The same test for Model 2 vs. Model 3, computed by hand:

ll2 <- as.numeric(logLik(mod2))
ll3 <- as.numeric(logLik(mod3))
chisq <- -2 * (ll2 - ll3)
chisq
## [1] 0.5172837
pchisq(chisq, df = 1, lower.tail = FALSE)   # p-value
## [1] 0.4720027

Model 1 improves on the null model, and Model 2 improves on Model 1, but Model 3 does not improve on Model 2. We prefer Model 2.

9.2 McFadden’s pseudo-\(R^2\)

ll_null <- as.numeric(logLik(mod0))
tibble(model  = c("Model 1", "Model 2", "Model 3"),
       logLik = c(logLik(mod1), logLik(mod2), logLik(mod3))) %>%
  mutate(McFadden_R2 = 1 - logLik / ll_null) %>%
  kable(digits = 3)
model logLik McFadden_R2
Model 1 -185.205 0.045
Model 2 -182.210 0.060
Model 3 -181.951 0.062

10 Logistic regression as a classifier

For classification we use a model with two continuous predictors, self-esteem and negative affect, and classify a person as “high depression” when the predicted probability is above 0.5.

modC <- glm(high_depression ~ self_esteem + neg_affect, data = d, family = binomial(link = "logit"))
summary(modC)
## 
## Call:
## glm(formula = high_depression ~ self_esteem + neg_affect, family = binomial(link = "logit"), 
##     data = d)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.28153    1.01998  -1.256   0.2090    
## self_esteem -0.05671    0.02735  -2.074   0.0381 *  
## neg_affect   0.12752    0.02238   5.698 1.22e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 387.81  on 279  degrees of freedom
## Residual deviance: 331.40  on 277  degrees of freedom
## AIC: 337.4
## 
## Number of Fisher Scoring iterations: 3

10.1 The decision boundary

The model predicts \(p = 0.5\) (log-odds \(= 0\)) along a straight line. Solving \(B_0 + B_1 \times \text{self-esteem} + B_2 \times \text{negative affect} = 0\) for negative affect gives:

bC <- coef(modC)
c(intercept = -bC[[1]] / bC[[3]], slope = -bC[[2]] / bC[[3]])
## intercept     slope 
## 10.049567  0.444689

10.2 The confusion matrix

predicted <- as.integer(predict(modC, type = "response") > 0.5)
cm <- confusionMatrix(factor(predicted), factor(d$high_depression), positive = "1")
cm
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 108  50
##          1  37  85
##                                          
##                Accuracy : 0.6893         
##                  95% CI : (0.6315, 0.743)
##     No Information Rate : 0.5179         
##     P-Value [Acc > NIR] : 4.201e-09      
##                                          
##                   Kappa : 0.3757         
##                                          
##  Mcnemar's Test P-Value : 0.1983         
##                                          
##             Sensitivity : 0.6296         
##             Specificity : 0.7448         
##          Pos Pred Value : 0.6967         
##          Neg Pred Value : 0.6835         
##              Prevalence : 0.4821         
##          Detection Rate : 0.3036         
##    Detection Prevalence : 0.4357         
##       Balanced Accuracy : 0.6872         
##                                          
##        'Positive' Class : 1              
## 

Rows are the model’s guesses, columns are the truth: 108 people were correctly cleared, 85 were correctly flagged, 37 were false alarms and 50 were missed. The accuracy is (108 + 85) / 280 = 0.689, clearly above the no-information rate of 0.518 (always guessing the larger group). These numbers are measured on the same 280 people the model was fitted to, so they are somewhat optimistic for new people.

11 Coefficient plot and marginal predictions

11.1 Coefficient plot

Each dot is a coefficient from Model 3 on the log-odds scale, with its 95% confidence interval. An interval that crosses 0 means no significant effect.

modelplot(mod3) + geom_vline(xintercept = 0, linetype = "dotted")

11.2 Marginal predictions

The predicted probability of high depression across the observed range of self-esteem, with a 95% confidence band: first for Model 1, then separately for men and women (Models 2 and 3). The band is widest where there are fewest people.

plot_predictions(mod1, condition = list("self_esteem")) +
  labs(y = "P(high depression)")

plot_predictions(mod2, condition = list("self_esteem", "gender")) +
  labs(y = "P(high depression)", colour = "gender\n(0 = man,\n1 = woman)", fill = "gender\n(0 = man,\n1 = woman)")

plot_predictions(mod3, condition = list("c_self_esteem", "gender")) +
  labs(y = "P(high depression)", colour = "gender\n(0 = man,\n1 = woman)", fill = "gender\n(0 = man,\n1 = woman)")

12 Session info

sessionInfo()
## R version 4.4.2 (2024-10-31)
## Platform: aarch64-apple-darwin20
## Running under: macOS 27.0
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: Europe/Zurich
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] knitr_1.51             marginaleffects_0.32.0 modelsummary_2.6.0     caret_7.0-1            lattice_0.22-9        
##  [6] lmtest_0.9-40          zoo_1.8-15             lubridate_1.9.5        forcats_1.0.1          stringr_1.6.0         
## [11] dplyr_1.2.0            purrr_1.2.1            readr_2.2.0            tidyr_1.3.2            tibble_3.3.1          
## [16] ggplot2_4.0.2          tidyverse_2.0.0        pacman_0.5.1          
## 
## loaded via a namespace (and not attached):
##  [1] pROC_1.19.0.1        sandwich_3.1-1       rlang_1.1.7          magrittr_2.0.4       multcomp_1.4-30     
##  [6] otel_0.2.0           e1071_1.7-17         compiler_4.4.2       vctrs_0.7.2          reshape2_1.4.5      
## [11] pkgconfig_2.0.3      fastmap_1.2.0        backports_1.5.0      labeling_0.4.3       effectsize_1.0.2    
## [16] rmarkdown_2.30       prodlim_2026.03.11   tzdb_0.5.0           xfun_0.57            cachem_1.1.0        
## [21] jsonlite_2.0.0       recipes_1.3.1        parallel_4.4.2       R6_2.6.1             bslib_0.10.0        
## [26] tables_0.9.33        stringi_1.8.7        RColorBrewer_1.1-3   parallelly_1.46.1    rpart_4.1.24        
## [31] jquerylib_0.1.4      estimability_1.5.1   Rcpp_1.1.1           iterators_1.0.14     future.apply_1.20.2 
## [36] parameters_0.28.3    Matrix_1.7-5         splines_4.4.2        nnet_7.3-20          timechange_0.4.0    
## [41] tidyselect_1.2.1     rstudioapi_0.18.0    yaml_2.3.12          timeDate_4052.112    codetools_0.2-20    
## [46] listenv_0.10.1       plyr_1.8.9           withr_3.0.2          bayestestR_0.17.0    S7_0.2.1            
## [51] coda_0.19-4.1        evaluate_1.0.5       future_1.70.0        survival_3.8-6       proxy_0.4-29        
## [56] pillar_1.11.1        checkmate_2.3.4      foreach_1.5.2        stats4_4.4.2         insight_1.4.6       
## [61] generics_0.1.4       hms_1.1.4            scales_1.4.0         xtable_1.8-8         globals_0.19.1      
## [66] class_7.3-23         glue_1.8.0           emmeans_2.0.2        tools_4.4.2          data.table_1.18.2.1 
## [71] ModelMetrics_1.2.2.2 gower_1.0.2          mvtnorm_1.3-6        grid_4.4.2           datawizard_1.3.0    
## [76] ipred_0.9-15         nlme_3.1-168         performance_0.16.0   cli_3.6.5            lava_1.8.2          
## [81] gtable_0.3.6         sass_0.4.10          digest_0.6.39        TH.data_1.1-5        farver_2.1.2        
## [86] htmltools_0.5.9      lifecycle_1.0.5      hardhat_1.4.2        MASS_7.3-65