diff --git a/01_Mesobromion.R b/01_Mesobromion.R
index ab81ea02ff8e23cc21945aaff3df7ab8956485f2..c83884714c4ed71178cff3e161833cf1dce85dfb 100644
--- a/01_Mesobromion.R
+++ b/01_Mesobromion.R
@@ -603,22 +603,28 @@ CWM.wide <- species %>%
#### 4.1 PCA of Y (Bealls) matrix ####
library(vegan)
W.beals <- as.data.frame(beals(species, include=T, type=2))
-#save(W.beals, file="_derived/Mesobromion/MatrixY_Beals.RData")
write.table(W.beals, sep="\t", file="_derived/Mesobromion/MatrixY_Beals.csv")
pca.out <- rda(W.beals)
varexpl <- round((pca.out$CA$eig)/sum(pca.out$CA$eig)*100,1)
cwms.envfit <- envfit(pca.out, CWM.wide, na.rm = T, choices = 1:3)
-
+env.envfit <- envfit(pca.out, env %>%
+ dplyr::select(Temp, Prec, pH=PHIPHOX, C.org=ORCDRC), choices=1:3)
PCA1_2 <- ggplot() +
geom_point(data=as.data.frame(pca.out$CA$u[,1:2]),
aes(x=PC1, y=PC2), pch="+", size=2, alpha=0.8) +
geom_segment(data=as.data.frame(cwms.envfit$vectors$arrows*.2),
aes(x=0, xend=PC1, y=0, yend=PC2), col="Dark blue", arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
- geom_label(data=as.data.frame(cwms.envfit$vectors$arrows*.23) %>%
- rownames_to_column("Trait"),
- aes(x=PC1, y=PC2, label=Trait), col="Dark blue", size=2,
+ geom_segment(data=as.data.frame(env.envfit$vectors$arrows*.2),
+ aes(x=0, xend=PC1, y=0, yend=PC2), col="Dark red", arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
+ geom_label(data=as.data.frame(env.envfit$vectors$arrows*.23) %>%
+ rownames_to_column("Env"),
+ aes(x=PC1, y=PC2, label=Env), col="Dark red", size=2,
position = position_dodge(1) ) +
+ geom_label(data=as.data.frame(cwms.envfit$vectors$arrows*.23) %>%
+ rownames_to_column("Trait"),
+ aes(x=PC1, y=PC2, label=Trait), col="Dark blue", size=2,
+ position = position_dodge(1) ) +
xlab(paste("PC1 (", varexpl[1], "%)", sep="")) +
ylab(paste("PC2 (", varexpl[2], "%)", sep="")) +
theme_bw() +
@@ -631,6 +637,12 @@ PCA1_2 <- ggplot() +
aes(x=PC1, y=PC3), pch="+", size=2, alpha=0.8) +
geom_segment(data=as.data.frame(cwms.envfit$vectors$arrows*.2),
aes(x=0, xend=PC1, y=0, yend=PC3), col="Dark blue", arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
+ geom_segment(data=as.data.frame(env.envfit$vectors$arrows*.2),
+ aes(x=0, xend=PC1, y=0, yend=PC3), col="Dark red", arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
+ geom_label(data=as.data.frame(env.envfit$vectors$arrows*.23) %>%
+ rownames_to_column("Env"),
+ aes(x=PC1, y=PC3, label=Env), col="Dark red", size=2,
+ position = position_dodge(1)) +
geom_label(data=as.data.frame(cwms.envfit$vectors$arrows*.23) %>%
rownames_to_column("Trait"),
aes(x=PC1, y=PC3, label=Trait), col="Dark blue", size=2,
@@ -673,29 +685,55 @@ ggplot() +
theme(panel.grid = element_blank()) +
title("PCA of CWMs")
-#### 4.3 PCA of individual traits #####
-trait.pca <- vegan::rda(as.matrix(traits %>%
- filter(!is.na(rowSums(.)))), na.action=na.omit)
-varexpl.t <- round(trait.pca$CA$eig/sum(trait.pca$CA$eig)*100,1)
-
-ggplot() +
- geom_point(data=as.data.frame(trait.pca$CA$u[,1:2]),
- aes(x=PC1, y=PC2), pch="+", size=2, alpha=0.8) +
- geom_segment(data=as.data.frame(trait.pca$CA$v*2) %>%
- rownames_to_column("Trait") %>%
- mutate(in.try=ifelse(Trait %in% colnames(alltry), "Dark red", "Dark blue")),
- aes(x=0, xend=PC1, y=0, yend=PC2, col=in.try), arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
- geom_text(data=as.data.frame(trait.pca$CA$v*2.1) %>%
- rownames_to_column("Trait") %>%
- mutate(in.try=ifelse(Trait %in% colnames(alltry), "Dark red", "Dark blue")),
- aes(x=PC1, y=PC2, label=Trait, col=in.try), size=3 ) +
+#### 4.3 PCoA of individual traits #####
+# assign weights to traits [Growth forms and root structures are coded as binary, but should be downweighted]
+myweights <- c(rep(1/9, 9), rep(1/17, 17), rep(1, 26))
+traits.gowdis <- FD::gowdis(traits, w=myweights, ord="podani")
+pcoa.species <- wcmdscale(traits.gowdis, eig=T)
+varexpl <- round((pcoa.species$eig)/sum(pcoa.species$eig)*100,1)
+traits.envfit <- envfit(pcoa.species, traits %>%
+ dplyr::select(any_of(traits.sign.alone)),
+ na.rm = T, choices = 1:3)
+
+set.seed(14)
+PCoA.t1 <- ggplot() +
+ geom_point(data=as.data.frame(pcoa.species$points[,1:2]),
+ aes(x=Dim1, y=Dim2), pch="+", size=2, alpha=0.8) +
+ geom_segment(data=as.data.frame(traits.envfit$vectors$arrows*0.2),
+ aes(x=0, xend=Dim1, y=0, yend=Dim2, col="Dark green"),
+ arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
+ geom_label(data=as.data.frame(jitter(traits.envfit$vectors$arrows*0.22, factor = 300)) %>%
+ rownames_to_column("Trait"),
+ aes(x=Dim1, y=Dim2, label=Trait, col="Dark green"), size=2,
+ position = position_dodge(2)) +
scale_color_identity() +
- xlab(paste("PC1 (", varexpl[1], "%)", sep="")) +
- ylab(paste("PC2 (", varexpl[2], "%)", sep="")) +
+ coord_equal(xlim = c(-.25, .25), ylim=c(-.25, .25)) +
+ xlab(paste("PCo1 (", varexpl[1], "%)", sep="")) +
+ ylab(paste("PCo2 (", varexpl[2], "%)", sep="")) +
theme_bw() +
- theme(panel.grid = element_blank()) +
- title("PCA of species-level traits")
+ theme(panel.grid = element_blank())# +
+ #ggtitle("PCoA of species-level traits")
+
+PCoA.t2 <- ggplot() +
+ geom_point(data=as.data.frame(pcoa.species$points[,1:3]),
+ aes(x=Dim1, y=Dim3), pch="+", size=2, alpha=0.8) +
+ geom_segment(data=as.data.frame(traits.envfit$vectors$arrows*0.2),
+ aes(x=0, xend=Dim1, y=0, yend=Dim3, col="Dark green"),
+ arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
+ geom_label(data=as.data.frame(traits.envfit$vectors$arrows*0.22) %>%
+ rownames_to_column("Trait"),
+ aes(x=Dim1, y=Dim3, label=Trait, col="Dark green"), size=2,
+ position = position_dodge(2)) +
+ scale_color_identity() +
+ coord_equal(xlim = c(-.25, .25), ylim=c(-.25, .25)) +
+ xlab(paste("PCo1 (", varexpl[1], "%)", sep="")) +
+ ylab(paste("PCo3 (", varexpl[3], "%)", sep="")) +
+ theme_bw() +
+ theme(panel.grid = element_blank())
+
+PC_traits <- cowplot::plot_grid(PCoA.t1, PCoA.t2, nrow=1)
+ggsave("_pics/PCoA_Traits.png", width=10, height=5, dpi=300, PC_traits)
#### 4.4 PCoA of X (Fuzzy weighted) matrix ####
@@ -770,69 +808,4 @@ ggsave(filename="_pics/PlotDistribution.png", height=8, width=6, dpi=300, basema
-### renaming variables
-"LF_Macroph = LEB_F_Makrophanerophyt
-LF_Nanoh = LEB_F_Nanophanerophyt
-LF_Hemicr = LEB_F_Hemikryptophyt
-LF_Geoph = LEB_F_Geophyt
-LF_Hemiph = LEB_F_Hemiphanerophyt
-LF_Theropf = LEB_F_Therophyt
-LF_Hydroph = LEB_F_Hydrophyt
-LF_Pseudoph = LEB_F_Pseudophanerophyt
-LF_Chamaeph = LEB_F_Chamaephyt
-LL_Pluri_poll = LEB_D_plurienn_pollakanth
-LL_Pluri_hapa = LEB_D_plurienn_hapaxanth
-LL_Annual = LEB_D_annuell
-LL_Bienn = LEB_D_bienn
-VS_Absent = V_V_VER_absent
-VS_RootSprout = V_VER_Wurzelspross
-VS_Runners = V_VER_Ausläufer
-VS_Rhizon = V_VER_Rhizom
-VS_Tuber = V_VER_Innovationsknopse.mit.Wurzelknolle
-VS_StorageRoot = V_VER_Innovationsknospe.mit.Speicherwurzel
-VS_RunningTuber = V_VER_Ausläuferknolle
-VS_Brutsprösschen = V_VER_Brutsprösschen ## How to translate?
-VS_Fragment = V_VER_Fragmentation
-VS_Turio = V_VER_Turio
-24 V_VER_Sprossknolle
-25 V_VER_phyllogener_Spross
-26 V_VER_Rhizompleiokorm
-27 V_VER_Zwiebel
-28 V_VER_Ausläuferrhizom
-29 V_VER_Ausläuferzwiebel
-30 V_VER_Bulbille
-31 ROS_T_rosettenlose.Pflanzen
-32 ROS_T_Halbrosettenpflanze
-33 ROS_T_Ganzrosettenpflanzen
-34 BL_AUSD_immergrün
-35 BL_AUSD_sommergrün
-36 BL_AUSD_vorsommergrün
-37 BL_AUSD_überwinternd_grün
-38 BL_ANAT_skleromorph
-39 BL_ANAT_mesomorph
-40 BL_ANAT_hygromorph
-41 BL_ANAT_hydromorph
-42 BL_ANAT_blattsukkulent
-43 BL_ANAT_helomorph
-44 BL_BEG
-45 BL_END
-46 BL_DAU
-47 LeafArea
-48 SLA
-49 LeafC.perdrymass
-50 LeafN
-51 LeafP
-52 PlantHeight
-53 SeedMass
-54 Seed.length
-55 LDMC
-56 LeafNperArea
-57 LeafNPratio
-58 Leaf.delta.15N
-59 Seed.num.rep.unit
-60 Leaffreshmass
-61 Disp.unit.leng
-62 BLU_KL
-63 REPR_T
-64 STRAT_T"