Homework 4 Answers

BSTA 513/613

Author

Nicky Wakim

Modified

May 29, 2026

Questions Part 1

Question 1

In this problem, we will practice performing model diagnostics in a logistic regression model.

This question is taken from the Hosmer and Lemeshow textbook. The ICU study data set consists of a sample of 200 subjects who were part of a much larger study on survival of patients following admission to an adult intensive care unit (ICU). The dataset should be available in our shared folder. The major goal of this study was to develop a logistic regression model to predict the probability of survival to hospital discharge of these patients. In this question, the primary outcome variable is vital (survival) status at hospital discharge, STA. Clinicians associated with the study felt that a key determinant of survival was the patient’s age at admission, AGE. We will build to a multivariable logistic regression model while adjusting for cancer part of the present problem (CAN), CPR prior to ICU admission (CPR), infection probable at ICU admission (INF), and level of consciousness at ICU admission (LOC).

A code sheet for the variables to be considered is displayed in Table 1.5 below (from the Hosmer and Lemeshow textbook, pg. 23). We refer to this data set as the ICU data.

You will need to use some of the mutations implemented in HW 2, Q2, Part d.

We will use the following model: \[\text{logit}(\pi(\textbf{X}))=\beta_0 + \beta_1 \cdot I(CAN=\text{``Yes"}) + \beta_2 \cdot I(CPR=\text{``Yes"}) + \\ \beta_3 \cdot I(INF=\text{``Yes"})\]

icu = read_csv(here("data", "icu.csv"))
icu1 = icu %>% mutate(STA = as.factor(STA) %>% relevel(ref = "Lived"))
icu2 = icu1 %>% mutate(CAN = as.factor(CAN) %>% relevel(ref = "No"), 
                     CPR = as.factor(CPR) %>% relevel(ref = "No"), 
                     INF = as.factor(INF) %>% relevel(ref = "No"), 
                     LOC = as.factor(LOC) %>% 
                       relevel(ref = "No Coma or Deep Stupor"))
term estimate std.error statistic p.value conf.low conf.high
(Intercept) −1.924 0.284 −6.764 0.000 −2.521 −1.399
CANYes 0.093 0.614 0.151 0.880 −1.246 1.219
CPRYes 1.514 0.602 2.516 0.012 0.326 2.730
INFYes 0.807 0.371 2.176 0.030 0.084 1.547

Part a

How many potential covariate patterns does the regression equation have?

Not given

Part b

Plot the standardized Deviance residual by predicted probability. Do you notice any potential outliers? Explain your reasoning.

ggplot(dx_m1, aes(x = .fitted, y = .std.resid)) +
  geom_point() +
  geom_hline(yintercept = c(-3, 3), col = "red", linetype="dashed") +
  labs(x = "Predicted Probability", y = "Std. Deviance Residual")

Part c

Plot the Cook’s distance by predicted probability. Do you notice any potential influential points? Explain your reasoning.

Not given

Part d

Plot the change in coefficient estimate by predicted probability. Do you notice any influential points? Explain your reasoning.

Not given

Question 2

This question stems from an example from an online textbook by Dr. Ramzi W. Nahhas. The dataset for this problem includes a subset of individuals from the 2019 National Survey on Drug Use and Health (NSDUH). Overall, our study aims included investigating potential risk factors for lifetime heroin use. Lifetime heroin use is a binary outcome, which we regress on age at first use of alcohol (alc_agefirst), age with 6 categories (demog_age_cat6), and sex assigned at birth (demog_sex).

load(here("data", "nsduh2019_adult_sub_rmph.RData"))
nsduh = nsduh_adult_sub %>% 
  dplyr::select(her_lifetime, alc_agefirst, demog_age_cat6, demog_sex) %>% 
  drop_na()

Part a

Using the nsduh dataset from the above chunk of code, please run a regression model and present the model summary using lifetime heroin use as our outcome, and age at first use of alcohol, categorical age, and sex assigned at birth as covariates in our model. No need to write out your model, you just need to write the R code to run it.

term estimate std.error statistic p.value conf.low conf.high
(Intercept) −15.061 1,024.209 −0.015 0.988 NA 43.320
alc_agefirst −0.244 0.064 −3.821 0.000 −0.373 −0.121
demog_age_cat626-34 15.539 1,024.209 0.015 0.988 −42.842 NA
demog_age_cat635-49 15.447 1,024.209 0.015 0.988 −47.994 NA
demog_age_cat650-64 15.602 1,024.209 0.015 0.988 −42.779 NA
demog_age_cat665+ 15.363 1,024.209 0.015 0.988 −43.018 NA
demog_sexFemale −1.238 0.653 −1.897 0.058 −2.727 −0.076

Part b

Are we encountering a numerical problem with our regression? If yes, please name the numerical issue. What first clued you into that issue? Provide conclusive evidence of this numerical issue (with a contingency table), and explain which variable(s) are causing this problem.

Yes

Part c

What would you do to “fix” this numerical issue? Please apply your “fix” and rerun the regression

Not given

Question 3

In this problem, we will practice fitting and interpreting a log-binomial regression.

This question is taken from the Hosmer and Lemeshow textbook. The ICU study data set consists of a sample of 200 subjects who were part of a much larger study on survival of patients following admission to an adult intensive care unit (ICU). The dataset should be available in our shared folder. The major goal of this study was to develop a logistic regression model to predict the probability of survival to hospital discharge of these patients. In this question, the primary outcome variable is vital (survival) status at hospital discharge, STA. Clinicians associated with the study felt that a key determinant of survival was the patient’s age at admission, AGE. We will build to a multivariable logistic regression model while adjusting for cancer part of the present problem (CAN), CPR prior to ICU admission (CPR), infection probable at ICU admission (INF), and level of consciousness at ICU admission (LOC).

A code sheet for the variables to be considered is displayed in Table 1.5 below (from the Hosmer and Lemeshow textbook, pg. 23). We refer to this data set as the ICU data.

You will need to use some of the mutations implemented in HW 2, Q2, Part d.

Part a

Write down the population equation for the log-binomial regression model of STA on AGE, CAN, CPR, and INF. How many parameters does this model contain?

This model contains 5 parameters.

Part b

Try using glm() to obtain the maximum likelihood estimates of the parameters of the log-binomial regression model in Part a. Do you run into any issues using glm()? If you try logbin(), does it fix the issues? Explain why and what warnings logbin() gives you.

Hint: keep the glm() function in its own code chunk so you can add #| eval: false. glm() may throw an error, so we want to show the work for glm() even though it’ll break your qmd rendering.

glm() does not work. logbin() has warnings.

Part c

Using logbin(), obtain the maximum likelihood estimates of the parameters of the log-binomial regression model in Part a, but now take out age. Using these estimates, write down the equation with the fitted values.


Call:
logbin(formula = STA ~ CAN + CPR + INF, data = icu2)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.4845  -0.7570  -0.5297  -0.5197   2.0886  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -2.0332     0.2374  -8.566  < 2e-16 ***
CANYes       -0.1478     0.4826  -0.306 0.759330    
CPRYes        0.9859     0.2874   3.430 0.000603 ***
INFYes        0.6434     0.2910   2.211 0.027045 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

    Null deviance: 200.16  on 199  degrees of freedom
Residual deviance: 186.90  on 196  degrees of freedom

  AIC: 194.9 
AIC_c: 195.1 

Number of iterations: 21093 (best: 480)

Part d

Interpret the exponential of the coefficient (risk ratio) estimate for CPR.

Not given

Part e

We were not able to fit a model with age, but let’s just entertain a scenario here. Let’s say we fit the model in Part a and got a coefficient estimate of 0.29 with a 95% confidence interval of 0.23 to 0.35. Using the model is Part a, interpret the exponential of the coefficient (risk ratio) estimate for AGE.

Not given

Questions Part 2

Question 4

For each of the following outcomes, what type of regression would you use? Explain your answer.

Part a

Number of minutes of moderate-to-vigorous physical activity per week (range: 0 to 600+)

Which regression model is most appropriate?

  1. Linear regression
  2. Logistic regression
  3. Log-binomial regression
  4. Poisson regression
  5. Multinomial logistic regression
  1. Linear regression

Part b

Whether the participant meets CDC recommendations for weekly physical activity (Yes/No)

Which regression model is most appropriate?

  1. Linear regression
  2. Logistic regression
  3. Log-binomial regression
  4. Poisson regression
  5. Multinomial logistic regression

Not given

Part c

Number of workouts the person completed in the past week (values: 0, 1, 2, …, 14)

Which regression model is most appropriate?

  1. Linear regression
  2. Logistic regression
  3. Log-binomial regression
  4. Poisson regression
  5. Multinomial logistic regression

Not given

Part d

Self-reported activity level: Sedentary, Moderately active, Highly active

Which regression model is most appropriate?

  1. Linear regression
  2. Logistic regression
  3. Log-binomial regression
  4. Poisson regression
  5. Multinomial logistic regression

Not given