permustats for anova.cca
Eduard Szöcs
July 1, 2015
What I did
See here for a diff to vegandevs/vegan@master.
Critical parts are maybe:
- added
anova.ccaclass to all anova.cca* functions - all return an attribute
F.perm anova.ccalist()returns pfvals as F.perm (not 100% sure about this)
Install branch
Install my branch from github (into `~/R-dev’ to avoid overwriting of current vegan version):
# install.packages(devtools)
require(devtools)
dev_mode()
install_github('EDiLD/vegan@d_perstats')
dev_mode()require(vegan, lib.loc = "~/R-dev")
vegan:::permustats.anova.cca## function (x, ...)
## {
## if (is.null(attr(x, "F.perm")))
## stop("no permutation data available")
## F.perm <- attr(x, "F.perm")
## k <- !is.na(x$F)
## F.0 <- x$F[k]
## structure(list(statistic = structure(F.0, names = rownames(x)[k]),
## permutations = F.perm, alternative = "greater"), class = "permustats")
## }
## <environment: namespace:vegan>Test data
# Test data
data(dune)
data(dune.env)Test function
#' test function
#' @param obj An anova.cca object
chk <- function(obj) {
print(obj)
print(ps <- permustats(obj))
print(summary(ps))
densityplot(ps)
} Test cases
anova.cca
chk(anova(rda(dune ~ Manure + Use, dune.env)))## Permutation test for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure + Use, data = dune.env)
## Df Variance F Pr(>F)
## Model 6 35.363 1.5713 0.016 *
## Residual 13 48.761
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num 1.57
## ..- attr(*, "names")= chr "Model"
## $ permutations: num [1:999] 1 1.066 1.058 0.999 1.067 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## Model 1.5713 2.3221 1.0272 1.0095 1.4323
##
## (Interval (Upper - Lower) = 0.95)
anova.ccabyaxis
chk(anova(rda(dune ~ Manure + Use, dune.env), by = 'axis'))## Permutation test for rda under reduced model
## Marginal tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure + Use, data = dune.env)
## Df Variance F Pr(>F)
## RDA1 1 16.509 4.4015 0.002 **
## RDA2 1 7.955 2.1209 0.049 *
## RDA3 1 5.182 1.3817 0.185
## RDA4 1 3.188 0.8500 0.559
## RDA5 1 1.366 0.3643 0.975
## RDA6 1 1.161 0.3096 0.984
## Residual 13 48.761
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num [1:6] 4.401 2.121 1.382 0.85 0.364 ...
## ..- attr(*, "names")= chr [1:6] "RDA1" "RDA2" "RDA3" "RDA4" ...
## $ permutations: num [1:999, 1:6] 2.087 0.626 0.772 1.229 1.101 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## RDA1 4.4015 5.9860 1.0171 0.8868 2.0895
## RDA2 2.1209 2.0012 1.0174 0.8724 2.1025
## RDA3 1.3817 0.7634 1.0035 0.9057 1.9271
## RDA4 0.8500 -0.3339 1.0312 0.9028 2.1556
## RDA5 0.3643 -1.2788 1.0166 0.9087 1.9594
## RDA6 0.3096 -1.2606 1.0116 0.8862 2.0391
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ Manure, dune.env), by = 'axis'))## Permutation test for rda under reduced model
## Marginal tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure, data = dune.env)
## Df Variance F Pr(>F)
## RDA1 1 16.003 4.3402 0.002 **
## RDA2 1 6.683 1.8125 0.087 .
## RDA3 1 4.642 1.2591 0.253
## RDA4 1 1.489 0.4037 0.972
## Residual 15 55.307
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num [1:4] 4.34 1.813 1.259 0.404
## ..- attr(*, "names")= chr [1:4] "RDA1" "RDA2" "RDA3" "RDA4"
## $ permutations: num [1:999, 1:4] 0.677 0.499 1.063 1.342 1.731 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## RDA1 4.3402 6.2839 1.0207 0.9027 2.0820
## RDA2 1.8125 1.3821 1.0282 0.8903 2.0512
## RDA3 1.2591 0.4535 1.0288 0.9228 2.0186
## RDA4 0.4037 -1.3441 0.9986 0.9010 1.8349
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ 1, dune.env), by = 'axis'))## No constrained component
##
## Model: rda(formula = dune ~ 1, data = dune.env)
## Df Var F N.Perm Pr(>F)
## Model 0 0.000 0 0
## Residual 19 84.124## Error in permustats.anova.cca(obj): no permutation data availableanova.ccabymargin
chk(anova(rda(dune ~ Manure + Use, dune.env), by = 'margin'))## Permutation test for rda under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure + Use, data = dune.env)
## Df Variance F Pr(>F)
## Manure 4 25.153 1.6765 0.020 *
## Use 2 6.546 0.8726 0.607
## Residual 13 48.761
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num [1:2] 1.676 0.873
## ..- attr(*, "names")= chr [1:2] "Manure" "Use"
## $ permutations: num [1:999, 1:2] 0.969 0.696 1.689 1.021 1.136 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## Manure 1.6765 2.4206 1.0192 0.9868 1.5204
## Use 0.8726 -0.4269 1.0199 0.9617 1.6618
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ Manure, dune.env), by = 'margin'))## Permutation test for rda under NA model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure, data = dune.env)
## Df Variance F Pr(>F)
## Manure 4 28.817 1.9539 0.009 **
## Residual 15 55.307
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num 1.95
## ..- attr(*, "names")= chr "Manure"
## $ permutations: num [1:999, 1] 1.043 0.999 0.871 1.15 0.788 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## Manure 1.9539 3.2873 1.0198 0.9785 1.5334
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ 1, dune.env), by = 'margin'))## No constrained component
##
## Model: rda(formula = dune ~ 1, data = dune.env)
## Df Var F N.Perm Pr(>F)
## Model 0 0.000 0 0
## Residual 19 84.124## Error in permustats.anova.cca(obj): no permutation data availableanova.ccabyterm
chk(anova(rda(dune ~ Manure + Use, dune.env), by = 'term'))## Permutation test for rda under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure + Use, data = dune.env)
## Df Variance F Pr(>F)
## Manure 4 28.817 1.9207 0.007 **
## Use 2 6.546 0.8726 0.606
## Residual 13 48.761
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num [1:2] 1.921 0.873
## ..- attr(*, "names")= chr [1:2] "Manure" "Use"
## $ permutations: num [1:999, 1:2] 1.135 0.862 0.746 1.315 0.834 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## Manure 1.9207 3.2141 1.0266 0.9925 1.5233
## Use 0.8726 -0.3777 0.9979 0.9468 1.6245
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ Manure, dune.env), by = 'term'))## Permutation test for rda under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = dune ~ Manure, data = dune.env)
## Df Variance F Pr(>F)
## Manure 4 28.817 1.9539 0.01 **
## Residual 15 55.307
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num 1.95
## ..- attr(*, "names")= chr "Manure"
## $ permutations: num [1:999, 1] 1.347 1.382 1.139 1.11 0.842 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## Manure 1.9539 3.2472 1.0253 0.9835 1.5528
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ 1, dune.env), by = 'term'))## No constrained component
##
## Model: rda(formula = dune ~ 1, data = dune.env)
## Df Var F N.Perm Pr(>F)
## Model 0 0.000 0 0
## Residual 19 84.124## Error in permustats.anova.cca(obj): no permutation data availableanova.ccalist
chk(anova(rda(dune ~ Manure + Use, dune.env),
rda(dune ~ Manure, dune.env)))## Permutation tests for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model 1: dune ~ Manure + Use
## Model 2: dune ~ Manure
## Res.Df Res.Variance Df Variance F Pr(>F)
## 1 13 48.761
## 2 15 55.307 -2 -6.5461 0.8726 0.599
## List of 3
## $ statistic : Named num 0.873
## ..- attr(*, "names")= chr "2"
## $ permutations: num [1:999, 1] 0.768 0.814 1.517 1.207 0.696 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## 2 0.8726 -0.3789 0.9963 0.9522 1.5947
##
## (Interval (Upper - Lower) = 0.95)
chk(anova(rda(dune ~ Manure + Use, dune.env),
rda(dune ~ Manure, dune.env),
rda(dune ~ 1, dune.env)))## Permutation tests for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model 1: dune ~ Manure + Use
## Model 2: dune ~ Manure
## Model 3: dune ~ 1
## Res.Df Res.Variance Df Variance F Pr(>F)
## 1 13 48.761
## 2 15 55.307 -2 -6.5461 0.8726 0.605
## 3 19 84.124 -4 -28.8167 1.9207 0.005 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## List of 3
## $ statistic : Named num [1:2] 0.873 1.921
## ..- attr(*, "names")= chr [1:2] "2" "3"
## $ permutations: num [1:999, 1:2] 0.652 0.844 0.544 0.672 0.786 ...
## $ alternative : chr "greater"
## - attr(*, "class")= chr "permustats"
## NULL
##
## statistic SES mean lower median upper
## 2 0.8726 -0.4241 1.0166 0.9617 1.6496
## 3 1.9207 3.3250 1.0176 0.9726 1.5371
##
## (Interval (Upper - Lower) = 0.95)