varespecvarechem_Indirect
FV
2023-10-12
packages
# install if required
if(!require(tidyverse)){install.packages("tidyverse")}
if(!require(googlesheets4)){install.packages("googlesheets4")}
if(!require(GGally)){install.packages("GGally")}
if(!require(reshape2)){install.packages("reshape2")}
if(!require(lme4)){install.packages("lme4")}
if(!require(compiler)){install.packages("compiler")}
if(!require(parallel)){install.packages("parallel")}
if(!require(boot)){install.packages("boot")}
if(!require(lattice)){install.packages("lattice")}
if(!require(car)){install.packages("car")}
if(!require(ggpubr)){install.packages("ggpubr")}
if(!require(vegan)){install.packages("vegan")}
if(!require(ggsci)){install.packages("ggsci")}
if(!require(ggrepel)){install.packages("ggrepel")}
if(!require(ggforce)){install.packages("ggforce")}
if(!require(plotly)){install.packages("plotly")}
if(!require(rgl)){install.packages("rgl")}
if(!require(mgcv)){install.packages("mgcv")}
if(!require(factoextra)){install.packages("factoextra")}
if(!require(devtools)){install.packages("devtools")}
if(!require(gridExtra)){install.packages("gridExtra")}
if(!require(BiodiversityR)){install.packages("BiodiversityR")}
# Load
library(tidyverse)
library(googlesheets4); gs4_deauth()
library(googledrive)
library(GGally)
library(reshape2)
library(lme4)
library(compiler)
library(parallel)
library(boot)
library(lattice)
library(car)
library(ggpubr)
library(vegan)
library(readxl)
library(ggsci)
library(ggrepel)
library(ggforce)
library(plotly)
library(rgl)
library(mgcv)
library(dplyr)
library(factoextra)
library(devtools)
library(gridExtra)
library(BiodiversityR) # also loads vegan
library(readxl)
library(ggsci)
library(ggrepel)
library(ggforce)
require(stats)
install_github("vqv/ggbiplot")
require(ggbiplot)Data
[Varespec-Varechem Datasets] ](https://www.researchgate.net/publication/227830523_Effects_of_reindeer_grazing_on_vegetation_in_dry_Pinus_sylvestris_forests)
data(varespec)
data(varechem)
M1 <- varespec
M2 <- varechemvegan
PCA [Vegan PCKG]
ord.model1 <- rda(M1, center=TRUE)
summary(ord.model1)##
## Call:
## rda(X = M1, center = TRUE)
##
## Partitioning of variance:
## Inertia Proportion
## Total 1826 1
## Unconstrained 1826 1
##
## Eigenvalues, and their contribution to the variance
##
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6
## Eigenvalue 982.9788 464.3040 132.25052 73.9337 48.41829 37.00937
## Proportion Explained 0.5384 0.2543 0.07244 0.0405 0.02652 0.02027
## Cumulative Proportion 0.5384 0.7927 0.86519 0.9057 0.93220 0.95247
## PC7 PC8 PC9 PC10 PC11 PC12
## Eigenvalue 25.72624 19.70557 12.274191 10.435361 9.350783 2.798400
## Proportion Explained 0.01409 0.01079 0.006723 0.005716 0.005122 0.001533
## Cumulative Proportion 0.96657 0.97736 0.984083 0.989799 0.994921 0.996454
## PC13 PC14 PC15 PC16 PC17
## Eigenvalue 2.327555 1.3917180 1.2057303 0.8147513 0.3312842
## Proportion Explained 0.001275 0.0007623 0.0006604 0.0004463 0.0001815
## Cumulative Proportion 0.997729 0.9984910 0.9991515 0.9995977 0.9997792
## PC18 PC19 PC20 PC21 PC22
## Eigenvalue 0.1866564 1.065e-01 6.362e-02 2.521e-02 1.652e-02
## Proportion Explained 0.0001022 5.835e-05 3.485e-05 1.381e-05 9.048e-06
## Cumulative Proportion 0.9998814 9.999e-01 1.000e+00 1.000e+00 1.000e+00
## PC23
## Eigenvalue 4.590e-03
## Proportion Explained 2.514e-06
## Cumulative Proportion 1.000e+00
##
## Scaling 2 for species and site scores
## * Species are scaled proportional to eigenvalues
## * Sites are unscaled: weighted dispersion equal on all dimensions
## * General scaling constant of scores: 14.31485
##
##
## Species scores
##
## PC1 PC2 PC3 PC4 PC5 PC6
## Callvulg -1.470e-01 0.3483315 -2.506e-01 0.6029965 0.3961208 -4.608e-01
## Empenigr 1.645e-01 -0.3731965 -5.543e-01 -0.7565217 -0.7432702 -3.484e-01
## Rhodtome -6.787e-02 -0.0810708 -6.464e-02 -0.0963929 -0.1749038 -4.008e-02
## Vaccmyrt -5.429e-01 -0.7446016 1.973e-01 -0.2791788 -0.5688937 -3.689e-01
## Vaccviti 9.013e-02 -0.7247759 -6.555e-01 -1.3366155 -0.7107105 -2.615e-01
## Pinusylv 3.404e-02 -0.0132527 4.465e-04 -0.0079399 -0.0012504 3.800e-03
## Descflex -7.359e-02 -0.0961364 6.485e-02 -0.0158203 -0.0741307 -5.328e-02
## Betupube -8.789e-04 -0.0008548 -6.637e-03 -0.0077726 -0.0072147 4.759e-04
## Vacculig -6.959e-02 0.2123669 6.410e-02 0.0136389 -0.0712061 -1.179e-01
## Diphcomp 3.807e-03 0.0333497 -1.026e-02 -0.0291455 -0.0086379 -2.263e-02
## Dicrsp -4.663e-01 -0.3675288 -1.143e-01 -0.2490898 0.7694280 1.225e+00
## Dicrfusc -1.176e+00 -0.4484656 -1.338e+00 2.1943059 -0.9187333 4.027e-02
## Dicrpoly -1.597e-02 -0.0326269 -6.610e-02 -0.1166483 -0.0369723 6.576e-02
## Hylosple -2.935e-01 -0.4217426 4.523e-01 -0.0656152 -0.0358715 -2.234e-01
## Pleuschr -3.890e+00 -4.3657447 2.343e+00 0.3747564 0.1476496 -3.060e-01
## Polypili -5.032e-04 0.0069583 7.660e-03 -0.0005528 -0.0034856 5.849e-03
## Polyjuni -9.888e-02 -0.0584117 -7.234e-02 -0.0518007 0.0438920 1.472e-01
## Polycomm -4.028e-03 -0.0057429 -5.248e-03 -0.0104641 -0.0058180 9.059e-05
## Pohlnuta 9.473e-03 -0.0116498 -1.010e-02 -0.0112498 -0.0009520 7.691e-03
## Ptilcili -2.749e-02 -0.0257213 -2.522e-01 -0.3745203 -0.2917694 -4.907e-03
## Barbhatc -6.264e-03 -0.0070227 -7.152e-02 -0.1029116 -0.0825601 -1.183e-03
## Cladarbu -7.281e-01 3.0180666 2.016e-01 0.1400429 0.8708266 -1.253e+00
## Cladrang 8.117e-01 4.2257879 2.338e+00 0.2831533 -0.7644829 4.473e-01
## Cladstel 9.580e+00 -2.0074825 6.078e-01 0.4155714 0.1401899 -2.112e-01
## Cladunci -3.796e-01 0.0706264 -6.867e-01 0.1413783 0.9384995 1.495e-01
## Cladcocc 5.652e-03 0.0096766 -3.946e-03 0.0117038 0.0033077 -4.020e-03
## Cladcorn -1.367e-02 -0.0047214 -1.714e-02 -0.0270900 0.0004555 1.197e-03
## Cladgrac -1.037e-02 0.0083562 -6.489e-03 -0.0132067 0.0018124 1.012e-02
## Cladfimb 3.163e-03 0.0028937 -1.436e-02 0.0004203 -0.0043761 -1.344e-02
## Cladcris -2.649e-02 0.0010675 -4.626e-02 -0.0128097 0.0150875 -3.677e-02
## Cladchlo 1.347e-02 -0.0054141 -7.871e-03 -0.0099859 -0.0028496 -4.324e-04
## Cladbotr -2.051e-03 -0.0009378 -5.815e-03 -0.0095913 -0.0057131 -1.141e-03
## Cladamau 3.733e-05 0.0020956 3.925e-04 -0.0009025 -0.0009981 2.615e-04
## Cladsp 6.254e-03 -0.0021214 -1.997e-03 0.0011793 0.0023250 -3.246e-03
## Cetreric -4.938e-03 0.0079861 -1.740e-02 0.0046582 0.0488265 1.964e-02
## Cetrisla 1.715e-02 -0.0103936 -1.530e-02 -0.0214109 -0.0140577 3.490e-03
## Flavniva 8.793e-02 0.0932642 -3.761e-03 0.0536027 0.1450323 1.915e-03
## Nepharct -5.549e-02 -0.0368786 -3.638e-02 0.0055964 0.0415127 1.027e-01
## Stersp -5.531e-02 0.3478548 2.425e-01 0.0088566 -0.1772134 2.711e-01
## Peltapht -3.661e-03 -0.0014739 1.729e-03 -0.0067221 -0.0043042 -5.434e-05
## Icmaeric -1.435e-03 0.0040087 1.296e-03 0.0023703 -0.0024341 2.959e-03
## Cladcerv 8.746e-04 0.0001554 2.667e-05 0.0004662 0.0006877 7.336e-04
## Claddefo -5.274e-02 0.0023107 -8.009e-02 -0.0278485 0.0198955 -1.490e-02
## Cladphyl 1.581e-02 -0.0057716 5.064e-04 0.0010606 0.0027075 -1.441e-03
##
##
## Site scores (weighted sums of species scores)
##
## PC1 PC2 PC3 PC4 PC5 PC6
## 18 -1.02674 2.59169 -1.53869 -3.23113 -1.104731 -2.20341
## 15 -2.64700 -1.62646 1.01889 1.88048 1.152767 -1.41888
## 24 -2.44595 -2.01412 0.12182 -2.44215 5.875564 8.03640
## 27 -3.02575 -4.32492 4.44163 -1.54364 -2.466762 -2.57141
## 23 -1.86899 -0.35380 -1.98920 -4.11912 -2.428561 -0.85734
## 19 -0.02298 -1.59265 -0.43949 -1.43686 1.058122 -1.18048
## 22 -2.53975 -1.70578 -4.10001 7.63091 -5.062274 -0.64415
## 16 -2.08714 0.06163 -2.32296 5.84188 -3.568851 1.63566
## 28 -3.77083 -5.80241 6.57611 0.02751 0.908046 -3.10610
## 13 -0.38714 2.82852 -0.32352 2.12507 3.631355 -3.88730
## 14 -1.75570 0.75366 -5.69429 1.68697 6.098213 -1.28171
## 20 -1.65653 0.32765 -1.93100 -2.50536 0.065456 1.72887
## 25 -2.39601 -1.83554 -1.65364 0.24230 1.756334 4.36881
## 7 -1.08594 5.78134 1.89095 -0.90494 -0.006891 -2.44158
## 5 -0.80203 6.28624 4.85360 0.70195 -3.707269 5.90752
## 6 -0.19764 5.11581 1.24710 -0.70641 1.876853 -4.68229
## 3 3.79496 1.11636 1.93079 1.75008 -0.164043 1.75417
## 4 1.24779 1.77832 -0.17925 1.13363 2.999648 0.07212
## 2 5.49157 -0.67853 1.58337 1.13337 -1.472922 1.97780
## 9 6.02767 -3.10958 -1.09190 -0.28517 0.894284 -1.20958
## 12 4.19915 -1.40814 -0.06645 -0.79569 -0.790560 -0.37424
## 10 6.18267 -2.32211 -0.76545 0.43683 0.430813 -0.78609
## 11 1.09782 0.55023 3.41127 0.38027 -0.325506 1.07113
## 21 -0.32554 -0.41738 -4.97967 -7.00076 -5.649086 0.09209.long formats
plot2 <- ordiplot(ord.model1, choices = c(1,2), display = "all" , scaling=2)
sites.long <- sites.long(plot2, env.data=M2)
species.long <- species.long(plot2)
axis.long <- axis.long(ord.model1, choices = c(1,2))env.fit
spec.envfit <- envfit(plot2, env=M1)
spec.data.envfit <- data.frame(r=spec.envfit$vectors$r, p=spec.envfit$vectors$pvals)
species.long2 <- species.long(plot2, spec.data=spec.data.envfit)
species.long2## r p axis1 axis2 labels
## Callvulg 0.051187552 0.558 -1.469634e-01 0.3483315194 Callvulg
## Empenigr 0.065059383 0.503 1.645106e-01 -0.3731965262 Empenigr
## Rhodtome 0.106778779 0.296 -6.787182e-02 -0.0810707906 Rhodtome
## Vaccmyrt 0.325836407 0.024 -5.428697e-01 -0.7446016365 Vaccmyrt
## Vaccviti 0.122428205 0.255 9.013033e-02 -0.7247759128 Vaccviti
## Pinusylv 0.158281736 0.144 3.403747e-02 -0.0132526677 Pinusylv
## Descflex 0.223857621 0.076 -7.358820e-02 -0.0961363923 Descflex
## Betupube 0.005151419 0.942 -8.789003e-04 -0.0008548446 Betupube
## Vacculig 0.153297953 0.217 -6.959288e-02 0.2123669228 Vacculig
## Diphcomp 0.054347092 0.524 3.807206e-03 0.0333497081 Diphcomp
## Dicrsp 0.106350827 0.313 -4.663296e-01 -0.3675287970 Dicrsp
## Dicrfusc 0.170332891 0.118 -1.176342e+00 -0.4484656402 Dicrfusc
## Dicrpoly 0.025553708 0.776 -1.597212e-02 -0.0326268996 Dicrpoly
## Hylosple 0.408173512 0.006 -2.935457e-01 -0.4217426333 Hylosple
## Pleuschr 0.853396295 0.001 -3.889920e+00 -4.3657447357 Pleuschr
## Polypili 0.114907841 0.217 -5.031813e-04 0.0069583267 Polypili
## Polyjuni 0.055502447 0.499 -9.887841e-02 -0.0584117017 Polyjuni
## Polycomm 0.093283226 0.333 -4.028316e-03 -0.0057429098 Polycomm
## Pohlnuta 0.194779169 0.086 9.472722e-03 -0.0116498455 Pohlnuta
## Ptilcili 0.003032132 0.951 -2.749459e-02 -0.0257213287 Ptilcili
## Barbhatc 0.002113016 0.963 -6.264138e-03 -0.0070226891 Barbhatc
## Cladarbu 0.763853471 0.001 -7.281023e-01 3.0180666192 Cladarbu
## Cladrang 0.742032353 0.001 8.117243e-01 4.2257878862 Cladrang
## Cladstel 0.993190176 0.001 9.580302e+00 -2.0074825204 Cladstel
## Cladunci 0.053808329 0.559 -3.795841e-01 0.0706264286 Cladunci
## Cladcocc 0.141894160 0.205 5.652253e-03 0.0096766300 Cladcocc
## Cladcorn 0.023541316 0.798 -1.366727e-02 -0.0047213924 Cladcorn
## Cladgrac 0.126072374 0.244 -1.037304e-02 0.0083561740 Cladgrac
## Cladfimb 0.025276835 0.779 3.163257e-03 0.0028936771 Cladfimb
## Cladcris 0.034362884 0.706 -2.649356e-02 0.0010674739 Cladcris
## Cladchlo 0.276083753 0.029 1.346539e-02 -0.0054141054 Cladchlo
## Cladbotr 0.016629369 0.814 -2.051457e-03 -0.0009377894 Cladbotr
## Cladamau 0.121916274 0.260 3.732507e-05 0.0020955658 Cladamau
## Cladsp 0.144320284 0.223 6.254496e-03 -0.0021214420 Cladsp
## Cetreric 0.017965951 0.864 -4.938379e-03 0.0079861459 Cetreric
## Cetrisla 0.149063015 0.179 1.715353e-02 -0.0103935843 Cetrisla
## Flavniva 0.035275985 0.630 8.792701e-02 0.0932642250 Flavniva
## Nepharct 0.040198652 0.606 -5.548726e-02 -0.0368785816 Nepharct
## Stersp 0.255585428 0.034 -5.531083e-02 0.3478547644 Stersp
## Peltapht 0.020204820 0.860 -3.660595e-03 -0.0014739143 Peltapht
## Icmaeric 0.261912394 0.047 -1.434714e-03 0.0040087454 Icmaeric
## Cladcerv 0.047790334 0.626 8.745883e-04 0.0001554259 Cladcerv
## Claddefo 0.095649612 0.364 -5.273539e-02 0.0023107405 Claddefo
## Cladphyl 0.355294514 0.012 1.580780e-02 -0.0057715835 Cladphyl## Vector Fit
vectors.envfit <- envfit(plot2, env= M2)
vectors.long3 <- vectorfit.long(vectors.envfit)
vectors.long3## vector axis1 axis2 r p
## N N -0.70851888 0.70569186 0.24520397 0.042
## P P -0.08868186 -0.99606000 0.38987120 0.006
## K K -0.31858192 -0.94789533 0.29058690 0.027
## Ca Ca -0.33361020 -0.94271111 0.35689387 0.011
## Mg Mg -0.42512286 -0.90513565 0.27366376 0.045
## S S 0.25137004 -0.96789106 0.19065688 0.124
## Al Al 0.85495235 0.51870655 0.42540274 0.007
## Fe Fe 0.79122906 0.61151989 0.36992184 0.009
## Mn Mn -0.65591734 -0.75483273 0.51521299 0.002
## Zn Zn -0.35219925 -0.93592505 0.21313172 0.069
## Mo Mo 0.40478509 0.91441185 0.02742188 0.752
## Baresoil Baresoil -0.99956100 -0.02962782 0.28058127 0.035
## Humdepth Humdepth -0.68995633 -0.72385100 0.46214509 0.003
## pH pH 0.89119076 0.45362874 0.24269465 0.066### Spec >= 40%
species.long3 <- species.long2[species.long2$r >= 0.4, ]
species.long3## r p axis1 axis2 labels
## Hylosple 0.4081735 0.006 -0.2935457 -0.4217426 Hylosple
## Pleuschr 0.8533963 0.001 -3.8899202 -4.3657447 Pleuschr
## Cladarbu 0.7638535 0.001 -0.7281023 3.0180666 Cladarbu
## Cladrang 0.7420324 0.001 0.8117243 4.2257879 Cladrang
## Cladstel 0.9931902 0.001 9.5803025 -2.0074825 Cladstel### Var Ambientales r >= 0.35
vectors.long3 <- vectors.long3[vectors.long3$r >= 0.35, ]
vectors.long3## vector axis1 axis2 r p
## P P -0.08868186 -0.9960600 0.3898712 0.006
## Ca Ca -0.33361020 -0.9427111 0.3568939 0.011
## Al Al 0.85495235 0.5187065 0.4254027 0.007
## Fe Fe 0.79122906 0.6115199 0.3699218 0.009
## Mn Mn -0.65591734 -0.7548327 0.5152130 0.002
## Humdepth Humdepth -0.68995633 -0.7238510 0.4621451 0.003plots [r >= 35% ]
f = 5
plotgg2 <- ggplot() +
geom_vline(xintercept = c(0), color = "grey70", linetype = 2) +
geom_hline(yintercept = c(0), color = "grey70", linetype = 2) +
xlab(axis.long[1, "label"]) +
ylab(axis.long[2, "label"]) +
scale_x_continuous(sec.axis = dup_axis(labels=NULL, name=NULL)) +
scale_y_continuous(sec.axis = dup_axis(labels=NULL, name=NULL)) +
geom_point(data=sites.long,
aes(x=axis1, y=axis2),
size=5, shape=1) +
geom_point(data=species.long2,aes(x=axis1*f, y=axis2*f)) +
geom_segment(data=species.long3,aes(x=0, y=0, xend=axis1*f, yend=axis2*f),colour="red", linewidth=0.7, arrow=arrow(angle=20)) +
geom_text_repel(data=species.long3,aes(x=axis1*f, y=axis2*f, label=labels),colour="red") +
#geom_segment(data=vectors.long3,aes(x=0, y=0, xend=axis1*f, yend=axis2*f),colour="blue", linewidth=0.7, arrow=arrow(angle=20)) +
#geom_text_repel(data=vectors.long3,aes(x=axis1*f, y=axis2*f, label=vector),colour="blue") +
ggsci::scale_colour_npg() +
coord_fixed(ratio=1)
plotgg2
Ordisurf
plot2 <- ordiplot(ord.model1, choices = c(1,2), display = "all" , scaling=2)
ordisurf(plot2, y=as.numeric(M2$P), add=TRUE, col="blue")