[R-sig-eco] split-split-plot design and adonis

[Steve Brewer]

Kay,

Seems like a clever and reasonable approach, but I have a couple of
comments/questions.

First, it seems that with this approach, you cannot evaluate the sediment
x hydrology, the sediment x depth or the sediment x hydrology x depth
interactions. I'm not sure if Arnaud is interested in these interactions,
but one generally is in most split-plots.

Second, I tried your approach with my own data (which is just a split
plot). Perhaps I did it wrong, but the resulting anova seems wrong to me.
The residual df are too high, I think because the analysis is not taking
the strata variation and its df (in Arnaud's case, the site effect) out of
the error term.

Here's a suggestion. First, simplify Arnaud's design to a split plot
(removing the depth effect. In the species file, this will entail summing
abundances across depths for each site x hydrology combination.

Use a design predictor file as originally described by Arnaud (except
without the depth factor); I.e., site1, site2, ... site9, three sediment
types - sed1, sed2, sed3, and two hydrology types - dry, wet.

To analyze the split-plot effects, I would suggest the following

Adonis(community ~ sediment*hydrology + site, strata = site)

You should get an anova table that will give you the sediment main effect,
the hydrology main effect, the sediment*hydrology interaction, the site
effect (probably not of interest), and residual error. ***Note that the
sediment main effect is garbage here because it is tested with the
residual error term, which is not the correct error term and must be
disregarded; Your approach to getting it seems reasonable to me.*** The
hydrology effect and the hydrology*treatment interaction should be
correct, however, because both are tested with the residual (I.e.,
split-plot) error term.

For the split-split plot effects, you need to make a term equivalent to
the site x hydrology term (as described by Arnaud). That will be the
strata effect for the split-split plot (but not for the split-plot
effects). In this case, you'll use a species dataset in which the
abundances are not summed across depths for each site x hydro combination.
This is identical to the original species data file. In this case, I
suggest the following analysis:

Adonis(community ~ sediment*hydrology*depth + site_hydro, strata =
site_hydro).

Here, you disregard all results except those of depth and its interactions
with sediment and/or hydrology.

Another possible way to test the sediment main effect would be to create a
new species data file in which abundances are summed across all
hydrologies and depths for each site. The adonis statement is then simply:

Adonis(community ~ sediment)

I'd appreciate any feedback.

Thanks,
Steve
 


J. Stephen Brewer 
Professor 
Department of Biology
PO Box 1848
 University of Mississippi
University, Mississippi 38677-1848
 Brewer web page - http://home.olemiss.edu/~jbrewer/
FAX - 662-915-5144
Phone - 662-915-1077




On 3/5/12 4:45 AM, "Kay Cichini" <kay.cichini at gmail.com> wrote:

>Hi Arnaud and list-members,
>
>what about the following:
>
>library(vegan)
>
>Sediment<-as.factor(rep(c("Sed1","Sed2","Sed3"),each=36))
>Site<-as.factor(rep(c("Site1","Site2","Site3","Site4","Site5","Site6","Sit
>e7","Site8","Site9"),each=12))
>Hydrology<-as.factor(rep(rep(c("Dry","Wet"),each=6),9))
>Depth<-as.factor(rep(rep(c("D1","D2"),each=3),18))
>
>nobs = length(Site)
>
># depth and hydrology are nested within sites:
>table(Site, Depth)
>table(Site, Hydrology)
>
># so the following call is ok as it restricts permutation on sites -
># that is, levels of depth and hydrology are shuffled only within and not
># across sites:
>set.seed(123)
>community <- matrix(rnorm(nobs * 3, 10, 1), nobs, 3)
>adonis(community ~ Hydrology * Depth, perm = 999, strata = Site)
>
># sediment is different:
>table(Site, Sediment)
>
># to allow for this structure you would need to shuffle sites as
># a whole. the below restriction will randomly assign levels of sediment
># to the sites - like this:
>ctrl = permControl(strata = Site, permute.strata = TRUE)
>table(Site, perm.Sediment = Sediment[permuted.index2(nobs, control =
>ctrl)])
>
>### to test sediment you'd need a customized code:
>(fit <- adonis(community ~ Sediment, permutations = 1))
>
>### number of perms
>B <- 999
>
>### setting up frame which will be populated by
>### random r2 values:
>pop <- rep(NA, B + 1)
>
>### the first entry will be the true r2:
>pop[1] <- fit$aov.tab[1, 5]
>
>for(i in 2:(B+1)){
>     idx <- permuted.index2(nobs, control = ctrl)
>     fit.rand <- adonis(community ~ Sediment[idx], permutations = 1)
>     pop[i] <- fit.rand$aov.tab[1,5]
>}
>
>### get the p-value (H0: Sediment has no effect on location of
>communities):
>print(pval <- sum(pop >= pop[1]) / (B + 1))
>
>### make a histogram to see random R2-values and the true one:
>hist(pop, xlab = "Population R2")
>abline(v = pop[1], col = 2, lty = 3)
>
>What do you think?
>Kay
>
>
>2012/3/4 Arnaud Foulquier <arnaud.foulquier at gmail.com>
>
>> I have run the two analysis. The only differences are the p-values that
>> appear in the anova tables. Other values (Df, SumsOfSqs, MeanSqs,
>>F.Model)
>> are not modified because they are calculated on the "original" (not
>> permuted) data.
>>
>> The p-value is calculated by comparing the value of F obtained with the
>> original data to the distribution obtained by permutations. If i
>>understand
>> well, the strata argument is used to constrain the permutations of the
>> original data within the groups specified in the strata argument. The
>> distribution of the F-statistic that is generated by permutations is not
>> the
>> same according to whether the strata argument is specified or not. This
>>is
>> why the p-values are modified when the strata argument is specified.
>>
>> However, I'm not sure that strata=Site_Hydro correctly specify the
>>nested
>> structure of my data.
>>
>>
>>
>>
>>
>> --
>> View this message in context:
>> 
>>http://r-sig-ecology.471788.n2.nabble.com/split-split-plot-design-and-ado
>>nis-tp7332029p7342690.html
>> Sent from the r-sig-ecology mailing list archive at Nabble.com.
>>
>> _______________________________________________
>> R-sig-ecology mailing list
>> R-sig-ecology at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
>>
>
>	[[alternative HTML version deleted]]
>
>_______________________________________________
>R-sig-ecology mailing list
>R-sig-ecology at r-project.org
>https://stat.ethz.ch/mailman/listinfo/r-sig-ecology