1
$\begingroup$

I'd like to test if the slopes of two linear regression models differ. However, the caveat is that one of the regressions fits a subset of the data, and the other fits the whole dataset. The two models use the same predictor variables and aren't nested. Likelihood ratio tests require the compared models to be nested and fitted to the same dataset size. AIC and BIC have been suggested to compare non-nested models but I'm not interested in which model fits best, just if the slopes differ significantly from one another.

Example (R):

library('lmtest')   
df <- data.frame(time=c(1,2,3,4,5,6,7,8,9,10), response=c(1,2,4,3,7,7,6,6,7,10))
lm1 <- lm(response ~ time, data=df)
lm2 <- lm(response ~ time, data=df[1:6,])
lrtest(lm1, lm2) # this approach does not work because the datasets are different sizes

What other approach might work here? Can I just look at the confidence intervals of the slopes and check if they overlap? Thank you for your help!

$\endgroup$

1 Answer 1

6
$\begingroup$

Welcome here.

A simple solution is to define a new binary variables which takes the value $1$ for the subset of the data and $0$ for all other observations not included in the subset. Then you can interact the slope you are interested in with the binary variable to see whether you find a statistical significant differences in the slope between the subset and the whole dataset:

df$sample = 0 # 'sample' takes the value 0 for the whole dataset
df[1:6, "sample"] = 1 # 'sample' takes the value 1 for the subset
linear_regression_model_with_interaction = lm(response ~ time + time:sample, data = df)
summary(linear_regression_model_with_interaction)

The results look as following. In the example data you presented, we find statistical differences in the slope estimated on the 10 percent level (ignoring potential heteroscedasticity etc. here, or whether it would be wise to also include a main effect for the 'sample' variable):

Call:
lm(formula = response ~ time + time:sample, data = df)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.61401 -0.63063  0.03804  0.58723  1.30236 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -0.5893     0.9292  -0.634 0.546128    
time          0.9287     0.1230   7.548 0.000132 ***
time:sample   0.3721     0.1632   2.281 0.056591 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.027 on 7 degrees of freedom
Multiple R-squared:  0.8915,    Adjusted R-squared:  0.8606 
F-statistic: 28.77 on 2 and 7 DF,  p-value: 0.0004201

PS: Potential heteroscedasticity can be taken into account via heteroscedasticity-robust standard error:

library('lmtest')
library('sandwich')
coeftest(linear_regression_model_with_interaction, vcov = vcovHC(linear_regression_model_with_interaction, type = "HC0"))

But keep in mind that robust standard errors might be not valid with small samples, this might be the reason why the p-value drops here (of course heteroscedasticity can bias your variance estimation in both directions but usually the standard error get larger).

t test of coefficients:

             Estimate Std. Error t value  Pr(>|t|)    
(Intercept) -0.589293   0.357539 -1.6482   0.14330    
time         0.928694   0.079517 11.6792 7.625e-06 ***
time:sample  0.372132   0.109051  3.4124   0.01125 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1  
$\endgroup$
2
  • 2
    $\begingroup$ You beat me by a minute or so, so I'll delete my answer, which is essentially the same! $\endgroup$
    – jbowman
    Dec 8, 2021 at 18:59
  • $\begingroup$ Thank you! I'm glad this has such a simple solution. Could someone please help me understand how this works conceptually, though? Or point me towards a resource that explains it? $\endgroup$
    – cas
    Dec 9, 2021 at 15:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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