[R] simple main effect.

439 views
Skip to first unread message

Or Duek

unread,
Dec 19, 2009, 12:27:35 PM12/19/09
to R-h...@r-project.org
Hi, I'm a bit new to R and I would like to know how can I compare simple
main effects when using the aov function.
I'm doing a mixed model ANOVA with two between subjects variables and one
within.
When I get an interaction of two of the variables I don't know how to check
for simple main effect of that interaction (A at B1 and A at B2 for
example).
The aov function is very simple but for some reason I can't find how to do
this.
Thank you very much.
Or Duek.

[[alternative HTML version deleted]]

______________________________________________
R-h...@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Tal Galili

unread,
Dec 20, 2009, 3:00:34 AM12/20/09
to Or Duek, R-h...@r-project.org
Try going through this:
http://www.personality-project.org/r/r.anova.html

<http://www.personality-project.org/r/r.anova.html>


----------------Contact
Details:-------------------------------------------------------
Contact me: Tal.G...@gmail.com | 972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com/ (English)
----------------------------------------------------------------------------------------------

Or Duek

unread,
Dec 20, 2009, 4:56:03 AM12/20/09
to Tal Galili, R-h...@r-project.org
Thanks.
But there is no simple main effect there.

Tal Galili

unread,
Dec 20, 2009, 8:35:03 AM12/20/09
to Or Duek, R-h...@r-project.org
Hi Or,

Maybe I didn't understand you.

Are you asking "how can I read the output of a fitted (complex within and
between) model for finding the simple main effect" ?


Tal

----------------Contact
Details:-------------------------------------------------------
Contact me: Tal.G...@gmail.com | 972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com/ (English)
----------------------------------------------------------------------------------------------

Or Duek

unread,
Dec 20, 2009, 8:36:52 AM12/20/09
to Tal Galili, R-h...@r-project.org
I don't think so.
I'm asking how can I see/analyse the simple main effect. I don't think it
shows in the summary report of the aov() function.

Tal Galili

unread,
Dec 20, 2009, 9:18:22 AM12/20/09
to Or Duek, R-h...@r-project.org
Could you please write the aov formula you are using ?

----------------Contact
Details:-------------------------------------------------------
Contact me: Tal.G...@gmail.com | 972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com/ (English)
----------------------------------------------------------------------------------------------

Or Duek

unread,
Dec 20, 2009, 9:30:47 AM12/20/09
to Tal Galili, R-h...@r-project.org
No problem.
If I have a mixed model with A as within subject (A is exposure, has 2
levels) and drug(2 levels between subject) and strain(2 levels between
subject).
I reshape the data to fit the R aov() function.
thus instead of looking like this:
subject drug strain exposure_1 exposure2
1 1 1 34 25
2 2 1 26 22
etc.

it looks like this:
subject drug strain exposure dependent
1 1 1 1 34
1 1 1 2 25
etc.

and the I do:
aov(dependent~(exposure*strain*drug) + Error(subject/exposure) +
(drug*strain), data)

Thank you very much for your help.
Or.

Tal Galili

unread,
Dec 20, 2009, 10:52:23 AM12/20/09
to Or Duek, R-h...@r-project.org
Thanks Or,
So I am failing to understand,
When you put that aov expression into the "summary()" Why isn't what you are
getting what you need ?

Tal


----------------Contact
Details:-------------------------------------------------------
Contact me: Tal.G...@gmail.com | 972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com/ (English)
----------------------------------------------------------------------------------------------

Tal Galili

unread,
Dec 20, 2009, 10:52:53 AM12/20/09
to Or Duek, R-h...@r-project.org
Also, do you have any missing data ?

Tal

----------------Contact
Details:-------------------------------------------------------
Contact me: Tal.G...@gmail.com | 972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com/ (English)
----------------------------------------------------------------------------------------------

Or Duek

unread,
Dec 20, 2009, 10:56:10 AM12/20/09
to Tal Galili, R-h...@r-project.org
I don't have missing data.
about what I need.
Lets say the drug*strain interaction is significant - now I want to check
for drug under the levels of strain - compare drug 1 and 2 only on strain 1
and then only on strain 2.
Or I'd like to compare the strains under levels of exposure.
This is the kind of data I fail to see in summary() but it is important to
understand the interactions.
thank you.

S Devriese

unread,
Dec 20, 2009, 11:32:54 AM12/20/09
to Or Duek, R-h...@r-project.org
On 12/20/2009 04:56 PM, Or Duek wrote:
> I don't have missing data.
> about what I need.
> Lets say the drug*strain interaction is significant - now I want to check
> for drug under the levels of strain - compare drug 1 and 2 only on strain 1
> and then only on strain 2.
> Or I'd like to compare the strains under levels of exposure.
> This is the kind of data I fail to see in summary() but it is important to
> understand the interactions.
> thank you.
>

