[Vegan-commits] r1467 - in pkg/vegan: . tests tests/Examples
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Jan 20 15:24:15 CET 2011
Author: jarioksa
Date: 2011-01-20 15:24:15 +0100 (Thu, 20 Jan 2011)
New Revision: 1467
Added:
pkg/vegan/tests/
pkg/vegan/tests/Examples/
pkg/vegan/tests/Examples/vegan-Ex.Rout.save
pkg/vegan/tests/vegan-tests.R
pkg/vegan/tests/vegan-tests.Rout.save
Log:
add tests/ with anova.cca unit tests and Examples
Added: pkg/vegan/tests/Examples/vegan-Ex.Rout.save
===================================================================
--- pkg/vegan/tests/Examples/vegan-Ex.Rout.save (rev 0)
+++ pkg/vegan/tests/Examples/vegan-Ex.Rout.save 2011-01-20 14:24:15 UTC (rev 1467)
@@ -0,0 +1,7971 @@
+
+R version 2.12.1 (2010-12-16)
+Copyright (C) 2010 The R Foundation for Statistical Computing
+ISBN 3-900051-07-0
+Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
+
+R is free software and comes with ABSOLUTELY NO WARRANTY.
+You are welcome to redistribute it under certain conditions.
+Type 'license()' or 'licence()' for distribution details.
+
+ Natural language support but running in an English locale
+
+R is a collaborative project with many contributors.
+Type 'contributors()' for more information and
+'citation()' on how to cite R or R packages in publications.
+
+Type 'demo()' for some demos, 'help()' for on-line help, or
+'help.start()' for an HTML browser interface to help.
+Type 'q()' to quit R.
+
+> pkgname <- "vegan"
+> source(file.path(R.home("share"), "R", "examples-header.R"))
+> options(warn = 1)
+> library('vegan')
+This is vegan 1.18-22
+>
+> assign(".oldSearch", search(), pos = 'CheckExEnv')
+> cleanEx()
+> nameEx("BCI")
+> ### * BCI
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: BCI
+> ### Title: Barro Colorado Island Tree Counts
+> ### Aliases: BCI
+> ### Keywords: datasets
+>
+> ### ** Examples
+>
+> data(BCI)
+> ## UTM Coordinates (in metres)
+> UTM.EW <- rep(seq(625754, 626654, by=100), each=5)
+> UTM.NS <- rep(seq(1011569, 1011969, by=100), len=50)
+>
+>
+>
+> cleanEx()
+> nameEx("CCorA")
+> ### * CCorA
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: CCorA
+> ### Title: Canonical Correlation Analysis
+> ### Aliases: CCorA print.CCorA biplot.CCorA
+> ### Keywords: multivariate
+>
+> ### ** Examples
+>
+> # Example using two mite groups. The mite data are available in vegan
+> data(mite)
+> # Two mite species associations (Legendre 2005, Fig. 4)
+> group.1 <- c(1,2,4:8,10:15,17,19:22,24,26:30)
+> group.2 <- c(3,9,16,18,23,25,31:35)
+> # Separate Hellinger transformations of the two groups of species
+> mite.hel.1 <- decostand(mite[,group.1], "hel")
+> mite.hel.2 <- decostand(mite[,group.2], "hel")
+> rownames(mite.hel.1) = paste("S",1:nrow(mite),sep="")
+> rownames(mite.hel.2) = paste("S",1:nrow(mite),sep="")
+> out <- CCorA(mite.hel.1, mite.hel.2)
+> out
+
+Canonical Correlation Analysis
+
+Call:
+CCorA(Y = mite.hel.1, X = mite.hel.2)
+
+ Y X
+Matrix Ranks 24 11
+
+Pillai's trace: 4.573009
+
+Significance of Pillai's trace:
+from F-distribution: 0.0032737
+
+ CanAxis1 CanAxis2 CanAxis3 CanAxis4 CanAxis5 CanAxis6
+Canonical Correlations 0.92810 0.82431 0.81209 0.74981 0.70795 0.65950
+ CanAxis7 CanAxis8 CanAxis9 CanAxis10 CanAxis11
+Canonical Correlations 0.50189 0.48179 0.41089 0.37823 0.28
+
+ Y | X X | Y
+RDA R squares 0.33224 0.5376
+adj. RDA R squares 0.20560 0.2910
+
+> summary(out)
+ Length Class Mode
+Pillai 1 -none- numeric
+Eigenvalues 11 -none- numeric
+CanCorr 11 -none- numeric
+Mat.ranks 2 -none- numeric
+RDA.Rsquares 2 -none- numeric
+RDA.adj.Rsq 2 -none- numeric
+nperm 1 -none- numeric
+p.Pillai 1 -none- numeric
+p.perm 1 -none- logical
+Cy 770 -none- numeric
+Cx 770 -none- numeric
+corr.Y.Cy 264 -none- numeric
+corr.X.Cx 121 -none- numeric
+corr.Y.Cx 264 -none- numeric
+corr.X.Cy 121 -none- numeric
+call 3 -none- call
+> biplot(out, "ob") # Two plots of objects
+> biplot(out, "v", cex=c(0.7,0.6)) # Two plots of variables
+> biplot(out, "ov", cex=c(0.7,0.6)) # Four plots (2 for objects, 2 for variables)
+> biplot(out, "b", cex=c(0.7,0.6)) # Two biplots
+> biplot(out, xlabs = NA, plot.axes = c(3,5)) # Plot axes 3, 5. No object names
+> biplot(out, plot.type="biplots", xlabs = NULL) # Replace object names by numbers
+>
+> # Example using random numbers. No significant relationship is expected
+> mat1 <- matrix(rnorm(60),20,3)
+> mat2 <- matrix(rnorm(100),20,5)
+> out2 = CCorA(mat1, mat2, nperm=99)
+> out2
+
+Canonical Correlation Analysis
+
+Call:
+CCorA(Y = mat1, X = mat2, nperm = 99)
+
+ Y X
+Matrix Ranks 3 5
+
+Pillai's trace: 0.6455578
+
+Significance of Pillai's trace:
+based on 99 permutations: 0.69
+from F-distribution: 0.70352
+
+ CanAxis1 CanAxis2 CanAxis3
+Canonical Correlations 0.69691 0.38140 0.12
+
+ Y | X X | Y
+RDA R squares 0.17066 0.1368
+adj. RDA R squares -0.12553 -0.0250
+
+> biplot(out2, "b")
+>
+>
+>
+> cleanEx()
+> nameEx("MOStest")
+> ### * MOStest
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: MOStest
+> ### Title: Mitchell-Olds & Shaw Test for the Location of Quadratic Extreme
+> ### Aliases: MOStest print.MOStest plot.MOStest fieller.MOStest
+> ### profile.MOStest confint.MOStest
+> ### Keywords: models regression
+>
+> ### ** Examples
+>
+> ## The Al-Mufti data analysed in humpfit():
+> mass <- c(140,230,310,310,400,510,610,670,860,900,1050,1160,1900,2480)
+> spno <- c(1, 4, 3, 9, 18, 30, 20, 14, 3, 2, 3, 2, 5, 2)
+> mod <- MOStest(mass, spno)
+> ## Insignificant
+> mod
+
+Mitchell-Olds and Shaw test
+Null: hump of a quadratic linear predictor is at min or max
+
+Family: gaussian
+Link function: identity
+
+ hump min max
+ 46.89749 140.00000 2480.00000
+***** Caution: hump/pit not bracketed by the data ******
+
+ min/max F Pr(>F)
+hump at min 140 0.0006 0.9816
+hump at max 2480 0.3161 0.5852
+Combined 0.9924
+> ## ... but inadequate shape of the curve
+> op <- par(mfrow=c(2,2), mar=c(4,4,1,1)+.1)
+> plot(mod)
+> ## Looks rather like log-link with Poisson error and logarithmic biomass
+> mod <- MOStest(log(mass), spno, family=quasipoisson)
+> mod
+
+Mitchell-Olds and Shaw test
+Null: hump of a quadratic linear predictor is at min or max
+
+Family: quasipoisson
+Link function: log
+
+ min hump max
+4.941642 6.243371 7.816014
+
+ min/max F Pr(>F)
+hump at min 4.9416 7.1367 0.02174 *
+hump at max 7.8160 9.0487 0.01191 *
+Combined 0.03338 *
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> plot(mod)
+> par(op)
+> ## Confidence Limits
+> fieller.MOStest(mod)
+ 2.5 % 97.5 %
+5.255827 6.782979
+> confint(mod)
+Loading required package: MASS
+ 2.5 % 97.5 %
+5.816021 6.574378
+> plot(profile(mod))
+>
+>
+>
+> graphics::par(get("par.postscript", pos = 'CheckExEnv'))
+> cleanEx()
+
+detaching ‘package:MASS’
+
+> nameEx("RsquareAdj")
+> ### * RsquareAdj
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: RsquareAdj
+> ### Title: Adjusted R-square
+> ### Aliases: RsquareAdj RsquareAdj.default RsquareAdj.rda RsquareAdj.cca
+> ### RsquareAdj.lm RsquareAdj.glm RsquareAdj
+> ### Keywords: univar multivariate
+>
+> ### ** Examples
+>
+> data(mite)
+> data(mite.env)
+> ## rda
+> m <- rda(decostand(mite, "hell") ~ ., mite.env)
+> RsquareAdj(m)
+$r.squared
+[1] 0.5265047
+
+$adj.r.squared
+[1] 0.4367038
+
+> ## default method
+> RsquareAdj(0.8, 20, 5)
+[1] 0.7285714
+>
+>
+>
+> cleanEx()
+> nameEx("SSarrhenius")
+> ### * SSarrhenius
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: SSarrhenius
+> ### Title: Self-Starting nls Species-Area Models
+> ### Aliases: SSarrhenius SSlomolino SSgitay SSgleason
+> ### Keywords: models
+>
+> ### ** Examples
+>
+> ## Get species area data: sipoo.area gives the areas of islands
+> example(sipoo)
+
+sipoo> data(sipoo)
+
+sipoo> ## Areas of the islands in hectares
+sipoo> sipoo.area <- c(1.1, 2.1, 2.2, 3.1, 3.5, 5.8, 6, 6.1, 6.5, 11.4, 13,
+sipoo+ 14.5, 16.1 ,17.5, 28.7, 40.5, 104.5, 233)
+> S <- specnumber(sipoo)
+> plot(S ~ sipoo.area, xlab = "Island Area (ha)", ylab = "Number of Species",
++ ylim = c(1, max(S)))
+> ## The Arrhenius model
+> marr <- nls(S ~ SSarrhenius(sipoo.area, k, z))
+> marr
+Nonlinear regression model
+ model: S ~ SSarrhenius(sipoo.area, k, z)
+ data: parent.frame()
+ k z
+3.4062 0.4364
+ residual sum-of-squares: 78.1
+
+Number of iterations to convergence: 5
+Achieved convergence tolerance: 1.056e-06
+> ## confidence limits from profile likelihood
+> confint(marr)
+Waiting for profiling to be done...
+ 2.5% 97.5%
+k 2.6220312 4.3033906
+z 0.3813576 0.4944693
+> ## draw a line
+> xtmp <- seq(min(sipoo.area), max(sipoo.area), len=51)
+> lines(xtmp, predict(marr, newdata=data.frame(sipoo.area = xtmp)), lwd=2)
+> ## The normal way is to use linear regression on log-log data,
+> ## but this will be different from the previous:
+> mloglog <- lm(log(S) ~ log(sipoo.area))
+> mloglog
+
+Call:
+lm(formula = log(S) ~ log(sipoo.area))
+
+Coefficients:
+ (Intercept) log(sipoo.area)
+ 1.0111 0.4925
+
+> lines(xtmp, exp(predict(mloglog, newdata=data.frame(sipoo.area=xtmp))),
++ lty=2)
+> ## Gleason: log-linear
+> mgle <- nls(S ~ SSgleason(sipoo.area, k, slope))
+> lines(xtmp, predict(mgle, newdata=data.frame(sipoo.area=xtmp)),
++ lwd=2, col=2)
+> ## Gitay: quadratic of log-linear
+> mgit <- nls(S ~ SSgitay(sipoo.area, k, slope))
+> lines(xtmp, predict(mgit, newdata=data.frame(sipoo.area=xtmp)),
++ lwd=2, col = 3)
+> ## Lomolino: using original names of the parameters (Lomolino 2000):
+> mlom <- nls(S ~ SSlomolino(sipoo.area, Smax, A50, Hill))
+> mlom
+Nonlinear regression model
+ model: S ~ SSlomolino(sipoo.area, Smax, A50, Hill)
+ data: parent.frame()
+ Smax A50 Hill
+53.493 94.697 2.018
+ residual sum-of-squares: 55.37
+
+Number of iterations to convergence: 6
+Achieved convergence tolerance: 9.715e-07
+> lines(xtmp, predict(mlom, newdata=data.frame(sipoo.area=xtmp)),
++ lwd=2, col = 4)
+> ## One canned model of standard R:
+> mmic <- nls(S ~ SSmicmen(sipoo.area, slope, Asym))
+> lines(xtmp, predict(mmic, newdata = data.frame(sipoo.area=xtmp)),
++ lwd =2, col = 5)
+> legend("bottomright", c("Arrhenius", "log-log linear", "Gleason", "Gitay",
++ "Lomolino", "Michaelis-Menten"), col=c(1,1,2,3,4,5), lwd=c(2,1,2,2,2,2),
++ lty=c(1,2,1,1,1,1))
+> ## compare models (AIC)
+> allmods <- list(Arrhenius = marr, Gleason = mgle, Gitay = mgit,
++ Lomolino = mlom, MicMen= mmic)
+>
+>
+>
+> cleanEx()
+> nameEx("add1.cca")
+> ### * add1.cca
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: add1.cca
+> ### Title: Add or Drop Single Terms to a Constrained Ordination Model
+> ### Aliases: add1.cca drop1.cca
+> ### Keywords: multivariate models
+>
+> ### ** Examples
+>
+> data(dune)
+> data(dune.env)
+> ## Automatic model building based on AIC but with permutation tests
+> step(cca(dune ~ 1, dune.env), reformulate(names(dune.env)), test="perm")
+Start: AIC=87.66
+dune ~ 1
+
+ Df AIC F N.Perm Pr(>F)
++ Moisture 3 86.608 2.2536 199 0.010 **
++ Management 3 86.935 2.1307 199 0.005 **
++ A1 1 87.411 2.1400 199 0.020 *
+<none> 87.657
++ Manure 4 88.832 1.5251 199 0.025 *
++ Use 2 89.134 1.1431 99 0.250
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+
+Step: AIC=86.61
+dune ~ Moisture
+
+ Df AIC F N.Perm Pr(>F)
+<none> 86.608
++ Management 3 86.813 1.4565 199 0.020 *
++ A1 1 86.992 1.2624 99 0.150
++ Use 2 87.259 1.2760 199 0.120
++ Manure 4 87.342 1.3143 199 0.090 .
+- Moisture 3 87.657 2.2536 199 0.005 **
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+Call: cca(formula = dune ~ Moisture, data = dune.env)
+
+ Inertia Proportion Rank
+Total 2.1153 1.0000
+Constrained 0.6283 0.2970 3
+Unconstrained 1.4870 0.7030 16
+Inertia is mean squared contingency coefficient
+
+Eigenvalues for constrained axes:
+ 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
+
+> ## 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
+
+> ## Manual model building
+> ## -- define the maximal model for scope
+> mbig <- rda(dune ~ ., dune.env)
+> ## -- define an empty model to start with
+> m0 <- rda(dune ~ 1, dune.env)
+> ## -- manual selection and updating
+> 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 .
+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
+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 .
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> m0 <- update(m0, . ~ . + Moisture)
+> ## -- included variables still significant?
+> 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 **
+---
+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
+>
+>
+>
+> cleanEx()
+> nameEx("adipart")
+> ### * adipart
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: adipart
+> ### Title: Additive Diversity Partitioning and Hierarchical Null Model
+> ### Testing
+> ### Aliases: adipart print.adipart hiersimu print.hiersimu
+> ### Keywords: multivariate
+>
+> ### ** Examples
+>
+> ## NOTE: 'nsimul' argument usually needs to be >= 99
+> ## here much lower value is used for demonstration
+>
+> data(mite)
+> data(mite.xy)
+> data(mite.env)
+> ## Function to get equal area partitions of the mite data
+> cutter <- function (x, cut = seq(0, 10, by = 2.5)) {
++ 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))}
+> ## The hierarchy of sample aggregation
+> levsm <- data.frame(
++ l1=as.factor(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)))
+> ## Let's see in a map
+> par(mfrow=c(1,3))
+> plot(mite.xy, main="l1", col=as.numeric(levsm$l1)+1)
+> plot(mite.xy, main="l2", col=as.numeric(levsm$l2)+1)
+> 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=9)
+adipart with 9 simulations
+with index richness, weights unif
+
+ statistic z 2.5% 50% 97.5% Pr(sim.)
+alpha.1 15.114 -44.931 22.217 22.386 22.594 0.1
+alpha.2 29.750 -28.284 34.500 34.750 35.000 0.1
+alpha.3 33.000 0.000 35.000 35.000 35.000 0.1
+gamma 35.000 0.000 35.000 35.000 35.000 1.0
+beta.1 14.636 11.309 12.156 12.286 12.710 0.1
+beta.2 3.250 16.971 0.000 0.250 0.500 0.1
+beta.3 2.000 0.000 0.000 0.000 0.000 0.1
+> ## Hierarchical null model testing
+> ## diversity analysis (similar to adipart)
+> hiersimu(mite ~., levsm, diversity, relative=TRUE, nsimul=9)
+hiersimu with 9 simulations
+
+ statistic z 2.5% 50% 97.5% Pr(sim.)
+l1 0.76064 -61.57993 0.93421 0.93718 0.9418 0.1
+l2 0.89736 -98.45713 0.99629 0.99723 0.9992 0.1
+l3 0.92791 -550.28201 0.99914 0.99931 0.9995 0.1
+l4 1.00000 0.00000 1.00000 1.00000 1.0000 1.0
+> ## Hierarchical testing with the Morisita index
+> morfun <- function(x) dispindmorisita(x)$imst
+> hiersimu(mite ~., levsm, morfun, drop.highest=TRUE, nsimul=9)
+hiersimu with 9 simulations
+
+ statistic z 2.5% 50% 97.5% Pr(sim.)
+l1 0.520702 9.445924 0.336022 0.363409 0.3815 0.1
+l2 0.602337 9.398210 0.098335 0.135400 0.2364 0.1
+l3 0.675085 18.347360 -0.293851 -0.189212 -0.1565 0.1
+>
+>
+>
+> graphics::par(get("par.postscript", pos = 'CheckExEnv'))
+> cleanEx()
+> nameEx("adonis")
+> ### * adonis
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: adonis
+> ### Title: Permutational Multivariate Analysis of Variance Using Distance
+> ### Matrices
+> ### Aliases: adonis print.adonis
+> ### Keywords: multivariate nonparametric
+>
+> ### ** Examples
+>
+> data(dune)
+> data(dune.env)
+> adonis(dune ~ Management*A1, data=dune.env, permutations=99)
+
+Call:
+adonis(formula = dune ~ Management * A1, data = dune.env, permutations = 99)
+
+ 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 *
+Management:A1 3 0.5892 0.19639 1.3090 0.13705 0.18
+Residuals 12 1.8004 0.15003 0.41878
+Total 19 4.2990 1.00000
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+>
+>
+> ### Example of use with strata, for nested (e.g., block) designs.
+>
+> dat <- expand.grid(rep=gl(2,1), NO3=factor(c(0,10)),field=gl(3,1) )
+> dat
+ rep NO3 field
+1 1 0 1
+2 2 0 1
+3 1 10 1
+4 2 10 1
+5 1 0 2
+6 2 0 2
+7 1 10 2
+8 2 10 2
+9 1 0 3
+10 2 0 3
+11 1 10 3
+12 2 10 3
+> Agropyron <- with(dat, as.numeric(field) + as.numeric(NO3)+2) +rnorm(12)/2
+> Schizachyrium <- with(dat, as.numeric(field) - as.numeric(NO3)+2) +rnorm(12)/2
+> total <- Agropyron + Schizachyrium
+> library(lattice)
+> dotplot(total ~ NO3, dat, jitter.x=TRUE, groups=field,
++ type=c('p','a'), xlab="NO3", auto.key=list(columns=3, lines=TRUE) )
+>
+> Y <- data.frame(Agropyron, Schizachyrium)
+> mod <- metaMDS(Y)
+Loading required package: MASS
+Run 0 stress 8.56469
+Run 1 stress 17.68491
+Run 2 stress 8.556588
+... New best solution
+... procrustes: rmse 0.001526914 max resid 0.003279392
+*** Solution reached
+
+> plot(mod)
+> ### Hulls show treatment
+> ordihull(mod, group=dat$NO3, show="0")
+> ordihull(mod, group=dat$NO3, show="10", col=3)
+> ### Spider shows fields
+> ordispider(mod, group=dat$field, lty=3, col="red")
+>
+> ### Correct hypothesis test (with strata)
+> adonis(Y ~ NO3, data=dat, strata=dat$field, perm=1e3)
+
+Call:
+adonis(formula = Y ~ NO3, data = dat, permutations = 1000, 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
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+>
+> ### Incorrect (no strata)
+> adonis(Y ~ NO3, data=dat, perm=1e3)
+
+Call:
+adonis(formula = Y ~ NO3, data = dat, permutations = 1000)
+
+ 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
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+>
+>
+>
+> cleanEx()
+
+detaching ‘package:MASS’, ‘package:lattice’
+
+> nameEx("anosim")
+> ### * anosim
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: anosim
+> ### Title: Analysis of Similarities
+> ### Aliases: anosim print.anosim summary.anosim plot.anosim
+> ### Keywords: multivariate nonparametric htest
+>
+> ### ** Examples
+>
+> data(dune)
+> data(dune.env)
+> dune.dist <- vegdist(dune)
+> attach(dune.env)
+> dune.ano <- anosim(dune.dist, Management)
+> summary(dune.ano)
+
+Call:
+anosim(dat = dune.dist, grouping = Management)
+Dissimilarity: bray
+
+ANOSIM statistic R: 0.2579
+ Significance: 0.006
+
+Based on 999 permutations
+
+Empirical upper confidence limits of R:
+ 90% 95% 97.5% 99%
+0.116 0.160 0.203 0.233
+
+Dissimilarity ranks between and within classes:
+ 0% 25% 50% 75% 100% N
+Between 4 58.50 104.00 145.500 188.0 147
+BF 5 15.25 25.50 41.250 57.0 3
+HF 1 7.25 46.25 68.125 89.5 10
+NM 6 64.75 124.50 156.250 181.0 15
+SF 3 32.75 53.50 99.250 184.0 15
+
+> plot(dune.ano)
+Warning in bxp(list(stats = c(4, 58.5, 104, 145.5, 188, 5, 15.25, 25.5, :
+ some notches went outside hinges ('box'): maybe set notch=FALSE
+>
+>
+>
+> cleanEx()
+
+detaching ‘dune.env’
+
+> nameEx("anova.cca")
+> ### * anova.cca
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: anova.cca
+> ### Title: Permutation Test for Constrained Correspondence Analysis,
+> ### Redundancy Analysis and Constrained Analysis of Principal Coordinates
+> ### Aliases: anova.cca anova.ccanull anova.ccabyaxis anova.ccabyterm
+> ### anova.ccabymargin permutest permutest.default permutest.cca
+> ### print.permutest.cca
+> ### Keywords: multivariate htest
+>
+> ### ** Examples
+>
+> data(varespec)
+> data(varechem)
+> vare.cca <- cca(varespec ~ Al + P + K, varechem)
+> ## overall test
+> anova(vare.cca)
+Permutation test for cca under reduced model
+
+Model: cca(formula = varespec ~ Al + P + K, data = varechem)
+ Df Chisq F N.Perm Pr(>F)
+Model 3 0.6441 2.9840 199 0.005 **
+Residual 20 1.4391
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> ## Test for axes
+> anova(vare.cca, by="axis", perm.max=500)
+Model: cca(formula = varespec ~ Al + P + K, data = varechem)
+ Df Chisq F N.Perm Pr(>F)
+CCA1 1 0.3616 5.0249 199 0.005 **
+CCA2 1 0.1700 2.3621 199 0.010 **
+CCA3 1 0.1126 1.5651 399 0.100 .
+Residual 20 1.4391
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> ## Sequential test for terms
+> anova(vare.cca, by="terms", permu=200)
+Permutation test for cca under reduced model
+Terms added sequentially (first to last)
+
+Model: cca(formula = varespec ~ Al + P + K, data = varechem)
+ Df Chisq F N.Perm Pr(>F)
+Al 1 0.2982 4.1440 199 0.005 **
+P 1 0.1899 2.6393 199 0.015 *
+K 1 0.1561 2.1688 199 0.030 *
+Residual 20 1.4391
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> ## Marginal or Type III effects
+> anova(vare.cca, by="margin")
+Permutation test for cca under reduced model
+Marginal effects of terms
+
+Model: cca(formula = varespec ~ Al + P + K, data = varechem)
+ Df Chisq F N.Perm Pr(>F)
+Al 1 0.3118 4.3340 199 0.00500 **
+P 1 0.1681 2.3362 199 0.01500 *
+K 1 0.1561 2.1688 599 0.02833 *
+Residual 20 1.4391
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+> ## Marginal test knows 'scope'
+> anova(vare.cca, by = "m", scope="P")
+Permutation test for cca under reduced model
+Marginal effects of terms
+
+Model: cca(formula = varespec ~ Al + P + K, data = varechem)
+ Df Chisq F N.Perm Pr(>F)
+P 1 0.1681 2.3362 199 0.015 *
+Residual 20 1.4391
+---
+Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+>
+>
+>
+> cleanEx()
+> nameEx("as.mlm")
+> ### * as.mlm
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: as.mlm.cca
+> ### Title: Refit Constrained Ordination as a Multiple Response Linear Model
+> ### Aliases: as.mlm as.mlm.cca as.mlm.rda
+> ### Keywords: models multivariate
+>
+> ### ** Examples
+>
+> data(varespec)
+> data(varechem)
+> mod <- cca(varespec ~ Al + P + K, data=varechem)
+> lmod <- as.mlm(mod)
+> ## Coefficients
+> lmod
+
+Call:
+lm(formula = wa ~ . - 1, data = as.data.frame(X))
+
+Coefficients:
+ CCA1 CCA2 CCA3
+Al 0.007479 -0.001884 0.003381
+P -0.006491 -0.102190 -0.022307
+K -0.006756 0.015344 0.017067
+
+> coef(mod)
+ CCA1 CCA2 CCA3
+Al 0.007478556 -0.001883637 0.003380774
+P -0.006491081 -0.102189737 -0.022306682
+K -0.006755568 0.015343662 0.017067351
+> ## Influential observations
+> influence.measures(lmod)
+Influence measures of
+ lm(formula = wa ~ . - 1, data = as.data.frame(X)) :
+
+ dfb.Al dfb.P dfb.K CCA1 CCA2 CCA3 cov.r CCA1.1
+18 -0.251387 0.00976 -0.06310 0.2740 0.1806 -0.118754 0.0265 7.38e-03
+15 0.099858 0.13864 -0.11781 -0.1654 -0.0935 0.006898 0.0319 2.86e-03
+24 -0.003448 -0.44078 0.20788 -0.4824 -0.1750 -0.260788 0.0307 2.33e-02
+27 0.071682 -0.01707 -0.03516 -0.1018 -0.1676 0.022271 0.0406 1.13e-03
+23 -0.116533 0.06900 -0.02545 0.1441 0.2918 -0.220457 0.0355 2.23e-03
+19 -0.007394 -0.01169 0.01080 0.0136 -0.2318 -0.000417 0.0359 2.02e-05
+22 0.150916 0.14845 -0.13091 -0.2047 0.3815 0.168914 0.0376 4.50e-03
+16 0.107456 0.17900 -0.09917 -0.2120 0.2250 0.194432 0.0338 4.75e-03
+28 0.332161 -0.34398 -0.05414 -0.6774 0.0742 0.620990 0.0364 4.65e-02
+13 0.366880 -1.00834 1.23685 1.3392 0.4102 0.277067 0.1124 1.89e-01
+14 0.024147 0.02512 -0.01161 -0.0361 0.1491 0.053638 0.0355 1.42e-04
+20 0.000747 -0.00560 0.00205 -0.0066 0.2935 -0.190351 0.0368 4.77e-06
+25 0.166736 -0.11049 0.09341 -0.2095 -0.1627 -0.070753 0.0346 4.66e-03
+7 -0.397145 0.15747 0.15662 -0.5912 0.5842 -0.838287 0.0327 3.51e-02
+5 -0.279996 -0.09119 -0.35616 0.7358 0.3694 -0.326563 0.0281 5.18e-02
+6 0.003191 -0.00168 -0.01550 0.0259 0.3447 0.201072 0.0400 7.34e-05
+3 -0.302851 -0.07889 0.25932 -0.4196 -0.2766 0.536017 0.0386 1.85e-02
+4 -0.058151 -0.02719 0.00870 -0.0664 0.8199 0.131003 0.0486 4.83e-04
+2 0.020380 0.00416 -0.00373 0.0206 -0.4158 -0.160401 0.0395 4.62e-05
+9 0.074217 0.09551 -0.10857 0.1271 -0.3481 0.644579 0.0383 1.75e-03
+12 -0.097825 -0.20830 0.04637 0.2864 -0.6601 0.270324 0.0280 8.19e-03
+10 0.149178 0.66594 -0.12975 0.8935 -0.2510 0.000571 0.0118 5.90e-02
+11 0.014687 0.00691 0.00105 0.0191 0.1838 -0.301086 0.0377 4.00e-05
+21 0.148213 0.15461 -0.02915 -0.2531 -0.1892 -0.318491 0.0361 6.81e-03
+ CCA2.1 CCA3.1 hat inf
+18 0.003207 1.39e-03 0.0321
+15 0.000915 4.98e-06 0.0295
+24 0.003071 6.82e-03 0.1135
+27 0.003062 5.41e-05 0.1375
+23 0.009151 5.22e-03 0.0555
+19 0.005873 1.90e-08 0.0176
+22 0.015648 3.07e-03 0.1077
+16 0.005348 3.99e-03 0.0594
+28 0.000559 3.91e-02 0.2256
+13 0.017773 8.11e-03 0.7168 *
+14 0.002424 3.13e-04 0.0158
+20 0.009421 3.96e-03 0.0393
+25 0.002810 5.31e-04 0.0670
+7 0.034274 7.06e-02 0.1635
+5 0.013057 1.02e-02 0.1588
+6 0.012996 4.42e-03 0.1169
+3 0.008056 3.02e-02 0.1833
+4 0.073491 1.88e-03 0.2741
+2 0.018909 2.81e-03 0.1063
+9 0.013158 4.51e-02 0.0978
+12 0.043487 7.29e-03 0.0409
+10 0.004658 2.41e-08 0.0768
+11 0.003694 9.91e-03 0.0632
+21 0.003807 1.08e-02 0.1009
+> plot(mod, type = "n")
+> points(mod, cex = 10*hatvalues(lmod), pch=16, xpd = TRUE)
+> text(mod, display = "bp", col = "blue")
+>
+>
+>
+> cleanEx()
+> nameEx("beals")
+> ### * beals
+>
+> flush(stderr()); flush(stdout())
+>
+> ### Name: beals
+> ### Title: Beals Smoothing and Degree of Absence
+> ### Aliases: beals swan
+> ### Keywords: manip smooth
+>
+> ### ** Examples
+>
+> data(dune)
+> ## Default
+> x <- beals(dune)
+> ## Remove target species
+> x <- beals(dune, include = FALSE)
+> ## Smoothed values against presence or absence of species
+> pa <- decostand(dune, "pa")
+> boxplot(as.vector(x) ~ unlist(pa), xlab="Presence", ylab="Beals")
+> ## Remove the bias of tarbet species: Yields lower values.
+> beals(dune, type =3, include = FALSE)
+ Belper Empnig Junbuf Junart Airpra Elepal Rumace
+2 0.4917791 0.01829337 0.2070478 0.13017869 0.03361524 0.08656209 0.3031318
+13 0.2939740 0.01840390 0.3339728 0.28204736 0.02191119 0.22437598 0.3107471
+4 0.4398775 0.02527165 0.2318440 0.21320122 0.02876960 0.16112450 0.2610006
+16 0.1516354 0.00742115 0.1776781 0.54304667 0.00742115 0.51523810 0.1508409
+6 0.3696613 0.03093489 0.2166255 0.11631772 0.05967771 0.08784329 0.4330705
+1 0.4935662 0.00000000 0.1987270 0.14794933 0.01157407 0.05592045 0.3160963
+8 0.2720122 0.02102222 0.2323566 0.33219882 0.02588333 0.29777644 0.2331043
+5 0.4322142 0.03029752 0.2164230 0.10406655 0.06632771 0.05419613 0.3979230
+17 0.3129358 0.21968254 0.1038156 0.03333333 0.29233391 0.02777778 0.2913631
+15 0.1420251 0.02215971 0.1370235 0.51561767 0.02538032 0.53973883 0.1459458
+10 0.4705094 0.03425596 0.1678040 0.09080902 0.07349127 0.07787078 0.3226025
+11 0.3713794 0.07310076 0.1656396 0.13478872 0.10259239 0.11146095 0.2753878
+9 0.3292320 0.01838883 0.2675624 0.23134366 0.02213662 0.19018562 0.3464870
+18 0.3676424 0.06763669 0.1376386 0.18956922 0.08535723 0.17341352 0.2668100
+3 0.4033301 0.01074444 0.2544828 0.22291082 0.01520046 0.16099072 0.2909068
+20 0.1450721 0.05905666 0.1192488 0.53056145 0.06144615 0.58476297 0.1098600
+14 0.1682755 0.01989756 0.1206518 0.46685537 0.02298398 0.55982776 0.1241545
+19 0.2237843 0.26011417 0.1127832 0.13543701 0.35137675 0.12592768 0.1831168
+12 0.2802946 0.03413656 0.3231724 0.27306206 0.03625297 0.19518888 0.3507140
+7 0.4051777 0.02304603 0.2110045 0.11293676 0.05329633 0.06694311 0.4113622
+ Viclat Brarut Ranfla Cirarv Hyprad Leoaut Potpal
+2 0.18494940 0.7415172 0.13042865 0.08570299 0.07454127 0.9349640 0.026548839
+13 0.06604026 0.7465509 0.33309468 0.07905124 0.04250245 0.9255312 0.056084315
+4 0.12550967 0.7919786 0.21909541 0.12920228 0.06985986 0.9204237 0.039859621
+16 0.03353260 0.8029953 0.66893881 0.05184498 0.01731602 0.8131968 0.201352298
+6 0.18035860 0.8387650 0.10909966 0.03223112 0.11463153 0.9590466 0.030347896
+1 0.17244420 0.7476589 0.08105273 0.09504980 0.07702746 0.8898317 0.000000000
+8 0.10213052 0.8109194 0.40134447 0.06688407 0.06550319 0.8755515 0.100703044
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/vegan -r 1467
More information about the Vegan-commits
mailing list