[R-sig-eco] Permanova with nested data
Steve Brewer
jbrewer at olemiss.edu
Mon Feb 25 16:13:18 CET 2013
Beth and others,
Given several recent queries regarding how to analyze repeated-measures
and split-plot perManova using adonis, I thought I would pass along what I
think is a reasonable solution.
I just saw the recent exchange over the use of BiodiversityR to do nested
perMANOVA. I was unaware of this function until today.
With that in mind, it is possible to do a simple repeated-measures
permanova using two different analyses, one for the within-subjects
effects and one for the between-subjects effects. The same approach
applies to a simple split-plot analysis.
For the within-subjects (or sub-plot) effects, you use adonis and the
strata function. The model formula could look something like:
Assume species responses are in "Speciesdata" and the treatment, time, and
plot effects are in "envfactors"
Adonis(Speciesdata ~ betweensubtrtmt * time + plot, data = envfactors,
strata = plot)
Where plot is nested within the betweensubtrtmnt
Strata restricts the permutation, and the residual error term will give
you the correct test for the time effect and the betweensubtrtmnt * time
interaction, but the test for the betweensubtrmnt main effect will be
wrong because plot, and not the residual error term, is the correct error
term for testing it.
To get a test for the betweensubtrtmnt main effect, load the BiodiversityR
package (I use the 1.6 version, but see the recent discussion about this)
and use the nested.npmanova function.
Nested.npmanova(speciesdata ~ betweensubtrtmnt + plot; data = envfactors)
In this case, the betweensubtrtmnt is tested with plot; plot is tested
with the residual error term but that latter test is not correct in this
instance and is usually not of interest anyway.
Note that the default distance is euclidean; you'll to use "method" to
specify a different distance, e.g.,
Nested.npmanova(speciesdata ~ betweensubtrtmnt + plot; data = envfactors;
method ="bray")
The same principles apply to a simple split-plot design, except that the
whole-plot treatment is treated like a between-subjects treatment and the
sub-plot treatment is treated like the time effect.
Hope this is of some help.
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 12/5/12 10:04 AM, "Beth Atkinson" <beth.atkinson at bristol.ac.uk> wrote:
>Hi,
>
>Apologies if this has been asked before, I looked looked through the
>archives and couldn't find a solution.
>
>
>I have plant community data from 32 woodland plots. Plots are grouped
>into
>sites (8 sites in total). Half the sites (4) are on one soil type, and
>half
>on another. At each site there are 4 plots, each under a different
>management regime . I want to know if the plant communities under
>different management regimes and on different soil types differ.
>
>My data is as follows:
>
>
>veg <- a data frame containg the community data (abundances of each
>species
>at each plot)
>
>
>environ<-read.csv(file.choose(),header=TRUE)
>
>> str(environ)
>'data.frame': 32 obs. of 4 variables:
> $ code: Factor w/ 32 levels "cf.1","cf.3",..: 1 9 17 25 2 10 18 26 3 11
>...
> $ site: Factor w/ 8 levels "eight","five",..: 5 5 5 5 8 8 8 8 3 3 ...
> $ type: Factor w/ 4 levels "clearfell","plantation",..: 1 2 3 4 1 2 3 4
>1
>2 ...
> $ nvc : Factor w/ 2 levels "W10","W8": 2 2 2 2 2 2 2 2 1 1 ...
>
># site is the name of each site
># type refers to the management regime
># nvc refers to the soil type (actually the NVC classification)
>
>Prior to PERMANOVA I used betadisper() to test for homogeneity of
>multivariate dispersion. No difference in dispersion was found either
>between plots on different site types, or plots on different NVC types.
>
>
>I carried out a PERMANOVA using adonis on this data as follows:
>
>adon.mod1.bray<-adonis(veg~ type*nvc, data=environ,strata=environ$site,
> method = "bray", permutations=999)
>
>I used strata=environ$site as management regime is nested within site.
>
>> adon.mod1.bray
>
>Call:
>adonis(formula = veg ~ type * nvc, data = environ, permutations = 999,
> method = "bray", strata = environ$site)
>
>Terms added sequentially (first to last)
>
> Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
>type 3 2.2358 0.74528 3.5742 0.24662 0.001 ***
>nvc 1 1.0905 1.09045 5.2296 0.12028 0.001 ***
>type:nvc 3 0.7353 0.24509 1.1754 0.08110 0.208
>Residuals 24 5.0044 0.20852 0.55200
>Total 31 9.0659 1.00000
>
>
>Is this acceptable. Is the strata=environ$site part correct?
>
>If data points are only being permutted within sites, and each site only
>occurs on one soil type, is there no permutation between NVC/soil types?
>
>I've also heard that adonis does not give the correct p values when data
>is
>nested. Is this correct, and a problem in the above example? If so is
>there
>a more suitable analysis? Would a db-RDA as follows be ok?:
>
>
>dbRDA<-capscale(veg~type*nvc+Condition(site), data=environ,
>distance="bray", add=TRUE)
>
>anova.cca(dbRDA, stop=999)
>
>anova.cca(dbRDA, by="terms", step=999)
>
>> anova.cca(dbRDA, by="terms", step=999)
>Permutation test for capscale under reduced model
>Terms added sequentially (first to last)
>
>Model: capscale(formula = veg ~ type + nvc + Condition(site) + type:nvc,
>data = environ, distance = "bray", add = TRUE)
> Df Var F N.Perm Pr(>F)
>type 3 2.9615 3.3215 999 <2e-16 ***
>type:nvc 3 1.1299 1.2673 999 0.1141
>Residual 18 5.3496
>
>
>However, this doesn't show a result for nvc. Why?
>
>I've also been looking into the manyglm() function in mvabund, but I
>don't
>think this can accept random effects.
>
>Any help would be greatly appreciated,
>
>Thanks,
>
>Beth
>
>--------------------------------------------------------
>Beth Atkinson
>
>PhD student
>Community Ecology Group
>School of Biological Sciences
>University of Bristol
>Woodland Road
>Bristol
>BS8 1UG
>
>_______________________________________________
>R-sig-ecology mailing list
>R-sig-ecology at r-project.org
>https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
More information about the R-sig-ecology
mailing list