In a study at the Ohio State University Medical Center, doctors investigated the association between ICU acquired paresis (ICUAP) (or muscular weakness) and subject’s admission SOFA score dichotomized at 11 (sofa11) (0 for less than 11 and 1 for 11 and greater).
Specifically the investigators wanted to determine if handgrip strength (max_grip) could be used as a surrogate for ICUAP. The investigators believed that the number of days the patient was mechanically ventilated (MVdays) could affect the relationship between weakness and dichotomized SOFA.
Download and read data
url <- "https://d396qusza40orc.cloudfront.net/appliedregression/Homeworks/Acquired_weakness.csv"
download.file(url, "icuap.csv")
icu <- read.csv("icuap.csv")
names(icu) <- tolower(names(icu))
Using a two-sample t-test, determine whether the two SOFA Score groups differ significantly with respect to ICUAP (as measured by max_grip).
t.test(max_grip~sofa11, icu, var.equal = TRUE)
##
## Two Sample t-test
##
## data: max_grip by sofa11
## t = 2.205, df = 135, p-value = 0.02915
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 0.924972 17.021695
## sample estimates:
## mean in group 0 mean in group 1
## 15.140000 6.166667
Differ significantly on \(\alpha = .05\)
Perform the same analysis from Exercise One (i.e., 2 sample t-test) using a regression analysis.
summary(lm(max_grip ~ sofa11, icu))
##
## Call:
## lm(formula = max_grip ~ sofa11, data = icu)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.14 -7.14 -1.14 4.86 68.86
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 15.140 1.204 12.570 <2e-16 ***
## sofa11 -8.973 4.070 -2.205 0.0291 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 13.47 on 135 degrees of freedom
## Multiple R-squared: 0.03476, Adjusted R-squared: 0.02761
## F-statistic: 4.862 on 1 and 135 DF, p-value: 0.02915
Same as before. Overall F-statistic is equal to \(t^2\) and has the same p-value on 1 and \((n-k-1)\) degrees of freedom.
Using a single model determine the regression of max_grip on MVdays for each of the two SOFA score groups.
We have to fit a model with two terms and an interaction. R kindly does the dummy coding for us. Also note the syntax.
summary(lm(max_grip ~ mvdays * sofa11, icu))
##
## Call:
## lm(formula = max_grip ~ mvdays * sofa11, data = icu)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.715 -7.022 -2.246 4.691 65.591
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.1842 1.9915 10.638 < 2e-16 ***
## mvdays -0.6938 0.1865 -3.721 0.000292 ***
## sofa11 -10.3081 6.9788 -1.477 0.142023
## mvdays:sofa11 0.4084 0.3851 1.060 0.290905
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.88 on 133 degrees of freedom
## Multiple R-squared: 0.13, Adjusted R-squared: 0.1104
## F-statistic: 6.625 on 3 and 133 DF, p-value: 0.0003318
The interaction term can be excluded from the model.
Test whether the two SOFA cohorts differ significantly, controlling for mechanically ventilated days.
summary(lm(max_grip ~ mvdays + sofa11, icu))
##
## Call:
## lm(formula = max_grip ~ mvdays + sofa11, data = icu)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.360 -7.193 -1.958 3.611 66.042
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.3503 1.8304 11.118 < 2e-16 ***
## mvdays -0.5981 0.1632 -3.664 0.000357 ***
## sofa11 -4.3157 4.0965 -1.053 0.294009
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.89 on 134 degrees of freedom
## Multiple R-squared: 0.1227, Adjusted R-squared: 0.1096
## F-statistic: 9.367 on 2 and 134 DF, p-value: 0.0001557
\(Pr(>|t|) = .294\), so they do not.
Test whether the two SOFA cohorts have the same slope. Do they have the same intercept? Are they coincident?
- \(\beta_3\) is equal to 0, so they have the same slope
- \(\beta_2\) is 0, so they have the same intercept, which means they are coincident.
Might as well be illustrated:
library(ggplot2)
ggplot(icu, aes(mvdays, max_grip, col = factor(sofa11))) + geom_point() +
geom_smooth(method = "lm", se=F) + theme_classic()
By the way, this scatter look awful (not really important for our homework after all).