User Tools

Site Tools


Sidebar

Multivariate Analysis of Vegetation Data in R


Link to website with data and examples:
http://bit.ly/AnaDatR


Link to this website:
http://bit.ly/AnaVegR


Contact:
David Zelený ([email protected])
Ching-Feng Li ([email protected])

iavs2015:scripts

Script from the course

All scripts and data used in the class were packed into one folder anavegr and can be downloaded here: anavegr.zip

first practice.r
a <- sqrt (25)
b <- a*a
?sqrt
 
example (sqrt)
 
a <- 5*5
d <- "Dobry den"
class (d)
class (a)
a + d
 
e <- "25"
class (e)
a + e
 
 
A <- seq (from = 1, to = 30)
length (A)
A + a
class (A)
B <- c("Dobry den", "jak se mate", "date si", "pivo")
B
length (B)
C <- c (A, B)
C
 
M <- matrix (data = A, nrow = 6, ncol = 5, byrow = TRUE)
M
 
M [4,4]
M [2, 5]
M [2, ]
M [, 5]
dim (M)
unconstrained ordination.r
# Example of unconstrained ordination analysis using Vltava data
 
# Import of data via clipboard:
# vltava.spe <- read.delim (file = "clipboard", row.names = 1)
 
# Import of data from the file in the computer:
#vltava.spe <- read.delim (file = "C:\\Users\\Zeleny\\Dropbox\\anavegr\\vltava.spe.txt", row.names = 1)
#vltava.env <- read.delim (file = "C:\\Users\\Zeleny\\Dropbox\\anavegr\\vltava.env.txt")
 
# Import of data stored online:
vltava.spe <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/vltava-spe.txt', row.names = 1)
vltava.env <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/vltava-env.txt')
 
vltava.env$GROUP
 
# install.packages ("vegan")
library (vegan)
 
# log transform the species composition data:
#vltava.spe.log <- log (vltava.spe + 1)
vltava.spe.log <- log1p (vltava.spe)
 
# Calculate detrended correspondence analysis:
DCA <- decorana (vltava.spe.log)
DCA
 
# Calculate correspondence analysis (to know the total inertia):
CA <- cca (vltava.spe.log)
CA
 
plot (DCA) # this is drawing DCA ordination diagram
# plot (CA) # this is drawing CA ordination diag.
 
#plot (DCA, choices = c(2,3), display = 'sites')
 
windows ()
plot (DCA, display = 'sites', type = 'none')
# ordispider (ord = DCA, groups = vltava.env$GROUP, label = T)
ordispider (ord = DCA, groups = vltava.env$GROUP, label = T, show.group = 1)
ordispider (ord = DCA, groups = vltava.env$GROUP, label = T, show.group = 2, col = 'red')
ordispider (ord = DCA, groups = vltava.env$GROUP, label = T, show.group = 3, col = 'green')
ordispider (ord = DCA, groups = vltava.env$GROUP, label = T, show.group = 4, col = 'blue')
 
# Project environmental variables as supplementary to ordination diagram:
env <- vltava.env[, c("ASPSSW", "pH", "SOILDPT")]
ef <- envfit (DCA, env)
plot (ef)
 
# project randomly generated environmental variables as supplementary to ordination diagram
rand.var <- matrix (rnorm (6*97), ncol = 6, nrow = 97)
ef.rand <- envfit (DCA, rand.var)
ef.rand
plot (ef.rand, col = 'navy')
constrained ordination.r
# Example of constrained ordination using seedlings data from Ohrazeni
 
# Upload the dataset stored online:
seedlings.spe <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/seedl-spe.txt', row.names = 1)
seedlings.env <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/seedl-env.txt', row.names = 1)
 
 
library (vegan)
DCA <- decorana (seedlings.spe)
DCA
 
#seedlings.spe <- log1p (seedlings.spe)
 
# Calculate redundancy analysis (linear method) with treatment as explanatory variable and bloc as covariable (coded as factor)
RDA <- rda (seedlings.spe ~ treatment + Condition (as.factor (block)), data = seedlings.env)
RDA
 
set.seed (234)  # this sets the random number generator (use this if you want the same results of permutation tests)
anova (RDA)
windows ()
plot (RDA)
set.seed (234)
anova (RDA, by = 'axis')
 
# How to calculate adjusted R2:
r2 <- RsquareAdj (RDA)
r2$adj.r.squared
forward selection.r
# This script calculates forward selection on data from Carpathian wetlands (Hajek et al.)
library (vegan)
 
# Loading data
vasc <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/vasc_plants.txt', row.names = 1)
chem <- read.delim ('http://www.davidzeleny.net/anadat-r/data-download/chemistry.txt', row.names = 1)
 
decorana (vasc)
# Call:
#   decorana(veg = vasc) 
# 
# Detrended correspondence analysis with 26 segments.
# Rescaling of axes with 4 iterations.
# 
# DCA1   DCA2    DCA3    DCA4
# Eigenvalues     0.4172 0.1657 0.11176 0.09538
# Decorana values 0.4244 0.1657 0.09826 0.08197
# Axis lengths    3.0884 2.5454 1.40922 1.64015
 
vasc.hell <- decostand (vasc, method = "hell")
# chem2 <- chem[,c(1,2,3,4,5,6,7,8,9,10,11,12,13,14)]
# chem2 <- chem[,1:14]
chem2 <- chem[,-15]
names (chem2)
 
# Run the global RDA model
rda.all <- rda (vasc.hell ~ ., data = chem2)
rda.all
 
# variation explained by individual axes
total.inertia <- rda.all$tot.chi
eig.constrained <- rda.all$CCA$eig
eig.unconstr <- rda.all$CA$eig
 
# percentage variation explained by axes
expl.constr <- eig.constrained/total.inertia
expl.unconstr <- eig.unconstr/total.inertia
 
# drawing this variation
barplot (c(expl.constr, expl.unconstr))
 
rda.0 <- rda (vasc.hell ~ 1, data = chem2)
rda.0
 
expl.pca <- rda.0$CA$eig/total.inertia
expl.pca
 
# Calculates forward selection:
fw <- ordistep (rda.0, scope = formula (rda.all), direction = 'forward')
fw
 
# Draw the results of RDA with variables selected by forward selection:
windows ()
plot (fw)
anova (fw)
twinspan in R.r
# TWINSPAN in R applied on Danube dataset
 
# Install the twinspanR library from GitHub repository:
install.packages ('devtools')
devtools::install_github("zdealveindy/twinspanR")
 
# initiate the library and data
library (twinspanR)
data (danube)
?danube
 
# calculates standard twispan, with default setting
tw <- twinspan (danube$spe, modif = F)
 
# cuts the results into vector of sample assignments to groups
tw.groups <- cut (tw, cluster = 4)
 
# Calculates and draws the DCA ordination diagram
DCA <- decorana (log1p (danube$spe))
plot (DCA, display = 'sites', type = 'n')
points (DCA, col = tw.groups)
ordihull (DCA, groups = tw.groups)
plot (envfit (DCA, danube$env[,1:3]))
 
# Prints the two-way ordered table from twinspan
print (tw, 'table')
iavs2015/scripts.txt · Last modified: 2017/10/11 20:36 (external edit)