split-split-plot design and adonis

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","Site7","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-adonis-tp7332029p7342690.html
Sent from the r-sig-ecology mailing list archive at Nabble.com.

	[[alternative HTML version deleted]]