[Vegan-commits] r2488 - pkg/vegan/tests/Examples
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Apr 9 10:43:27 CEST 2013
Author: jarioksa
Date: 2013-04-09 10:43:27 +0200 (Tue, 09 Apr 2013)
New Revision: 2488
Modified:
pkg/vegan/tests/Examples/vegan-Ex.Rout.save
Log:
correct (= up-to-date) Examples in tests/
Modified: pkg/vegan/tests/Examples/vegan-Ex.Rout.save
===================================================================
--- pkg/vegan/tests/Examples/vegan-Ex.Rout.save 2013-04-08 17:13:29 UTC (rev 2487)
+++ pkg/vegan/tests/Examples/vegan-Ex.Rout.save 2013-04-09 08:43:27 UTC (rev 2488)
@@ -1,6 +1,6 @@
-R version 2.13.1 (2011-07-08)
-Copyright (C) 2011 The R Foundation for Statistical Computing
+R version 2.15.3 (2013-03-01) -- "Security Blanket"
+Copyright (C) 2013 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
@@ -23,7 +23,7 @@
> options(warn = 1)
> library('vegan')
Loading required package: permute
-This is vegan 2.1-1
+This is vegan 2.1-28
>
> assign(".oldSearch", search(), pos = 'CheckExEnv')
> cleanEx()
@@ -154,17 +154,17 @@
> plot(ef)
> ordisurf(mod ~ pH, varechem, knots = 1, add = TRUE)
Loading required package: mgcv
-This is mgcv 1.7-6. For overview type 'help("mgcv-package")'.
+This is mgcv 1.7-22. For overview type 'help("mgcv-package")'.
Family: gaussian
Link function: identity
Formula:
y ~ poly(x1, 1) + poly(x2, 1)
-<environment: 0x102618388>
+<environment: 0x1024f0310>
Total model degrees of freedom 3
-GCV score: 0.0427924
+GCV score: 0.04278782
>
>
>
@@ -425,77 +425,19 @@
Inertia is mean squared contingency coefficient
Eigenvalues for constrained axes:
- CCA1 CCA2 CCA3
-0.41868 0.13304 0.07659
+ CCA1 CCA2 CCA3
+0.4187 0.1330 0.0766
Eigenvalues for unconstrained axes:
- CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
-0.409782 0.225913 0.176062 0.123389 0.108171 0.090751 0.085878 0.060894
- CA9 CA10 CA11 CA12 CA13 CA14 CA15 CA16
-0.056606 0.046688 0.041926 0.020103 0.014335 0.009917 0.008505 0.008033
+ CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8 CA9 CA10 CA11
+0.4098 0.2259 0.1761 0.1234 0.1082 0.0908 0.0859 0.0609 0.0566 0.0467 0.0419
+ CA12 CA13 CA14 CA15 CA16
+0.0201 0.0143 0.0099 0.0085 0.0080
-> ## The same, but based on permutation P-values
-> ordistep(cca(dune ~ 1, dune.env), reformulate(names(dune.env)), perm.max=200)
-
-Start: dune ~ 1
-
- Df AIC F N.Perm Pr(>F)
-+ Moisture 3 86.608 2.2536 199 0.005 **
-+ Management 3 86.935 2.1307 199 0.005 **
-+ Manure 4 88.832 1.5251 199 0.025 *
-+ A1 1 87.411 2.1400 199 0.035 *
-+ Use 2 89.134 1.1431 99 0.130
----
-Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-
-Step: dune ~ Moisture
-
- Df AIC F N.Perm Pr(>F)
-- Moisture 3 87.657 2.2536 99 0.01 **
----
-Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-
- Df AIC F N.Perm Pr(>F)
-+ Management 3 86.813 1.4565 199 0.035 *
-+ Use 2 87.259 1.2760 199 0.095 .
-+ Manure 4 87.342 1.3143 199 0.095 .
-+ A1 1 86.992 1.2624 99 0.170
----
-Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-
-Step: dune ~ Moisture + Management
-
- Df AIC F N.Perm Pr(>F)
-- Management 3 86.608 1.4565 199 0.035 *
-- Moisture 3 86.935 1.5518 99 0.020 *
----
-Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-
- Df AIC F N.Perm Pr(>F)
-+ A1 1 86.190 1.6817 199 0.09 .
-+ Manure 3 88.430 0.8167 99 0.58
-+ Use 2 88.245 0.7534 99 0.65
----
-Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-
-Call: cca(formula = dune ~ Moisture + Management, data = dune.env)
-
- Inertia Proportion Rank
-Total 2.1153 1.0000
-Constrained 1.0024 0.4739 6
-Unconstrained 1.1129 0.5261 13
-Inertia is mean squared contingency coefficient
-
-Eigenvalues for constrained axes:
- CCA1 CCA2 CCA3 CCA4 CCA5 CCA6
-0.44583 0.28869 0.11239 0.07166 0.04937 0.03444
-
-Eigenvalues for unconstrained axes:
- CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
-0.350396 0.152057 0.125084 0.109838 0.092209 0.077107 0.059441 0.047755
- CA9 CA10 CA11 CA12 CA13
-0.036958 0.022266 0.020700 0.010827 0.008252
-
+> ## see ?ordistep to do the same, but based on permutation P-values
+> ## Not run:
+> ##D ordistep(cca(dune ~ 1, dune.env), reformulate(names(dune.env)), perm.max=200)
+> ## End(Not run)
> ## Manual model building
> ## -- define the maximal model for scope
> mbig <- rda(dune ~ ., dune.env)
@@ -505,21 +447,21 @@
> add1(m0, scope=formula(mbig), test="perm")
Df AIC F N.Perm Pr(>F)
<none> 89.620
-A1 1 89.591 1.9217 199 0.055 .
+A1 1 89.591 1.9217 199 0.070 .
Moisture 3 87.707 2.5883 199 0.005 **
Management 3 87.082 2.8400 199 0.005 **
-Use 2 91.032 1.1741 99 0.270
+Use 2 91.032 1.1741 99 0.180
Manure 4 89.232 1.9539 199 0.010 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> m0 <- update(m0, . ~ . + Management)
> add1(m0, scope=formula(mbig), test="perm")
- Df AIC F N.Perm Pr(>F)
-<none> 87.082
-A1 1 87.424 1.2965 99 0.21
-Moisture 3 85.567 1.9764 199 0.03 *
-Use 2 88.284 1.0510 99 0.41
-Manure 3 87.517 1.3902 199 0.07 .
+ Df AIC F N.Perm Pr(>F)
+<none> 87.082
+A1 1 87.424 1.2965 99 0.240
+Moisture 3 85.567 1.9764 199 0.005 **
+Use 2 88.284 1.0510 99 0.430
+Manure 3 87.517 1.3902 199 0.130
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> m0 <- update(m0, . ~ . + Moisture)
@@ -527,16 +469,16 @@
> drop1(m0, test="perm")
Df AIC F N.Perm Pr(>F)
<none> 85.567
-Management 3 87.707 2.1769 199 0.015 *
-Moisture 3 87.082 1.9764 199 0.005 **
+Management 3 87.707 2.1769 199 0.010 **
+Moisture 3 87.082 1.9764 199 0.015 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> add1(m0, scope=formula(mbig), test="perm")
Df AIC F N.Perm Pr(>F)
<none> 85.567
-A1 1 86.220 0.8359 99 0.66
-Use 2 86.842 0.8027 99 0.66
-Manure 3 85.762 1.1225 99 0.31
+A1 1 86.220 0.8359 99 0.72
+Use 2 86.842 0.8027 99 0.77
+Manure 3 85.762 1.1225 99 0.26
>
>
>
@@ -549,7 +491,8 @@
> ### Name: adipart
> ### Title: Additive Diversity Partitioning and Hierarchical Null Model
> ### Testing
-> ### Aliases: adipart hiersimu print.hiersimu
+> ### Aliases: adipart adipart.default adipart.formula hiersimu
+> ### hiersimu.default hiersimu.formula
> ### Keywords: multivariate
>
> ### ** Examples
@@ -565,10 +508,10 @@
+ out <- rep(1, length(x))
+ for (i in 2:(length(cut) - 1))
+ out[which(x > cut[i] & x <= cut[(i + 1)])] <- i
-+ return(as.factor(out))}
++ return(out)}
> ## The hierarchy of sample aggregation
> levsm <- data.frame(
-+ l1=as.factor(1:nrow(mite)),
++ l1=1:nrow(mite),
+ l2=cutter(mite.xy$y, cut = seq(0, 10, by = 2.5)),
+ l3=cutter(mite.xy$y, cut = seq(0, 10, by = 5)),
+ l4=cutter(mite.xy$y, cut = seq(0, 10, by = 10)))
@@ -579,41 +522,90 @@
> plot(mite.xy, main="l3", col=as.numeric(levsm$l3)+1)
> par(mfrow=c(1,1))
> ## Additive diversity partitioning
-> adipart(mite ~., levsm, index="richness", nsimul=19)
-adipart with 19 simulations
-with index richness, weights unif
+> adipart(mite, index="richness", nsimul=19)
+adipart object
- statistic z 2.5% 50% 97.5% Pr(sim.)
-alpha.1 15.1143 -38.7550 22.0321 22.3000 22.608 0.05 *
-alpha.2 29.7500 -27.1142 34.5000 34.7500 35.000 0.05 *
-alpha.3 33.0000 0.0000 35.0000 35.0000 35.000 0.05 *
-gamma 35.0000 0.0000 35.0000 35.0000 35.000 1.00
-beta.1 14.6357 9.0433 12.1629 12.4500 12.955 0.05 *
-beta.2 3.2500 16.4371 0.0000 0.2500 0.500 0.05 *
-beta.3 2.0000 0.0000 0.0000 0.0000 0.000 0.05 *
+Call: adipart(y = mite, index = "richness", nsimul = 19)
+
+nullmodel method ‘r2dtable’ with 19 simulations
+options: index richness, weights unif
+alternative hypothesis: simulated median is not equal to the statistic
+
+ statistic z mean 2.5% 50% 97.5% Pr(sim.)
+alpha.1 15.114 -38.43 22.344 22.032 22.300 22.608 0.05 *
+gamma 35.000 0.00 35.000 35.000 35.000 35.000 1.00
+beta.1 19.886 38.43 12.656 12.392 12.700 12.968 0.05 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> adipart(mite ~ ., levsm, index="richness", nsimul=19)
+adipart object
+
+Call: adipart(formula = mite ~ ., data = levsm, index = "richness",
+nsimul = 19)
+
+nullmodel method ‘r2dtable’ with 19 simulations
+options: index richness, weights unif
+alternative hypothesis: simulated median is not equal to the statistic
+
+ statistic z mean 2.5% 50% 97.5% Pr(sim.)
+alpha.1 15.114 -46.2370 22.39624 22.12571 22.44286 22.6236 0.05 *
+alpha.2 29.750 -21.7076 34.81579 34.36250 35.00000 35.0000 0.05 *
+alpha.3 33.000 0.0000 35.00000 35.00000 35.00000 35.0000 0.05 *
+gamma 35.000 0.0000 35.00000 35.00000 35.00000 35.0000 1.00
+beta.1 14.636 9.0407 12.41955 12.00750 12.42857 12.8743 0.05 *
+beta.2 3.250 13.1373 0.18421 0.00000 0.00000 0.6375 0.05 *
+beta.3 2.000 0.0000 0.00000 0.00000 0.00000 0.0000 0.05 *
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> ## Hierarchical null model testing
> ## diversity analysis (similar to adipart)
-> hiersimu(mite ~., levsm, diversity, relative=TRUE, nsimul=19)
-hiersimu with 19 simulations
+> hiersimu(mite, FUN=diversity, relative=TRUE, nsimul=19)
+hiersimu object
- statistic z 2.5% 50% 97.5% Pr(sim.)
-l1 0.76064 -65.47286 0.93511 0.93959 0.9437 0.05 *
-l2 0.89736 -127.77766 0.99635 0.99815 0.9989 0.05 *
-l3 0.92791 -516.33891 0.99921 0.99948 0.9997 0.05 *
-l4 1.00000 0.00000 1.00000 1.00000 1.0000 1.00
+Call: hiersimu(y = mite, FUN = diversity, relative = TRUE, nsimul = 19)
+
+nullmodel method ‘r2dtable’ with 19 simulations
+
+alternative hypothesis: simulated median is not equal to the statistic
+
+ statistic z mean 2.5% 50% 97.5% Pr(sim.)
+level_1 0.76064 -71.195 0.93904 0.93487 0.93856 0.9444 0.05 *
+leve_2 1.00000 0.000 1.00000 1.00000 1.00000 1.0000 1.00
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> hiersimu(mite ~., levsm, FUN=diversity, relative=TRUE, nsimul=19)
+hiersimu object
+
+Call: hiersimu(formula = mite ~ ., data = levsm, FUN = diversity,
+relative = TRUE, nsimul = 19)
+
+nullmodel method ‘r2dtable’ with 19 simulations
+
+alternative hypothesis: simulated median is not equal to the statistic
+
+ statistic z mean 2.5% 50% 97.5% Pr(sim.)
+l1 0.76064 -75.139 0.93833 0.93389 0.93819 0.9427 0.05 *
+l2 0.89736 -110.968 0.99811 0.99699 0.99814 0.9999 0.05 *
+l3 0.92791 -417.338 0.99940 0.99904 0.99943 0.9996 0.05 *
+l4 1.00000 0.000 1.00000 1.00000 1.00000 1.0000 1.00
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> ## Hierarchical testing with the Morisita index
> morfun <- function(x) dispindmorisita(x)$imst
> hiersimu(mite ~., levsm, morfun, drop.highest=TRUE, nsimul=19)
-hiersimu with 19 simulations
+hiersimu object
- statistic z 2.5% 50% 97.5% Pr(sim.)
-l1 0.52070 4.98527 0.31016 0.36570 0.4227 0.05 *
-l2 0.60234 12.33099 0.11979 0.17096 0.2283 0.05 *
-l3 0.67509 19.37352 -0.24164 -0.16761 -0.0895 0.05 *
+Call: hiersimu(formula = mite ~ ., data = levsm, FUN = morfun,
+drop.highest = TRUE, nsimul = 19)
+
+nullmodel method ‘r2dtable’ with 19 simulations
+
+alternative hypothesis: simulated median is not equal to the statistic
+
+ statistic z mean 2.5% 50% 97.5% Pr(sim.)
+l1 0.52070 8.5216 0.353253 0.322624 0.351073 0.3848 0.05 *
+l2 0.60234 14.3854 0.153047 0.096700 0.150434 0.1969 0.05 *
+l3 0.67509 20.3162 -0.182473 -0.234793 -0.195937 -0.0988 0.05 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
@@ -641,6 +633,8 @@
Call:
adonis(formula = dune ~ Management * A1, data = dune.env, permutations = 99)
+Terms added sequentially (first to last)
+
Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
Management 3 1.4686 0.48953 3.2629 0.34161 0.01 **
A1 1 0.4409 0.44089 2.9387 0.10256 0.02 *
@@ -677,12 +671,12 @@
>
> Y <- data.frame(Agropyron, Schizachyrium)
> mod <- metaMDS(Y)
-Run 0 stress 0.08556588
-Run 1 stress 0.1560545
-Run 2 stress 0.08556612
-... procrustes: rmse 0.0001630123 max resid 0.0003642025
+Run 0 stress 0.08556586
+Run 1 stress 0.1560544
+Run 2 stress 0.08556586
+... New best solution
+... procrustes: rmse 1.094365e-06 max resid 1.88838e-06
*** Solution reached
-
> plot(mod)
> ### Hulls show treatment
> ordihull(mod, group=dat$NO3, show="0")
@@ -691,28 +685,32 @@
> ordispider(mod, group=dat$field, lty=3, col="red")
>
> ### Correct hypothesis test (with strata)
-> adonis(Y ~ NO3, data=dat, strata=dat$field, perm=1e3)
+> adonis(Y ~ NO3, data=dat, strata=dat$field, perm=999)
Call:
-adonis(formula = Y ~ NO3, data = dat, permutations = 1000, strata = dat$field)
+adonis(formula = Y ~ NO3, data = dat, permutations = 999, strata = dat$field)
- Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
-NO3 1 0.055856 0.055856 4.0281 0.28714 0.008991 **
-Residuals 10 0.138667 0.013867 0.71286
-Total 11 0.194524 1.00000
+Terms added sequentially (first to last)
+
+ Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
+NO3 1 0.055856 0.055856 4.0281 0.28714 0.009 **
+Residuals 10 0.138667 0.013867 0.71286
+Total 11 0.194524 1.00000
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> ### Incorrect (no strata)
-> adonis(Y ~ NO3, data=dat, perm=1e3)
+> adonis(Y ~ NO3, data=dat, perm=999)
Call:
-adonis(formula = Y ~ NO3, data = dat, permutations = 1000)
+adonis(formula = Y ~ NO3, data = dat, permutations = 999)
- Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
-NO3 1 0.055856 0.055856 4.0281 0.28714 0.004995 **
-Residuals 10 0.138667 0.013867 0.71286
-Total 11 0.194524 1.00000
+Terms added sequentially (first to last)
+
+ Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
+NO3 1 0.055856 0.055856 4.0281 0.28714 0.005 **
+Residuals 10 0.138667 0.013867 0.71286
+Total 11 0.194524 1.00000
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
@@ -750,7 +748,7 @@
Based on 999 permutations
-Empirical upper confidence limits of R:
+Upper quantiles of permutations (null model):
90% 95% 97.5% 99%
0.116 0.160 0.203 0.233
@@ -1094,7 +1092,8 @@
> ### Name: betadisper
> ### Title: Multivariate homogeneity of groups dispersions (variances)
> ### Aliases: betadisper scores.betadisper anova.betadisper plot.betadisper
-> ### boxplot.betadisper TukeyHSD.betadisper ordimedian
+> ### boxplot.betadisper TukeyHSD.betadisper eigenvals.betadisper
+> ### ordimedian
> ### Keywords: methods multivariate hplot
>
> ### ** Examples
@@ -1118,17 +1117,13 @@
No. of Positive Eigenvalues: 15
No. of Negative Eigenvalues: 8
-Average distance to centroid:
+Average distance to medoid:
grazed ungrazed
0.3926 0.2706
Eigenvalues for PCoA axes:
- PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8 PCoA9 PCoA10
- 1.7552 1.1334 0.4429 0.3698 0.2454 0.1961 0.1751 0.1284 0.0972 0.0760
- PCoA11 PCoA12 PCoA13 PCoA14 PCoA15 PCoA16 PCoA17 PCoA18 PCoA19 PCoA20
- 0.0637 0.0583 0.0395 0.0173 0.0051 -0.0004 -0.0065 -0.0133 -0.0254 -0.0375
- PCoA21 PCoA22 PCoA23
--0.0480 -0.0537 -0.0741
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.7552 1.1334 0.4429 0.3698 0.2454 0.1961 0.1751 0.1284
>
> ## Perform test
> anova(mod)
@@ -1191,7 +1186,98 @@
> ## Draw a boxplot of the distances to centroid for each group
> boxplot(mod)
>
+> ## `scores` and `eigenvals` also work
+> scrs <- scores(mod)
+> str(scrs)
+List of 2
+ $ sites : num [1:24, 1:2] 0.0946 -0.3125 -0.3511 -0.3291 -0.1926 ...
+ ..- attr(*, "dimnames")=List of 2
+ .. ..$ : chr [1:24] "18" "15" "24" "27" ...
+ .. ..$ : chr [1:2] "PCoA1" "PCoA2"
+ $ centroids: num [1:2, 1:2] -0.1455 0.2786 0.0758 -0.2111
+ ..- attr(*, "dimnames")=List of 2
+ .. ..$ : chr [1:2] "grazed" "ungrazed"
+ .. ..$ : chr [1:2] "PCoA1" "PCoA2"
+> head(scores(mod, 1:4, display = "sites"))
+ PCoA1 PCoA2 PCoA3 PCoA4
+18 0.09459373 0.15914576 0.074400844 -0.202466025
+15 -0.31248809 0.10032751 -0.062243360 0.110844864
+24 -0.35106507 -0.05954096 -0.038079447 0.095060928
+27 -0.32914546 -0.17019348 0.231623720 0.019110623
+23 -0.19259443 -0.01459250 -0.005679372 -0.209718312
+19 -0.06794575 -0.14501690 -0.085645653 0.002431355
+> # group centroids/medoids
+> scores(mod, 1:4, display = "centroids")
+ PCoA1 PCoA2 PCoA3 PCoA4
+grazed -0.1455200 0.07584572 -0.01366220 -0.0178990
+ungrazed 0.2786095 -0.21114993 -0.03475586 0.0220129
+> # eigenvalues from the underlying principal coordinates analysis
+> eigenvals(mod)
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7
+ 1.7552165 1.1334455 0.4429018 0.3698054 0.2453532 0.1960921 0.1751131
+ PCoA8 PCoA9 PCoA10 PCoA11 PCoA12 PCoA13 PCoA14
+ 0.1284467 0.0971594 0.0759601 0.0637178 0.0583225 0.0394934 0.0172699
+ PCoA15 PCoA16 PCoA17 PCoA18 PCoA19 PCoA20 PCoA21
+ 0.0051011 -0.0004131 -0.0064654 -0.0133147 -0.0253944 -0.0375105 -0.0480069
+ PCoA22 PCoA23
+-0.0537146 -0.0741390
+>
+> ## try out bias correction; compare with mod3
+> (mod3B <- betadisper(dis, groups, type = "median", bias.adjust=TRUE))
+
+ Homogeneity of multivariate dispersions
+
+Call: betadisper(d = dis, group = groups, type = "median", bias.adjust
+= TRUE)
+
+No. of Positive Eigenvalues: 15
+No. of Negative Eigenvalues: 8
+
+Average distance to medoid:
+ grazed ungrazed
+ 0.4055 0.2893
+
+Eigenvalues for PCoA axes:
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.7552 1.1334 0.4429 0.3698 0.2454 0.1961 0.1751 0.1284
+>
+> ## should always work for a single group
+> group <- factor(rep("grazed", NROW(varespec)))
+> (tmp <- betadisper(dis, group, type = "median"))
+
+ Homogeneity of multivariate dispersions
+
+Call: betadisper(d = dis, group = group, type = "median")
+
+No. of Positive Eigenvalues: 15
+No. of Negative Eigenvalues: 8
+
+Average distance to medoid:
+grazed
+0.4255
+
+Eigenvalues for PCoA axes:
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.7552 1.1334 0.4429 0.3698 0.2454 0.1961 0.1751 0.1284
+> (tmp <- betadisper(dis, group, type = "centroid"))
+
+ Homogeneity of multivariate dispersions
+
+Call: betadisper(d = dis, group = group, type = "centroid")
+
+No. of Positive Eigenvalues: 15
+No. of Negative Eigenvalues: 8
+
+Average distance to centroid:
+grazed
+0.4261
+
+Eigenvalues for PCoA axes:
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.7552 1.1334 0.4429 0.3698 0.2454 0.1961 0.1751 0.1284
+>
> ## simulate missing values in 'd' and 'group'
+> ## using spatial medians
> groups[c(2,20)] <- NA
> dis[c(2, 20)] <- NA
> mod2 <- betadisper(dis, groups) ## warnings
@@ -1208,15 +1294,13 @@
No. of Positive Eigenvalues: 14
No. of Negative Eigenvalues: 5
-Average distance to centroid:
+Average distance to medoid:
grazed ungrazed
0.3984 0.3008
Eigenvalues for PCoA axes:
- PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8 PCoA9 PCoA10
- 1.4755 0.8245 0.4218 0.3456 0.2159 0.1688 0.1150 0.1060 0.0912 0.0639
- PCoA11 PCoA12 PCoA13 PCoA14 PCoA15 PCoA16 PCoA17 PCoA18 PCoA19
- 0.0420 0.0267 0.0157 0.0020 -0.0025 -0.0215 -0.0221 -0.0486 -0.0592
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.4755 0.8245 0.4218 0.3456 0.2159 0.1688 0.1150 0.1060
> permutest(mod2, control = permControl(nperm = 100))
Permutation test for homogeneity of multivariate dispersions
@@ -1245,30 +1329,28 @@
> boxplot(mod2)
> plot(TukeyHSD(mod2))
>
-> ## Using spatial median
-> mod3 <- betadisper(dis, groups, type = "median")
-Warning in betadisper(dis, groups, type = "median") :
+> ## Using group centroids
+> mod3 <- betadisper(dis, groups, type = "centroid")
+Warning in betadisper(dis, groups, type = "centroid") :
Missing observations due to 'group' removed.
-Warning in betadisper(dis, groups, type = "median") :
+Warning in betadisper(dis, groups, type = "centroid") :
Missing observations due to 'd' removed.
> mod3
Homogeneity of multivariate dispersions
-Call: betadisper(d = dis, group = groups, type = "median")
+Call: betadisper(d = dis, group = groups, type = "centroid")
No. of Positive Eigenvalues: 14
No. of Negative Eigenvalues: 5
Average distance to centroid:
grazed ungrazed
- 0.3984 0.3008
+ 0.4001 0.3108
Eigenvalues for PCoA axes:
- PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8 PCoA9 PCoA10
- 1.4755 0.8245 0.4218 0.3456 0.2159 0.1688 0.1150 0.1060 0.0912 0.0639
- PCoA11 PCoA12 PCoA13 PCoA14 PCoA15 PCoA16 PCoA17 PCoA18 PCoA19
- 0.0420 0.0267 0.0157 0.0020 -0.0025 -0.0215 -0.0221 -0.0486 -0.0592
+ PCoA1 PCoA2 PCoA3 PCoA4 PCoA5 PCoA6 PCoA7 PCoA8
+1.4755 0.8245 0.4218 0.3456 0.2159 0.1688 0.1150 0.1060
> permutest(mod3, control = permControl(nperm = 100))
Permutation test for homogeneity of multivariate dispersions
@@ -1283,22 +1365,27 @@
Mirrored permutations for Samples?: No
Response: Distances
- Df Sum Sq Mean Sq F N.Perm Pr(>F)
-Groups 1 0.039979 0.039979 2.4237 100 0.1287
-Residuals 18 0.296910 0.016495
+ Df Sum Sq Mean Sq F N.Perm Pr(>F)
+Groups 1 0.033468 0.033468 3.1749 100 0.06931 .
+Residuals 18 0.189749 0.010542
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> anova(mod3)
Analysis of Variance Table
Response: Distances
- Df Sum Sq Mean Sq F value Pr(>F)
-Groups 1 0.039979 0.039979 2.4237 0.1369
-Residuals 18 0.296910 0.016495
+ Df Sum Sq Mean Sq F value Pr(>F)
+Groups 1 0.033468 0.033468 3.1749 0.09166 .
+Residuals 18 0.189749 0.010542
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> plot(mod3)
> boxplot(mod3)
> plot(TukeyHSD(mod3))
>
>
>
+>
> cleanEx()
> nameEx("betadiver")
> ### * betadiver
@@ -1516,10 +1603,10 @@
0.5413 0.3265 0.1293
Eigenvalues for unconstrained axes:
- MDS1 MDS2 MDS3 MDS4 MDS5 MDS6 MDS7 MDS8
-0.906518 0.512743 0.337915 0.262598 0.203220 0.161762 0.124174 0.085570
- MDS9 MDS10 MDS11 MDS12 MDS13 MDS14 MDS15
-0.068881 0.058346 0.050083 0.027738 0.020839 0.007306 0.001345
+ MDS1 MDS2 MDS3 MDS4 MDS5 MDS6 MDS7 MDS8 MDS9 MDS10 MDS11
+0.9065 0.5127 0.3379 0.2626 0.2032 0.1618 0.1242 0.0856 0.0689 0.0583 0.0501
+ MDS12 MDS13 MDS14 MDS15
+0.0277 0.0208 0.0073 0.0013
> plot(vare.cap)
> anova(vare.cap)
@@ -1553,6 +1640,8 @@
1.4408 0.8523 0.6015 0.4888 0.4187 0.3538 0.2877 0.2160
(Showed only 8 of all 19 unconstrained eigenvalues)
+Constant added to distances: 0.2614286
+
> ## Avoid negative eigenvalues by taking square roots of dissimilarities
> capscale(varespec ~ N + P + K + Condition(Al), varechem,
+ dist = "bray", sqrt.dist= TRUE)
@@ -1660,16 +1749,14 @@
Inertia is mean squared contingency coefficient
Eigenvalues for constrained axes:
- CCA1 CCA2 CCA3 CCA4 CCA5 CCA6 CCA7 CCA8
-0.438870 0.291775 0.162847 0.142130 0.117952 0.089029 0.070295 0.058359
- CCA9 CCA10 CCA11 CCA12 CCA13 CCA14
-0.031141 0.013294 0.008364 0.006538 0.006156 0.004733
+ CCA1 CCA2 CCA3 CCA4 CCA5 CCA6 CCA7 CCA8 CCA9 CCA10 CCA11
+0.4389 0.2918 0.1628 0.1421 0.1180 0.0890 0.0703 0.0584 0.0311 0.0133 0.0084
+ CCA12 CCA13 CCA14
+0.0065 0.0062 0.0047
Eigenvalues for unconstrained axes:
- CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
-0.197765 0.141926 0.101174 0.070787 0.053303 0.033299 0.018868 0.015104
- CA9
-0.009488
+ CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8 CA9
+0.19776 0.14193 0.10117 0.07079 0.05330 0.03330 0.01887 0.01510 0.00949
> plot(vare.cca)
> ## Formula interface and a better model
@@ -1685,8 +1772,8 @@
Inertia is mean squared contingency coefficient
Eigenvalues for constrained axes:
- CCA1 CCA2 CCA3 CCA4 CCA5 CCA6
-0.37563 0.23419 0.14067 0.13229 0.10675 0.05614
+ CCA1 CCA2 CCA3 CCA4 CCA5 CCA6
+0.3756 0.2342 0.1407 0.1323 0.1068 0.0561
Eigenvalues for unconstrained axes:
CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
@@ -1705,12 +1792,12 @@
Inertia is mean squared contingency coefficient
Eigenvalues for constrained axes:
- CCA1
-0.1572
+ CCA1
+0.15722
Eigenvalues for unconstrained axes:
- CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
-0.47455 0.29389 0.21403 0.19541 0.17482 0.11711 0.11207 0.08797
+ CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
+0.4745 0.2939 0.2140 0.1954 0.1748 0.1171 0.1121 0.0880
(Showed only 8 of all 22 unconstrained eigenvalues)
> cca(varespec ~ Ca + Condition(pH), varechem)
@@ -1724,12 +1811,12 @@
Inertia is mean squared contingency coefficient
Eigenvalues for constrained axes:
- CCA1
-0.1827
+ CCA1
+0.18269
Eigenvalues for unconstrained axes:
- CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
-0.38343 0.27487 0.21233 0.17599 0.17013 0.11613 0.10892 0.08797
+ CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
+0.3834 0.2749 0.2123 0.1760 0.1701 0.1161 0.1089 0.0880
(Showed only 8 of all 21 unconstrained eigenvalues)
> ## RDA
@@ -1790,18 +1877,18 @@
>
> ### Name: clamtest
> ### Title: Multinomial Species Classification Method (CLAM)
-> ### Aliases: clamtest summary.clamtest print.summary.clamtest plot.clamtest
+> ### Aliases: clamtest summary.clamtest plot.clamtest
> ### Keywords: htest
>
> ### ** Examples
>
> data(mite)
> data(mite.env)
-> x <- clamtest(mite, mite.env$Shrub=="None", alpha=0.005, specialization = 0.667)
-> summary(x)
+> sol <- clamtest(mite, mite.env$Shrub=="None", alpha=0.005)
+> summary(sol)
Two Groups Species Classification Method (CLAM)
-Specialization threshold = 0.667
+Specialization threshold = 0.6666667
Alpha level = 0.005
Estimated sample coverage:
@@ -1817,7 +1904,7 @@
Specialist_FALSE 14 0.400
Specialist_TRUE 4 0.114
Too_rare 7 0.200
-> head(x)
+> head(sol)
Species Total_FALSE Total_TRUE Classes
1 Brachy 534 77 Generalist
2 PHTH 89 0 Specialist_FALSE
@@ -1825,7 +1912,7 @@
4 RARD 85 0 Specialist_FALSE
5 SSTR 22 0 Too_rare
6 Protopl 26 0 Too_rare
-> plot(x)
+> plot(sol)
>
>
>
@@ -1847,20 +1934,20 @@
+ array(replicate(n, sample(x)), c(dim(x), n))
> (cs <- commsim("r00", fun=f, binary=TRUE,
+ isSeq=FALSE, mode="integer"))
-An object of class "commsim"
-"r00" method (binary, non-sequential, integer mode)
+An object of class “commsim”
+‘r00’ method (binary, non-sequential, integer mode)
>
> ## retrieving the sequential swap algorithm
> (cs <- make.commsim("swap"))
-An object of class "commsim"
-"swap" method (binary, sequential, integer mode)
+An object of class “commsim”
+‘swap’ method (binary, sequential, integer mode)
>
> ## feeding a commsim object as argument
> make.commsim(cs)
-An object of class "commsim"
-"swap" method (binary, sequential, integer mode)
+An object of class “commsim”
+‘swap’ method (binary, sequential, integer mode)
>
> ## structural constraints
@@ -1877,46 +1964,46 @@
+ y <- simulate(m, nsim=n)
+ out <- rowMeans(sapply(1:dim(y)[3],
+ function(i) diagfun(attr(y, "data"), y[,,i])))
-+ z <- as.numeric(c(attr(y, "binary"), attr(y, "isSeq")))
-+ names(z) <- c("binary", "isSeq")
++ z <- as.numeric(c(attr(y, "binary"), attr(y, "isSeq"),
++ attr(y, "mode") == "double"))
++ names(z) <- c("binary", "isSeq", "double")
+ c(z, out)
+ }
> x <- matrix(rbinom(10*12, 1, 0.5)*rpois(10*12, 3), 12, 10)
> algos <- make.commsim()
> a <- t(sapply(algos, evalfun, x=x, n=10))
-Warning in storage.mode(state) <- object$commsim$mode :
- NAs introduced by coercion
> print(as.table(ifelse(a==1,1,0)), zero.print = ".")
- binary isSeq sum fill rowSums colSums rowFreq colFreq
-r00 1 . 1 1 . . . .
-c0 1 . 1 1 . 1 . 1
-r0 1 . 1 1 1 . 1 .
-r1 1 . 1 1 1 . 1 .
-r2 1 . 1 1 1 . 1 .
-quasiswap 1 . 1 1 1 1 1 1
-swap 1 1 1 1 1 1 1 1
-tswap 1 1 1 1 1 1 1 1
-backtrack 1 . 1 1 1 1 1 1
-r2dtable . . 1 . 1 1 . .
-swap_count . 1 . . . . . .
-quasiswap_count . . 1 1 1 1 . .
-swsh_samp . . 1 1 . . 1 1
-swsh_both . . 1 1 . . 1 1
-swsh_samp_r . . 1 1 1 . 1 1
-swsh_samp_c . . 1 1 . 1 1 1
-swsh_both_r . . 1 1 1 . 1 1
-swsh_both_c . . 1 1 . 1 1 1
-abuswap_r . 1 1 1 1 . 1 1
-abuswap_c . 1 1 1 . 1 1 1
-r00_samp . . 1 1 . . . .
-c0_samp . . 1 1 . 1 . 1
-r0_samp . . 1 1 1 . 1 .
-r00_ind . . 1 . . . . .
-c0_ind . . 1 . . 1 . .
-r0_ind . . 1 . 1 . . .
-r00_both . . 1 1 . . . .
-c0_both . . 1 1 . 1 . 1
-r0_both . . 1 1 1 . 1 .
+ binary isSeq double sum fill rowSums colSums rowFreq colFreq
+r00 1 . . 1 1 . . . .
+c0 1 . . 1 1 . 1 . 1
+r0 1 . . 1 1 1 . 1 .
+r0_old 1 . . 1 1 1 . 1 .
+r1 1 . . 1 1 1 . 1 .
+r2 1 . . 1 1 1 . 1 .
+quasiswap 1 . . 1 1 1 1 1 1
+swap 1 1 . 1 1 1 1 1 1
+tswap 1 1 . 1 1 1 1 1 1
+backtrack 1 . . 1 1 1 1 1 1
+r2dtable . . . 1 . 1 1 . .
+swap_count . 1 . 1 1 1 1 . .
+quasiswap_count . . . 1 1 1 1 . .
+swsh_samp . . . 1 1 . . 1 1
+swsh_both . . . 1 1 . . 1 1
+swsh_samp_r . . . 1 1 1 . 1 1
+swsh_samp_c . . . 1 1 . 1 1 1
+swsh_both_r . . . 1 1 1 . 1 1
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/vegan -r 2488
More information about the Vegan-commits
mailing list