[R-sig-eco] Output for interactions in models that do not include all main effects

Bob O'Hara bohara at senckenberg.de
Wed Apr 4 10:31:41 CEST 2012


On 04/03/2012 11:31 PM, Kristen Gorman wrote:
> Dear all,
> I have R code to run AIC including multi-model inference. I am running into a problem in calling the output from models where both parameters in an interaction are not included as main effects.
Why would you want to do that? Why would you (for example) expect the 
average of the Rlipid slope to be zero if the slope varies with the 
value of RFGinit? Does this make sense?

(this is the sort of thing that makes statisticians splutter into their 
tea when they see someone do it: it rarely makes sense. Well, unless you 
have nested effects - which you don't have here- where the interaction 
is the nested effect)

if you respect marginality, there won't be a problem because the main 
effect is always included. If you really want to include interactions 
without main effects, you can either write the formula "by hand", using 
paste():

something=Rlipid
form = paste("Slipid ~ ", something, " + RFGinit:", something, sep="")
lm(form, data = DataSet)

and then work out how to get the order. Or you could try using update():

mod1 = lm(formula = Slipid ~ RFGinit*Rlipid, data = DataSet)
mod2=update(mod1, . ~ . -RFGinit)

HTH

Bob

> In R, the interaction will be called depending on the parameter that was used as the only main effect in the model. So, I end up generating 2 different interactions (e.g., Rlipid:RFGinit vs RFGinit:Rlipid) that are actually the same. This becomes a problem in the remaining R code that requires weighted and summed values for the parameter and SE estimates. Thus, I would like to call the interaction consistently across models. See the following code:
>
> --
> lm(formula = Slipid ~ Rlipid + RFGinit:Rlipid, data = DataSet)
>
> Residuals:
>      Min      1Q  Median      3Q     Max
> -74.075 -19.047   7.233  20.445  45.391
>
> Coefficients:
>                  Estimate Std. Error t value Pr(>|t|)
> (Intercept)    120.33847    5.30405  22.688<2e-16 ***
> Rlipid           0.30493    0.23615   1.291    0.202
> Rlipid:RFGinit  -0.02099    0.01773  -1.184    0.241
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Residual standard error: 30.88 on 60 degrees of freedom
> Multiple R-squared: 0.02721,        Adjusted R-squared: -0.005221
> F-statistic: 0.839 on 2 and 60 DF,  p-value: 0.4372
>
>
> lm(formula = Slipid ~ RFGinit + Rlipid:RFGinit, data = DataSet)
>
> Residuals:
>     Min     1Q Median     3Q    Max
> -76.35 -21.63   7.09  22.46  45.71
>
> Coefficients:
>                   Estimate Std. Error t value Pr(>|t|)
> (Intercept)    131.028546   8.717104  15.031<2e-16 ***
> RFGinit         -0.933483   0.742083  -1.258    0.213
> RFGinit:Rlipid   0.003926   0.009283   0.423    0.674
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Residual standard error: 30.9 on 60 degrees of freedom
> Multiple R-squared: 0.02586,        Adjusted R-squared: -0.00661
> F-statistic: 0.7964 on 2 and 60 DF,  p-value: 0.4556
> --
>
>
> Is there a way to tell R to call the interaction based on alphabetical order of the 2 interaction terms and not based on the term that was used as a main effect?
>
> Thanks very much for any insight.
>
> Kristen Gorman
>
> _______________________________________________
> R-sig-ecology mailing list
> R-sig-ecology at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology


-- 
Bob O'Hara

Biodiversity and Climate Research Centre
Senckenberganlage 25
D-60325 Frankfurt am Main,
Germany

Tel: +49 69 798 40226
Mobile: +49 1515 888 5440
WWW:   http://www.bik-f.de/root/index.php?page_id=219
Blog: http://blogs.nature.com/boboh
Journal of Negative Results - EEB: www.jnr-eeb.org



More information about the R-sig-ecology mailing list