[R] Omnibus test for main effects in the faceofaninteraction containing the main effects.

Daniel Malter daniel at umd.edu
Tue Sep 8 06:36:31 CEST 2009


John, as I wrote in the post sciptum, an anova on ML (but not REML) fitted
models seems permissible (Faraway 2006, "Extending the linear model with R",
p. 158). I am certainly not an expert on this and there are better "sources"
of information on why and when (e.g., Deepayan Sarkar, Julian Faraway, Doug
Bates, and certainly many others).  

As for your second problem, I think my statement remains valid, too. Because
this may yield biased estimates, you should not allow for such a model.

set.seed(1750)
x1=rbinom(120,2,c(1/3,1/2))
x2=rbinom(120,1,1/2)
e=rnorm(120,0,0.5)
x1=factor(x1)
x2=factor(x2)
y=1*(x1==1)+(-1)*(x1==2)+1.5*(x2==1)+(x1==1)*(x2==1)*(-2)+(x1==2)*(x2==1)*(3
)+e

reg=lm(y~x1+I((x1==1)*(x2==1))+I((x1==2)*(x2==1)))
summary(reg) 
##Gives biased estimates

reg1=lm(y~x1+x2+x1:x2)
summary(reg1)
##Gives unbiased estimates

Therefore, I would not use the model with interaction effects if not all
direct effects are included. And in fact, I never see this done in analyses
in my field (business/economics).

Daniel

-------------------------
cuncta stricte discussurus
-------------------------

-----Ursprüngliche Nachricht-----
Von: John Sorkin [mailto:Jsorkin at grecc.umaryland.edu] 
Gesendet: Monday, September 07, 2009 9:42 PM
An: r-help at r-project.org; daniel at umd.edu
Betreff: Re: [R] Omnibus test for main effects in the faceofaninteraction
containing the main effects.

Daniel,
When Group is entered as a factor, and the factor has two levels, the ANOVA
table gives a p value for each level of the factor. What I am looking for is
the omnibus p value for the factor, i.e. the test that the factor (with all
its levels) improves the prediction of the outcome.

You are correct that normally one could rely on the fact that the model
Post<-Time+as.factor(Group)+as.factor(Group)*Time
contains the model
Post<-Time+as.factor(Group)
and compare the two models using anova(model1,model2). However, my model is
has a random effect, the comparison is not so easy. The REML comparions of
nested random effects models is not valid when the fixed effects are not the
same in the models, which is the essence of the problem in my case. 

In addition to the REML problem if one wants to perform an omnibus test for
Group, one would want to compare nested models, one containing Group, and
the other not containing group. This would suggest comparing
Post<-Time+      as.factor(Group)*Time to
Post<-Time+Group+as.factor(Group)*Time
The quandry here is whether one should or not "allow" the first model as it
is poorly specified - one term of the interaction, as.factor(Group)*Time,
as.factor(Group) does not appear as a main effect
- a no-no in model building. 
John


John Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
Baltimore VA Medical Center GRECC,
University of Maryland School of Medicine Claude D. Pepper OAIC, University
of Maryland Clinical Nutrition Research Unit, and Baltimore VA Center Stroke
of Excellence

University of Maryland School of Medicine Division of Gerontology Baltimore
VA Medical Center 10 North Greene Street GRECC (BT/18/GR) Baltimore, MD
21201-1524

(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)
jsorkin at grecc.umaryland.edu
>>> "Daniel Malter" <daniel at umd.edu> 09/07/09 9:23 PM >>>
John, your question is confusing. After reading it twice, I still cannot
figure out what exactly you want to compare.

Your model "a" is the unrestricted model, and model "b" is a restricted
version of model "a" (i.e., b is a hiearchically reduced version of a, or
put differently, all coefficients of b are in a with a having additional
coefficients). Thus, it is appropriate to compare the models (also called
nested models).

Comparing c with a and d with a is also appropriate for the same reason.
However, note that depedent on discipline, it may be highly unconventional
to fit an interaction without all direct effects of the interacted variables
(the reason for this being that you may get biased estimates).

What you might consider is:
1. Run an intercept only model
2. Run a model with group and time
3. Run a model with group, time, and the interaction

Then compare 2 to 1, and 3 to 2. This tells you whether including more
variables (hierarchically) makes your model better.

HTH,
Daniel

On a different note, if lme fits with "restricted maximum likelihood," I
think I remember that you cannot compare them. You have to fit them with
"maximum likelihood." I am pointing this out because lmer with restricted
maximum likelihood by standard, so lme might too.

-------------------------
cuncta stricte discussurus
-------------------------

-----Ursprüngliche Nachricht-----
Von: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
Im
Auftrag von John Sorkin
Gesendet: Monday, September 07, 2009 4:00 PM
An: r-help at r-project.org
Betreff: [R] Omnibus test for main effects in the face of aninteraction
containing the main effects.

R 2.9.1
Windows XP

UPDATE,
Even my first suggestion
anova(fita,fitb) is probably not appropriate as the fixed effects are
different in the two model, so I don't even know how to perform the ombnibus
test for the interaction!



I am fitting a random effects ANOVA with two factors Group which has two
levels and Time which has three levels:
 fita<-lme(Post~Time+factor(Group)+factor(Group)*Time,
random=~1|SS,data=blah$alldata)

I wantinteraction. I believe I can get the omnibus test for the interaction
by running the model:

fitb<-lme(Post~Time+factor(Group), random=~1|SS,data=blah$alldata) followed
by anova(fita,fitb).

How do I get the omnibus test for the main effects i.e. for Time and
factor(Group)? I could drop each from the model, i.e.
fitc<-lme(Post~          factor(Group)+factor(Group)*Time,
random=~1|SS,data=blah$alldata)
fitd<-lme(Post~Time+                        factor(Group)*Time,
random=~1|SS,data=blah$alldata)

and then run
anova(fita,fitc)
anova(fita,fitd)
but I don't like this option as it will have in interaction that contains a
factor that is not included in the model as a main effect. How then do I get
the omnibus test for Time and factor(Group)?

Thanks
John




John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology Baltimore
VA Medical Center 10 North Greene Street GRECC (BT/18/GR) Baltimore, MD
21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)

Confidentiality Statement:
This email message, including any attachments, is for th...{{dropped:17}}




More information about the R-help mailing list