[R-sig-eco] Issue with BiodiversityR::nested.npmanova
Baldwin, Jim -FS
jbaldwin at fs.fed.us
Sat Feb 23 16:29:25 CET 2013
I suspect it must be a version problem. Running Windows XP (yes, still XP) with R version 2.13.1 (2011-07-08) and just now loading vegan and BiodiversityR I get the following:
> library(BiodiversityR)
Loading required package: vegan
Loading required package: permute
This is vegan 2.0-2
Warning messages:
1: package 'BiodiversityR' was built under R version 2.13.2
2: package 'vegan' was built under R version 2.13.2
3: package 'permute' was built under R version 2.13.2
> library(vegan)
> data(mite)
>
> mite <- mite[1:60,]
> env <- data.frame(Fac = rep(LETTERS[1:2], each = 30), NFac = rep(letters[1:4], each = 15))
>
> nested.npmanova(mite ~ Fac + NFac, data = env, method = "jac", permutations = 1000)
Total sum of squares of distance matrix: 16.77718
Total sum of squares for non-parametric manova: 16.7771827228089
Nested anova for NFac nested within Fac
Df SumsofSquares F N.Perm Pr(>F)
Fac 1 3.6153 7.7348 1000 0.347652
NFac 2 0.9348 2.1408 1000 0.001998 **
Residual 56 12.2270
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Jim
-----Original Message-----
From: r-sig-ecology-bounces at r-project.org [mailto:r-sig-ecology-bounces at r-project.org] On Behalf Of Kay Cichini
Sent: Saturday, February 23, 2013 5:56 AM
To: Jari Oksanen; r-sig-ecology at r-project.org
Subject: Re: [R-sig-eco] Issue with BiodiversityR::nested.npmanova
I also tried the following:
library(BiodiversityR)
library(vegan)
data(mite)
mite <- mite[1:60,]
env <- data.frame(Fac = rep(LETTERS[1:2], each = 30), NFac = rep(letters[1:4], each = 15))
nested.npmanova(mite ~ Fac + NFac, data = env, method = "jac", permutations = 1000)
throwing this error:
Total sum of squares of distance matrix: 16.77718 Fehler in model.frame(data = data1, na.action = function (object, ...) :
Objekt 'data1' nicht gefunden
I inspected the source of the function but can not figure out why data1 is not assigned correctly:
function (formula, data, method = "euc", permutations = 100,
warnings = FALSE)
{
randomize = function(data, toplev, lowlev) {
newdata <- data
orig.levs <- levels(data[, lowlev])
nl <- length(orig.levs)
new.levs <- orig.levs[sample(nl)]
for (i in 1:nl) {
subs1 <- data[, lowlev] == orig.levs[i]
subs2 <- data[, lowlev] == new.levs[i]
newtops <- data[subs2, toplev]
newtops <- newtops[1]
newtops <- rep(newtops, sum(subs1))
newdata[subs1, toplev] <- newtops
}
return(newdata)
}
randomize2 = function(data, strata) {
newdata <- data
orig.levs <- levels(data[, strata])
nl <- length(orig.levs)
for (i in 1:nl) {
subs <- data[, strata] == orig.levs[i]
nsub <- sum(subs == T)
subdata <- data[subs, ]
newdata[subs, ] <- subdata[sample(nsub), ]
}
return(newdata)
}
ow <- options("warn")
if (warnings == FALSE) {
options(warn = 0)
}
formula <- as.formula(formula)
.BiodiversityR <- new.env()
environment(formula) <- .BiodiversityR
if (length(all.vars(formula)) > 3)
stop(paste("function only works with one main and one nested
factor"))
x <- eval(as.name((all.vars(formula)[1])))
if (inherits(x, "dist")) {
distmatrix <- as.matrix(x)
}
else {
distmatrix <- as.matrix(vegdist(x, method = method))
}
SStot <- sum(distmatrix^2)/(2 * nrow(distmatrix))
cat("Total sum of squares of distance matrix:", SStot, "\n")
resp <- all.vars(formula)[1]
toplev <- all.vars(formula)[2]
lowlev <- all.vars(formula)[3]
environment(formula) <- .GlobalEnv
data1 <- data
assign("data1", data, envir = .BiodiversityR)
model <- capscale(formula, data1, distance = "euclidean")
anovares <- anova(model, perm.max = 1, by = "terms")
anovadat <- data.frame(anovares)
df1 <- anovadat[1, 1]
df2 <- anovadat[2, 1]
df3 <- nrow(distmatrix) - df1 - df2 - 1
adonis1 <- adonis(formula, data1, permutations = 2, method = method)
adonis1 <- data.frame(adonis1$aov.tab)
sstop <- anovadat[1, 2] <- adonis1[1, 2]
sslow <- anovadat[2, 2] <- adonis1[2, 2]
ssres <- anovadat[3, 2] <- adonis1[3, 2]
vartot <- adonis1[4, 2]
Ftop <- anovadat[1, 3] <- (sstop/df1)/(sslow/df2)
Flow <- anovadat[2, 3] <- (sslow/df2)/(ssres/df3)
counter <- 1
for (i in 1:permutations) {
data2 <- randomize(data, toplev, lowlev)
assign("data2", data2, envir = .BiodiversityR)
adonis2 <- adonis(formula, data = data2, method = method,
permutations = 2)
adonis2 <- data.frame(adonis2$aov.tab)
Frand <- (adonis2[1, 2]/df1)/(adonis2[2, 2]/df2)
if (Frand >= Ftop) {
counter <- counter + 1
}
}
signi <- counter/(permutations + 1)
anovadat[1, 4] <- anovadat[2, 4] <- permutations
anovadat[1, 5] <- signi
counter <- 1
for (i in 1:permutations) {
data2 <- randomize2(data, toplev)
assign("data2", data2, envir = .BiodiversityR)
adonis2 <- adonis(formula, data = data2, method = method,
permutations = 2)
adonis2 <- data.frame(adonis2$aov.tab)
Frand <- (adonis2[2, 2]/df2)/(adonis2[3, 2]/df3)
if (Frand >= Flow) {
counter <- counter + 1
}
}
remove("data1", envir = .BiodiversityR)
remove("data2", envir = .BiodiversityR)
signi <- counter/(permutations + 1)
anovadat[2, 5] <- signi
colnames(anovadat) <- c("Df", "SumsofSquares", "F", "N.Perm",
"Pr(>F)")
mod <- paste("Nested anova for", lowlev, "nested within",
toplev, "\n")
head <- paste("Total sum of squares for non-parametric manova:",
vartot, "\n")
options(ow)
structure(anovadat, heading = c(head, mod), Random.seed = NA,
class = c("anova.cca", "anova", "data.frame")) }
<environment: namespace:BiodiversityR>
2013/2/23 Jari Oksanen <jari.oksanen at oulu.fi>
> Kay,
>
> This should not work if the function is correctly written. You say
> that the terms in your formula are in data=env, but there are no
> variables env$NFac and env$Fac in env -- but there are NFac and Fac.
>
> Cheers, Jari
>
> Sent from my iPad
>
> On 23.2.2013, at 12.15, "Kay Cichini" <kay.cichini at gmail.com> wrote:
>
> > Hi all,
> >
> > Does anyone have a clue why this is not working:
> >
> > library(BiodiversityR)
> > library(vegan)
> >
> > data(mite)
> > mite <- mite[1:60,]
> >
> > env <- data.frame(Fac = rep(LETTERS[1:2], each = 30), NFac =
> > rep(letters[1:4], each = 15))
> >
> > nested.npmanova(mite ~ env$Fac + env$NFac, data = env, method =
> > "jac", permutations = 1000)
> >
> > Regards,
> > Kay
> >
> > --
> >
> > Kay Cichini, MSc Biol
> >
> > Grubenweg 22, 6071 Aldrans
> >
> > Tel.: 0650 9359101
> >
> > E-Mail: kay.cichini at gmail.com
> >
> > Web: www.theBioBucket.blogspot.co.at<
> http://www.thebiobucket.blogspot.co.at/><
> http://www.theBioBucket.blogspot.co.at>
> > --
> >
> > [[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
>
--
Kay Cichini, MSc Biol
Grubenweg 22, 6071 Aldrans
Tel.: 0650 9359101
E-Mail: kay.cichini at gmail.com
Web: www.theBioBucket.blogspot.co.at<http://www.thebiobucket.blogspot.co.at/><http://www.theBioBucket.blogspot.co.at>
--
[[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
This electronic message contains information generated by the USDA solely for the intended recipients. Any unauthorized interception of this message or the use or disclosure of the information it contains may violate the law and subject the violator to civil or criminal penalties. If you believe you have received this message in error, please notify the sender and delete the email immediately.
More information about the R-sig-ecology
mailing list