0

I try to create more then three independent contrasts in PERMANOVA with 4 factors without success. I need to use all possible pairwise combinations of factor levels in my contr2df object. There are any way for make this possible? In my code:

#1st factor
treat <- gl(4, 15, labels = paste("t", 1:4, sep="")); treat 

#Variables
set.seed(124)
sp  <- cbind(c(rnorm(10,  5, 0.25), rnorm(50, 2.5, 0.25)), rnorm(60, 2.5, 0.25),
             c(rnorm(10, 12, 0.25), rnorm(50, 2.5, 0.25)), rnorm(60, 2.5, 0.25))
colnames(sp) <- c("sp1", "sp2", "sp3", "sp4")
head(sp))


#create a design matrix of the contrasts for "treat" 
Treat_Imp<-model.matrix(~treat-1)

require(vegan)

fullModel <- adonis(sp ~ treat, method = "euclidean", permutations = 9999)

fullModel

#Comparisons
TI    <- model.matrix(~ treat-1)
head(TI)

f <- nlevels(treat)
comb <- t(combn(1:f, 2))
n    <- nrow(comb)

contr2 <- NULL
for (x in 1:n) {
      i <- comb[x, 1]
      j <- comb[x, 2]
      tmp <- list(TI[,i] - TI[,j]); names(tmp) <- paste0("TI",i, "_", j)
       contr2 <- c(contr2, tmp) }
contr2


contr2df <- as.data.frame(contr2)

adonis(
     sp ~ ., data = contr2df,
     method = "euclidean",
     permutations = 9999)
# 

Thanks,

Alexandre

Leprechault
  • 1,531
  • 12
  • 28
  • _"I need to use all possible pairwise combinations of factor levels in my contr2df object"_ - e.g. `AlgDesign::gen.factorial(3, 6, TRUE, varNames=c("TI1_2", "TI1_3", "TI1_4", "TI2_3", "TI2_4", "TI3_4"))` creates a full factorial design. – lukeA Mar 11 '16 at 12:36
  • But doen't work in adonis() function, I have the error: Error in G * t(hat) – Leprechault Mar 11 '16 at 13:08
  • ... because `sp` has only 60 rows, and `contrdf` has now 729. – lukeA Mar 11 '16 at 13:21
  • And what the solution for this, because a need to inclued sp variables? – Leprechault Mar 11 '16 at 13:25
  • Well, the solution is to make `sp` have 729 rows :-) : e.g. `sp <- cbind(c(rnorm(29, 5, 0.25), rnorm(700, 2.5, 0.25)), rnorm(729, 2.5, 0.25), c(rnorm(29, 12, 0.25), rnorm(700, 2.5, 0.25)), rnorm(729, 2.5, 0.25))`. I dunno what you are analysing. Perhaps the question is better suited for Stackexchange than for Stackoverflow? – lukeA Mar 11 '16 at 13:36
  • Thank you very much, lukeA!!! – Leprechault Mar 11 '16 at 13:45
  • You're welcome. (BTW: I meant Cross Validated, not Stackexchange.) – lukeA Mar 11 '16 at 13:52

0 Answers0