Do you main pairwise multiple comparison tests like Tukey Honest
Significant Difference tests? Then you could use TukeyHSD in the stats
package or see the DTK package (Dunnett-Tukey-Kramer Pairwise Multiple
Comparison Test Adjusted for Unequal Variances and Unequal Sample Sizes)

Stephan

Or Duek

unread,
Dec 20, 2009, 11:35:10 AM12/20/09
to S Devriese, R-h...@r-project.org
For some reasion I wasn't able to use TukeyHSD - I think because I need to
set the different levels under a second variable.
Tukey only helps me when I have more than 2 levels of same variable.
Thanl you.

On Sun, Dec 20, 2009 at 6:32 PM, S Devriese <sdmai...@gmail.com> wrote:

> On 12/20/2009 04:56 PM, Or Duek wrote:
> > I don't have missing data.
> > about what I need.
> > Lets say the drug*strain interaction is significant - now I want to check
> > for drug under the levels of strain - compare drug 1 and 2 only on strain
> 1
> > and then only on strain 2.
> > Or I'd like to compare the strains under levels of exposure.
> > This is the kind of data I fail to see in summary() but it is important
> to
> > understand the interactions.
> > thank you.
> >
>
> Do you main pairwise multiple comparison tests like Tukey Honest
> Significant Difference tests? Then you could use TukeyHSD in the stats
> package or see the DTK package (Dunnett-Tukey-Kramer Pairwise Multiple
> Comparison Test Adjusted for Unequal Variances and Unequal Sample Sizes)
>
> Stephan
>

[[alternative HTML version deleted]]

Ista Zahn

unread,
Dec 20, 2009, 12:09:22 PM12/20/09
to Or Duek, R-h...@r-project.org
Hi Or,
I understand your question (and am not sure what the confusion is
about actually). Just like you said, you want to know the effect of A
at B=1, and the effect of A at B=2. In this case, you want to know if
drug has a significant effect for those with strain one, and whether
drug has a significant effect for those with strain two. So the good
part is, I understand the question.


The bad part is that unfortunately I don't know how to do it with
aov(). Also, I'm not a stats guru and could easily be wrong about the
following advice. You've been warned!

If you didn't have a mixed model, I would tell you to run the model
using lm() twice, setting drug = 1 as the reference group the first
time, and drug = 2 as the reference group the second time, but this
won't work in your case. The best I can offer is a suggestion to run
the model separately for each level of drug. Note that you can and
should use the error term from the overall model though -- you will
have to do this "by hand" (just divide MS_effect from the subset
model by MS_error from the full model, and evaluate using df_error
from the overall model). So basically, I'm suggesting that you do

Data.drug1 <- subset(Data, drug == "1")
aov.model.drug.1 <- aov(dependent~(exposure*strain) +
Error(subject/exposure) + (strain), data=Data.drug1)
summary(aov.model.drug.1)

Data.drug2 <- subset(Data, drug == "2")
aov.model.drug.2 <- aov(dependent~(exposure*strain) +
Error(subject/exposure) + (strain), data=Data.drug2)
summary(aov.model.drug.2)

Good luck!

-Ista

--
Ista Zahn
Graduate student
University of Rochester
Department of Clinical and Social Psychology
http://yourpsyche.org

Or Duek

unread,
Dec 20, 2009, 12:28:43 PM12/20/09
to Ista Zahn, R-h...@r-project.org
Thank you very much Ista.
Can you please be a bit more specific as to using the overall error. I don't
know how to actually do it in R.
thank you,
Or.

Ista Zahn

unread,
Dec 20, 2009, 12:45:57 PM12/20/09
to Or Duek, R-h...@r-project.org
Well, the F value is just MS_effect/MS_error. So you run the original
model, run the subset models, but disregard the F values. To calculate
the correct F values, divide MS_effect from the subset model by the
MS_error term from the original (full) model. It's not really a matter
of learning how to do it in R -- you're just doing simple division.

-Ista

Richard M. Heiberger

unread,
Dec 20, 2009, 5:16:24 PM12/20/09
to Or Duek, R-h...@r-project.org
For simple effects in the presence of interaction there are several
options included in the HH package. If you don't already have the HH
package, you can get it with
install.packages("HH")

Graphically, you can plot them with the function
interaction2wt(..., simple=TRUE)
See the examples in
?HH::interaction2wt

For tests on the simple effect of A conditional on a level of B, you
can use the model formula B/A and look at the partition of the sums of
squares using the split= argument
summary(mymodel.aov, split=<put your details here>)

For multiple comparisons from designs with Error() terms, you need to
specify the same sums of squares with an equivalent formula that doesn't
use the Error() function. See the maiz example in
?HH::MMC
Read the example all the way to the end of the help file.

Rich

Or Duek

unread,
Feb 27, 2010, 9:18:48 AM2/27/10
to Richard M. Heiberger, R-h...@r-project.org
I am very new to R and thus find those examples a bit confusing although I
believe the solution to my problems lies there.
Lets take for example an experiment in which I had two between subject
variables - Strain and treatment, and one within - exposure. all the
variables had 2 levels each.

