CODEWITHSUNDEEP
ranova
user3170629
user3170629

Reputation: 499

Two-way ANOVA with means and standard deviation

Is there a way to use a two-way ANOVA in R that also provides the means and standard deviation for the tested variables?

Currently I am using the code

summary(aov(Tortuosity~Site*Species,data=BF.C))

But this only provides the F value, probability, residuals and so on.

Thanks

Upvotes: 1

Views: 5420

Answers (1)

lord.garbage
lord.garbage

Reputation: 5960

Yes, you can use model.tables(). Try ?model.tables to get acquainted with it. What you want is probably (I’ll use an example):

aov_out <- aov(yield ~ N*P, npk)

model.tables(aov_out, "means") # This will give you the means.
model.tables(aov_out, se = TRUE) # Will give you standard errors for the effects.
model.tables(aov_out, "means", se = TRUE) # Will give you the standard error for the differences of means.

Note however, that the last command only works as long as you don’t have any random effects in your model. Hence, for a model like:

aov_out_ran <- aov(yield ~ N*P + Error(block), npk)

the last command will not work as it is not yet implemented. But you will see that in a warning message anyway.

You can also just compute the means easily by hand. Refit your model with contrast coding:

aov_out_contr <- aov(yield ~ N*P, npk, contrasts = list (N = "contr.sum", P = "contr.sum"))

Use coef() to get the coefficients of the model:

coef_aov_out <- coef(aov_out)

Because we used contrast coding the (Intercept) at position coef_aov_out[1] will be the grand mean and the further coefficients in coef_aov_out will be the effects of the main effects or interaction effects that need to be substracted or added to the grand mean in order to get the group specific means. You can compute them by hand like this:

# Compute mean of N0:
N0 <- coef_aov_out[1]+coef_aov_out[2]

# Compute mean of N1:
N1 <- coef_aov_out[1]-coef_aov_out[2]

# Compute mean of P0:
P0 <- coef_aov_out[1]+coef_aov_out[3]

# Compute mean of P1:
P1 <- coef_aov_out[1]-coef_aov_out[3]

# Compute mean of N0xP0:
NOP0 <- coef_aov_out[1]+coef_aov_out[2]+coef_aov_out[3]+coef_aov_out[4]

# Compute mean of N0xP1:
N0P1 <- coef_aov_out[1]+coef_aov_out[2]-coef_aov_out[3]-coef_aov_out[4]

# Compute mean of N1xP0:
N1P0 <- coef_aov_out[1]-coef_aov_out[2]+coef_aov_out[3]-coef_aov_out[4]

# Compute mean of N1xP1:
N1P1 <- coef_aov_out[1]-coef_aov_out[2]-coef_aov_out[3]+coef_aov_out[4]

You can compare the results with model.tables(aov_out_contr, "means"). It’s a good exercise to understand what’s going on.

Upvotes: 3

Related Questions