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

[Steve Brewer]

I just realized a mistake in my suggested approach. See below in CAPS.
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/8/12 1:07 PM, "Steve Brewer" <jbrewer at olemiss.edu> wrote:

>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.

UNFORTUNATELY, A SIMPLE SUMMING OF THE RAW ABUNDANCES WON'T WORK HERE
(UNLESS YOU'RE USING EUCLIDEAN DISTANCES, BUT MOST FOLKS ARE USING
BRAY-CURTIS DISTANCES FOR SPECIES COMPOSITION DATA).
>
>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.

THIS IS ONLY CORRECT IF IT IS INDEED A SPLIT-PLOT DESIGN AND NOT
SPLIT-SPLIT PLOT.

>
>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.

THIS WILL WORK. SO, IF YOU'RE ONLY INTERESTED IN THE EFFECTS OF DEPTH AND
ITS INTERACTIONS WITH SEDIMENT AND/OR HYDROLOGY, THEN YOU COULD DO THIS.
>
>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)

NOPE. THIS WON'T WORK. CANNOT SUM THE ABUNDANCES.

EITHER MUST BE DONE THE WAY KAY SUGGESTED OR ONE MUST DO A PCO AND
CALCULATE DISTANCES FROM CENTROIDS.

SORRY.
>
>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","Si
>>t
>>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-ad
>>>o
>>>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
>
>_______________________________________________
>R-sig-ecology mailing list
>R-sig-ecology at r-project.org
>https://stat.ethz.ch/mailman/listinfo/r-sig-ecology