I found an interaction between exposure and Strain and I want to compare
Strain A and B under every exposure (first and second).
The general model was with that function:
aov(duration~(Strain*exposure*treatment)+Error(subject/exposure),data)

in summary(aovmodel) there was a significant interaction between exposure
and strain.
how (using those HH packages) can I compare Strains under the conditions of
exposure?


BTW - I don't have to use aov (although its seems to be the simplest one).

Thank you very much.


On Mon, Dec 21, 2009 at 12:16 AM, Richard M. Heiberger <r...@temple.edu>wrote:

> For simple effects in the presence of interaction there are several
> options included in the HH package. If you don't already have the HH
> package, you can get it with
> install.packages("HH")
>
> Graphically, you can plot them with the function
> interaction2wt(..., simple=TRUE)
> See the examples in
> ?HH::interaction2wt
>
> For tests on the simple effect of A conditional on a level of B, you
> can use the model formula B/A and look at the partition of the sums of
> squares using the split= argument
> summary(mymodel.aov, split=<put your details here>)
>
> For multiple comparisons from designs with Error() terms, you need to
> specify the same sums of squares with an equivalent formula that doesn't
> use the Error() function. See the maiz example in
> ?HH::MMC
> Read the example all the way to the end of the help file.
>
> Rich
>
>

[[alternative HTML version deleted]]

DrorD

unread,
Feb 27, 2010, 11:11:53 AM2/27/10
to r-h...@r-project.org
I tried to implement Ista's procedure and would like to provide it as
a working example, with the intention to get feedback from the R
community:

The data contains three variables:
One dependent var: t.total
and two independent vars: group (between: D2C2, C2D2) and present.type
(within: C2, D2).

# First I do the "overall" ANOVA:
m.full=aov(t.total ~ group * present.type + Error(subj/present.type),
data=dat.net)
summary(m.full)

Error: subj
Df Sum Sq Mean Sq F value Pr(>F)
group 1 1430 1430 0.4224 0.528
Residuals 12 40634 3386
Error: subj:present.type
Df Sum Sq Mean Sq F value Pr(>F)
present.type 1 603.1 603.1 0.7988 0.3890145
group:present.type 1 22775.8 22775.8 30.1661 0.0001379 ***
Residuals 12 9060.1 755.0
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Error: Within
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 840 148493 177
-------------------------------

Now, since the interaction is significant, I want to compute two
simple main effects: to find out if there is a significant difference
between C2 and D2 (present.type var) (i) in group D2C2 and then also
(ii) in group C2D2 (won't be shown to avoid redundancy). To achieve
that:

(1) I run the model separately for each level of group:

dat.g1 = subset(dat.net, group=="D2C2")
m.g1 = aov(t.total ~ present.type + Error(subj/present.type),
data=dat.g1)
summary(m.g1)

Error: subj
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 6 22788 3798
Error: subj:present.type
Df Sum Sq Mean Sq F value Pr(>F)
present.type 1 15395.8 15395.8 18.694 0.004963 **
Residuals 6 4941.4 823.6
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Error: Within
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 420 80658 192
-------------------------------

(2) I use the error term from the overall model (dat.net) to calculate
the MS-Error term:

MS-Effect(from model m.g1) for present.type = 15395.8 with df = 1
MS-Error(from model m.full) for present.type = 755.0 with df = 12
(from Error: subj:present.type)

so we have F(1,12) = 15395.8 / 755.0
which means F = 20.4 and to calculate p-sig:
1 - pf(20.4,1,12)
->> p=0.0007070375

Well, is this the way to do it?
Is it equivalent to or different from using the HH package?

Thanks in advance and best to all,
dror

----------------------

> R-h...@r-project.org mailing listhttps://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guidehttp://www.R-project.org/posting-guide.html

RICHARD M. HEIBERGER

unread,
Feb 27, 2010, 5:10:10 PM2/27/10
to Or Duek, R-h...@r-project.org
> Lets take for example an experiment in which I had two between subject
> variables - Strain and treatment, and one within - exposure. all the
> variables had 2 levels each.
>
> I found an interaction between exposure and Strain and I want to compare
> Strain A and B under every exposure (first and second).
> The general model was with that function:
> aov(duration~(Strain*exposure*treatment)+Error(subject/exposure),data)
>
> in summary(aovmodel) there was a significant interaction between exposure
> and strain.
> how (using those HH packages) can I compare Strains under the conditions of
> exposure?

Your example is structurally similar to the maiz example in ?MMC.
Therefore the answer will also be similar. It is not possible to do
any more without the exact
structure of your data. As I indicated before, the duration variable
can be random numbers.
I need the full dataset with the actual values for Strain, exposure,
treatment, and subject.
You are welcome to use A,B,C,D for treatment levels and 1,2,...,n for
subject ID.

Rich

Reply all
Reply to author
Forward
0 new messages