1
$\begingroup$

There is a similar, if not duplicate, question from 9 years ago on this but it never got answered so I'm hoping to get some clarity on this.

Let's say we have N treatments, such as medicines or fertilizers, that can be applied at varying intensities (continuous) and we want to fit an ANCOVA model to make inference about the difference in response to increasing intensity among the N treatments. So we want to know if the slope is different between treatments. The issue is that we have some individuals with intensity 0, where it is not meaningful to assign those individuals to any one of the treatments.

I just want to confirm that the following approach in R is correct. I assigned the individuals with intensity 0 to a treatment of "none" or something similar, so we now have N+1 unique values in the treatment column. Then fit the ANCOVA model as normal and also do post-hoc contrasts comparing the slopes between all pairs of treatments, but exclude the "none" treatment from those contrasts. Does this properly take account of the information from the individuals where intensity is 0 and treatment is "none"?

Example with two treatments

exdat <- structure(list(treatment = c("none", "none", "none", "none", 
"none", "none", "trtA", "trtA", "trtA", "trtB", "trtB", "trtB", 
"trtA", "trtA", "trtA", "trtB", "trtB", "trtB", "trtA", "trtA", 
"trtA", "trtA", "trtB", "trtB", "trtB", "trtA", "trtA", "trtA", 
"trtA", "trtB", "trtB", "trtB"), y = c(1.887069649, 2.2721258855, 
1.9459101491, 1.8405496334, 2.2300144002, 2.2192034841, 2.9123506646, 
2.8094026954, 1.974081026, 1.9169226122, 2.1041341543, 1.960094784, 
2.7472709143, 2.9391619221, 2.6532419646, 2.4680995315, 2.0149030205, 
2.7850112422, 2.2925347571, 2.3513752572, 3.1945831323, 2.6532419646, 
3.3068867022, 2.9285235239, 2.5649493575, 2.9231615807, 3.0155349009, 
2.9391619221, 2.7536607124, 2.9444389792, 2.8033603809, 2.9069010598
), intensity = c(0L, 0L, 0L, 0L, 0L, 0L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 
7L, 7L, 7L, 7L, 7L)), row.names = c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 
8L, 9L, 11L, 12L, 13L, 15L, 16L, 17L, 19L, 20L, 21L, 23L, 24L, 
25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L), class = "data.frame")

library(emmeans)

model <- lm(y ~ intensity + treatment + intensity:treatment, data = exdat)
summary(model)

emt <- emtrends(model, ~ treatment, var = 'intensity', at = list(treatment = c('trtA', 'trtB')))
contrast(emt, 'pairwise')

In this case because there are only two treatments other than the "none" treatment, the contrast gives us the same t-ratio and p-value as we get from looking at the interaction coefficient in the model summary:

 contrast    estimate    SE df t.ratio p.value
 trtA - trtB   -0.203 0.079 27  -2.567  0.0161
$\endgroup$
1
  • 1
    $\begingroup$ Yes, I believe this yields a correct analysis. $\endgroup$
    – Russ Lenth
    Mar 2 at 20:47

1 Answer 1

1
$\begingroup$

The "None" group sounds like a control group.

Thus, you would have something like the following:

Control

Treatment 1 at intensity levels $1, 2, 3, 4, 5,\dots$

Treatment 2 at intensity levels $1, 2, 3, 4, 5,\dots$

Treatment 3 at intensity levels $1, 2, 3, 4, 5,\dots$

$\vdots$

Then the "covariate" for the ANCOVA model would be the intensity, and you would multiply the covariate by each treatment indicator variable (the interactions that indicate differing slopes). Since there is no intensity for the control group, I would omit the covariate on its own. For the groups above:

$$ Y_i = \beta_0 + \beta_1X_{treatment1} + \beta_2X_{treatment2} + \beta_3X_{treatment3} + \beta_4 X_{treatment1}X_{intensity}+ \beta_5 X_{treatment2}X_{intensity}+ \beta_6 X_{treatment3}X_{intensity} +\epsilon_i $$

A regular ANCOVA would have an $X_{intensity}$ on its own, but you correctly identify that not to be meaningful, so omit it, fit your regression model, and interpret the coefficients.

$\endgroup$
2
  • 1
    $\begingroup$ You might want to interact the three treatments with polynomials or splines to investigate nonlinear behavior, but that does go beyond that standard ANCOVA setup. $\endgroup$
    – Dave
    Mar 1 at 17:43
  • $\begingroup$ Yes, I do have a situation with nonlinearities and a common zero level but I didn't want to add irrelevant detail to the question. I should be able to extend this to that situation pretty easily. Thanks for the sanity check! $\endgroup$
    – qdread
    Mar 1 at 17:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.