{"id":768,"date":"2019-01-21T12:41:49","date_gmt":"2019-01-21T17:41:49","guid":{"rendered":"http:\/\/www.clayford.net\/statistics\/?p=768"},"modified":"2020-09-13T10:17:54","modified_gmt":"2020-09-13T14:17:54","slug":"most-useful-tests-for-an-ancova-model","status":"publish","type":"post","link":"https:\/\/www.clayford.net\/statistics\/most-useful-tests-for-an-ancova-model\/","title":{"rendered":"Most useful tests for an ANCOVA model"},"content":{"rendered":"

In his book Regression Modeling Strategies, 2nd Ed<\/em>, Frank Harrell provides a list of what he calls the \u201cmost useful tests\u201d for a 2-level factor \\(\\times\\) numeric model (Table 2.2, p. 19). This is often called an Analysis of Covariance, or ANCOVA. The basic idea is we have a numeric response with substantial variability and we seek to understand the variability by modeling the mean of the response as a function of a categorical variable and a numeric variable.<\/p>\n

Let’s simulate some data for such a model and then see how we can use R to carry out these tests.<\/p>\n

n <- 400\nset.seed(1)\nsex <- factor(sample(x = c(\"f\", \"m\"), size = n, replace = TRUE))\nage <- round(runif(n = n, min = 18, max = 65))\ny <- 1 + 0.8*age + 0.4*(sex == \"m\") - 0.7*age*(sex == \"m\") + rnorm(n, mean = 0, sd = 8)\ndat <- data.frame(y, age, sex)\n<\/code><\/pre>\n
The data contain a numeric response, y<\/code>, that is a function of age<\/code> and sex<\/code>. I set the \u201ctrue\u201d coefficient values to 1, 0.8, 0.4, and -0.7. They correspond to \\(\\beta_0\\) through \\(\\beta_3\\) in the following model:<\/p>\n
\\[y = \\beta_0 + \\beta_1 age + \\beta_2 sex + \\beta_3 age \\times sex\\]<\/p>\n
In addition the error component is a Normal distribution with a standard deviation of 8.<\/p>\n
Now let’s model the data and see how close we get to recovering the true parameter values.<\/p>\n
mod <- lm(y ~ age * sex, dat)\nsummary(mod)\n<\/code><\/pre>\n
## \n## Call:\n## lm(formula = y ~ age * sex, data = dat)\n## \n## Residuals:\n##      Min       1Q   Median       3Q      Max \n## -23.8986  -5.8552  -0.2503   6.0507  30.6188 \n## \n## Coefficients:\n##             Estimate Std. Error t value Pr(>|t|)    \n## (Intercept)  0.27268    1.93776   0.141    0.888    \n## age          0.79781    0.04316  18.484   <2e-16 ***\n## sexm         2.07143    2.84931   0.727    0.468    \n## age:sexm    -0.72702    0.06462 -11.251   <2e-16 ***\n## ---\n## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\n## \n## Residual standard error: 8.661 on 396 degrees of freedom\n## Multiple R-squared:  0.7874, Adjusted R-squared:  0.7858 \n## F-statistic:   489 on 3 and 396 DF,  p-value: < 2.2e-16\n<\/pre>\n
While the coefficient estimates for age and the age \\(\\times\\) sex interaction are pretty close to the true values, the same cannot be said for the intercept and sex coefficients. The residual standard error of 8.661 is close to the true value of 8.<\/p>\n
We can see in the summary output of the model that four hypothesis tests, one for each coefficient, are carried out for us. Each are testing if the coefficient is equal to 0. Of those four, only one qualifies as one of the most useful tests: the last one for age:sexm<\/code>. This tests if the effect of age is independent of sex and vice versa. Stated two other ways, it tests if age and sex are additive, or if the age effect is the same for both sexes. To get a better understanding of what we’re testing, let’s plot the data with fitted age slopes for each sex.<\/p>\n
library(ggplot2)\nggplot(dat, aes(x = age, y = y, color = sex)) +\n  geom_point() +\n  geom_smooth(method=\"lm\")\n<\/code><\/pre>\n
\"plot<\/p>\n
Visually it appears the effect of age is not<\/em> independent of sex. It seems more pronounced for females. Is this effect real or maybe due to chance? The hypothesis test in the summary output for age:sexm<\/code> evaluates this. Obviously the effect seems very real. We are not likely to see such a difference in slopes this large if there truly was no difference. It does appear the effect of age is different for each sex. The estimate of -0.72 estimates the difference in slopes (or age effect) for the males and females.<\/p>\n
The other three hypothesis tests are not very useful.<\/p>\n