From d02a6bc308c03b15f70fe1c61b3d329c3f0b0310 Mon Sep 17 00:00:00 2001
From: Francesco Sabatini <francesco.sabatini@idiv.de>
Date: Wed, 2 Sep 2020 18:49:10 +0200
Subject: [PATCH] Minor fix to pictures and ordinations
---
02_Mesobromion_ExamineOutput.R | 176 +++++++++++++++++++--------------
1 file changed, 103 insertions(+), 73 deletions(-)
diff --git a/02_Mesobromion_ExamineOutput.R b/02_Mesobromion_ExamineOutput.R
index a2e1fd3..1dfd91f 100644
--- a/02_Mesobromion_ExamineOutput.R
+++ b/02_Mesobromion_ExamineOutput.R
@@ -72,6 +72,9 @@ is.binary <- function(x){(all(na.omit(x) %in% 0:1))}
####1.1 Table S3 - Trait recap #####
trait.recap <- traits %>%
+ ##back transform log-transformed traits from TRY
+ mutate_at(.vars=vars(LA:DispUnitLength),
+ .funs=~exp(.)) %>%
dplyr::select(-ends_with("1")) %>% ## exclude two traits that shouldn't be there V_VER_present1 & V_VER_absent1
mutate_if(.predicate = ~is.numeric(.), .funs=~round(.,3)) %>%
summarize_all(.funs = list(xxxType.of.variable=~ifelse(is.binary(.), "binary",
@@ -82,7 +85,7 @@ trait.recap <- traits %>%
ifelse(is.binary(.),
paste(names(table(.)), "=", table(.), collapse="; "),
ifelse(is.numeric(.),
- paste(range(., na.rm=T), collapse=" - "),
+ paste(round(quantile(.,probs = c(0.05,.95), na.rm=T),1), collapse=" - "),
ifelse(is.factor(.),
paste(levels(.), collapse=", "),
NA)))))) %>%
@@ -96,6 +99,21 @@ trait.recap <- traits %>%
dplyr::select(trait.name, Long_English_name),
by=c("Variable" = "trait.name")) %>%
dplyr::select(Code=Variable, Variable=Long_English_name, everything())
+
+trait.recap <- trait.recap %>%
+ mutate(Variable=ifelse(Code %in% c("FP_Beg", "FP_End", "FP_Dur"), paste0(Variable, " (month)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="LA", paste0(Variable, " (mm²)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="SLA", paste0(Variable, " (m2kg-1)"), Variable)) %>%
+ mutate(Variable=ifelse(Code %in% c("LeafNPratio", "LDMC"), paste0(Variable, " (g g-1)"), Variable)) %>%
+ mutate(Variable=ifelse(Code %in% c("LeafN","LeafP", "LeafCdrymass"), paste0(Variable, " (mg g-1)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="Height", paste0(Variable, " (m)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="SeedMass", paste0(Variable, " (mg)"), Variable)) %>%
+ mutate(Variable=ifelse(Code %in% c("SeedLength","DispUnitLength"), paste0(Variable, " (mm)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="LeafNarea", paste0(Variable, " (g m-2)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="Leaf.delta.15N", paste0(Variable, " (ppm)"), Variable)) %>%
+ mutate(Variable=ifelse(Code=="LeafFreshMass", paste0(Variable, " (g)"), Variable))
+
+
## calculate trait coverage across plots, both on p\a and cover
trait.coverage <- species.cov %>%
@@ -128,29 +146,29 @@ trait.recap <- trait.recap %>%
write_csv(trait.recap, path="_output/S3_Traits.recap.csv")
-
+###### ___________________ ######
###### 2. Import output from Cluster ####
##### 2.1 Cover values ######
# ## traits onebyone
- myfilelist0 <- list.files(path="_derived/Mesobromion/Cover/onebyone/", pattern="HIDDENcov-nona2_[0-9]+_.RData", full.names = T)
- dataFiles0 = purrr::map(myfilelist0, function(x){get(load(x))})
- #load("_derived/Mesobromion/PresAbs/HIDDEN_round_11.RData")
- corXY.alone = bind_rows(dataFiles0) %>%
- as_tibble() %>%
- distinct()
- corXY.alone.ci <- get.ci(corXY.alone)
- corXY.alone.ci <- corXY.alone.ci %>%
- mutate(trait1=Trait.comb) %>%
- mutate_at(.vars=vars(trait1),
- .funs=~factor(.,
- levels=trait.labs$Trait.comb,
- labels=trait.labs$trait.name)) %>%
- arrange(ntraits, desc(Coef.obs)) %>%
- dplyr::select(Trait.comb, Test, n, ntraits, everything()) %>%
- #mutate(Trait.comb=NA) %>%
- mutate(run="alone") %>%
- arrange(desc(sign_plus, Coef.obs))
-
+# myfilelist0 <- list.files(path="_derived/Mesobromion/Cover/onebyone/", pattern="HIDDENcov-nona2_[0-9]+_.RData", full.names = T)
+# dataFiles0 = purrr::map(myfilelist0, function(x){get(load(x))})
+# #load("_derived/Mesobromion/PresAbs/HIDDEN_round_11.RData")
+# corXY.alone = bind_rows(dataFiles0) %>%
+# as_tibble() %>%
+# distinct()
+# corXY.alone.ci <- get.ci(corXY.alone)
+# corXY.alone.ci <- corXY.alone.ci %>%
+# mutate(trait1=Trait.comb) %>%
+# mutate_at(.vars=vars(trait1),
+# .funs=~factor(.,
+# levels=trait.labs$Trait.comb,
+# labels=trait.labs$trait.name)) %>%
+# arrange(ntraits, desc(Coef.obs)) %>%
+# dplyr::select(Trait.comb, Test, n, ntraits, everything()) %>%
+# #mutate(Trait.comb=NA) %>%
+# mutate(run="alone") %>%
+# arrange(desc(sign_plus, Coef.obs))
+#
### sequential trait combo
myfilelist1 <- list.files(path="_derived/Mesobromion/Cover/", pattern="HIDDENcov_round_[0-9]+.RData", full.names = T)
@@ -181,7 +199,7 @@ corXY.ci <- corXY.all.ci # %>%
#### 2.1.2 One by one - Graph of r(XY) ####
### list traits being significant when taken one by one
-traits.sign.alone.cov <- corXY.all.ci %>%
+traits.sign.alone.cov <- corXY.ci %>%
#filter(run=="alone") %>%
filter(ntraits==1) %>%
filter(sign_plus==T) %>%
@@ -194,12 +212,13 @@ mydata.byone <- corXY.ci %>%
#filter(run=="alone") %>%
filter(ntraits==1)%>%
arrange(Coef.obs) %>%
- mutate(seq=1:n())
+ mutate(seq=1:n()) %>%
+ left_join(trait.labs, by=c("trait1"="trait.name"))
write_csv(mydata.byone %>%
dplyr::select(Trait.comb:sign_plus) %>%
arrange(desc(sign_plus), desc(Coef.obs)),
- path = "_output/S992_SolutionsOneByOIne.cov.csv")
+ path = "_output/S992_SolutionsOneByOne.pa.csv")
(top.first <- ggplot(data=mydata.byone %>%
mutate(sign_plus=fct_rev(as.factor(sign_plus)))) +
@@ -207,7 +226,7 @@ write_csv(mydata.byone %>%
col=sign_plus)) +
geom_point(aes(x=Coef.obs, y=seq), pch=15) +
scale_y_continuous(breaks=mydata.byone$seq,
- labels=mydata.byone$trait1, name=NULL) +
+ labels=mydata.byone$Long_English_name , name=NULL) +
scale_x_continuous(name="RD correlation") +
labs(color = "Significant") +
theme_bw() +
@@ -218,7 +237,7 @@ write_csv(mydata.byone %>%
ggsave(filename = "_pics/Fig4_TopFirstPredictors_CI_cov.png", dpi=400,
- width=4, height=6, plot = top.first)
+ width=5, height=6, plot = top.first)
mydata <- corXY.ci %>%
filter(run=="seq") %>%
@@ -348,30 +367,40 @@ tt2 <- ttheme_minimal(
colhead = list(fg_params=list(cex = .7)),
rowhead = list(fg_params=list(cex = .7)))
-(topall.leg <- cowplot::plot_grid(top.all,
- tableGrob(trait.labs %>%
- mutate(Code=1:n()) %>%
- dplyr::select(Code, Trait=trait.name) %>%
- filter(Trait %in% traits.sign.alone.cov),
- theme=tt2, rows = NULL),
- #tableGrob(trait.labs %>%
- # mutate(Code=1:n()) %>%
- # dplyr::select(Code, Trait=trait.name) %>%
- # slice(33:n()), theme=tt2),
- nrow=1, rel_widths=c(0.70, 0.3)))
+
+ttlabs <- trait.labs %>%
+ mutate(Code=1:n()) %>%
+ filter(trait.name %in% traits.sign.alone.cov)
+
+tobold <- which(ttlabs$trait.name %in% best.traits.cov)
+tg <- tableGrob(ttlabs %>%
+ dplyr::select(Code, Trait=Long_English_name) %>%
+ mutate(Trait=replace(x = Trait,
+ list = Trait=="Vegetative Propagation - Fragmentation",
+ values = "Veg. Propag. - Fragmentation")),
+ theme=tt2, rows = NULL)
+
+## this needs to be triple checked
+for (i in (11 + tobold)) {
+ tg$grobs[[i]] <- editGrob(tg$grobs[[i]], gp=gpar(fontface="bold"))
+}
+
+
+(topall.leg <- cowplot::plot_grid(top.all, tg,
+ nrow=1, rel_widths=c(0.60, 0.4)))
ggsave(filename = "_pics/Fig5_Best_AllCombinations_CI_cov.png", dpi=400,
- width=6, height=3, topall.leg)
+ width=6, height=2, topall.leg)
break()
-
+###### ___________________ ######
###### 2.2 Presence/Absence data ####
####### 2.2.1 Import and adjust dataset ####
### merge output from one-by-one and all combo of sign traits
# ### all traits taken alone
-# myfilelist0 <- list.files(path="_derived/Mesobromion/PresAbs", pattern="HIDDENpa_[0-9]+_.RData", full.names = T)
+# myfilelist0 <- list.files(path="_derived/Mesobromion/PresAbs/onebyone", pattern="HIDDENpa-nona2_[0-9]+_.RData", full.names = T)
# dataFiles0 = purrr::map(myfilelist0, function(x){get(load(x))})
# #load("_derived/Mesobromion/PresAbs/HIDDEN_round_11.RData")
# corXY.alone = bind_rows(dataFiles0) %>%
@@ -383,10 +412,11 @@ break()
# mutate_at(.vars=vars(trait1),
# .funs=~factor(.,
# levels=trait.labs$Trait.comb,
-# labels=trait.labs$Short_english_name)) %>%
+# labels=trait.labs$trait.name)) %>%
# arrange(ntraits, Coef.obs) %>%
# dplyr::select(Trait.comb, Test, n, ntraits, everything()) %>%
-# mutate(Trait.comb=NA)
+# arrange(desc(sign_plus), desc(Coef.obs))
+# #mutate(Trait.comb=NA)
### combo of sign traits
@@ -421,7 +451,7 @@ traits.sign.alone.pa <- corXY.ci.pa %>%
### OLD ###[1] Leaf_Scleroph Height GF_Nanophan BL_Dur SLA LeafP LeafCdrymass CSR
#[9] LeafN SeedLength LDMC GF_Hemiphan
-# [1] Leaf_Scleroph Height GF_Nanophan SLA LeafCdrymass FP_Dur GF_Hemiphan LeafN
+### NEW ### [1] Leaf_Scleroph Height GF_Nanophan SLA LeafCdrymass FP_Dur GF_Hemiphan LeafN
# LeafP SeedLength GF_Macrophan
### Plot combinations one by one
@@ -577,7 +607,7 @@ ggsave(filename = "_pics/SXXX_Best_AllCombinations_CI_pa.png", dpi=400,
-
+###### ___________________ ######
##### 3. Calculate CWM ####
#species <- read_delim("_data/Mesobromion/species.out.10perc.txt", delim="\t")
species.cov <- read_delim("_data/Mesobromion/species.v2.10perc.cov.txt", delim="\t")
@@ -646,6 +676,7 @@ circleFun <- function(center = c(0,0),diameter = 1, npoints = 100){
return(data.frame(x = xx, y = yy))
}
dat <- circleFun(diameter = 2, npoints = 100)
+###### ___________________ ######
### 4.0 PCA of X (Fuzzy-weighted) matrix + CWMs ####
#best.traits.cov <- c("Height", "BL_Dur", "VP_Fragm", "CSR")
source("01b_MesobromionCluster.R")
@@ -654,7 +685,7 @@ W.fuzzy <- matrix.x(comm=species.cov %>%
traits=traits %>%
filter(species0 %in% colnames(species.cov)) %>%
column_to_rownames("species0") %>%
- dplyr::select(all_of(traits.sign.alone.cov)),
+ dplyr::select(all_of(best.traits.cov)),
scale = T)$matrix.X
pca.fuzz <- rda(W.fuzzy, scale = F)
varexpl <- round((pca.fuzz$CA$eig)/sum(pca.fuzz$CA$eig)*100,1)
@@ -690,9 +721,9 @@ basemap0 <- ggplot(data=as.data.frame(scores.pca*5)) +
aes(x=0, xend=PC1, y=0, yend=PC2, col=mycol),
arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
geom_path(data=dat,aes(x,y), col=gray(0.8), lwd=0.5) + #add correlation circle
- geom_point(data=myvectors %>%
- filter(categorical==1),
- aes(x=PC1, y=PC2, col=mycol), pch=16, size=1, alpha=0.8) +
+# geom_point(data=myvectors %>%
+# filter(categorical==1),
+# aes(x=PC1, y=PC2, col=mycol), pch=16, size=1, alpha=0.8) +
geom_label_repel(data=myvectors,# %>%
#mutate_if(~is.numeric(.), ~(.)*1.2),
aes(x=PC1, y=PC2, label=mylab, col=mycol, fontface=fontface0), size=2,
@@ -701,26 +732,26 @@ basemap0 <- ggplot(data=as.data.frame(scores.pca*5)) +
ylab(paste("PC2 (", varexpl[2], "%)", sep=""))
)
-PCA.fuzz2_3 <- basemap0 +
+PCA.fuzz1_3 <- basemap0 +
geom_point(data=as.data.frame(pca.fuzz$CA$u[,1:3]*5),
- aes(x=PC2, y=PC3), pch="+", size=2, alpha=0.8) +
+ aes(x=PC1, y=PC3), pch="+", size=2, alpha=0.8) +
geom_segment(data=myvectors%>%
filter(categorical==0),
- aes(x=0, xend=PC2, y=0, yend=PC3, col=mycol),
+ aes(x=0, xend=PC1, y=0, yend=PC3, col=mycol),
arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
geom_path(data=dat,aes(x,y), col=gray(0.8), lwd=0.5) + #add correlation circle
- geom_point(data=myvectors %>%
- filter(categorical==1),
- aes(x=PC2, y=PC3, col=mycol), pch=16, size=1, alpha=0.8) +
+# geom_point(data=myvectors %>%
+# filter(categorical==1),
+# aes(x=PC1, y=PC3, col=mycol), pch=16, size=1, alpha=0.8) +
geom_label_repel(data=myvectors,# %>%
#mutate_if(~is.numeric(.), ~(.)*1.2),
- aes(x=PC2, y=PC3, label=mylab, col=mycol, fontface=fontface0), size=2,
+ aes(x=PC1, y=PC3, label=mylab, col=mycol, fontface=fontface0), size=2,
position = position_dodge(1), segment.alpha=0.8, segment.colour=gray(0.7)) +
- xlab(paste("PC2 (", varexpl[2], "%)", sep="")) +
+ xlab(paste("PC1 (", varexpl[1], "%)", sep="")) +
ylab(paste("PC3 (", varexpl[3], "%)", sep=""))
-PC_fuzzy <- cowplot::plot_grid(PCA.fuzz1_2,PCA.fuzz2_3, nrow=1)
-ggsave("_pics/FigS991_PC_Fuzzy_1-3.png", width=10, height=5, dpi=300, last_plot())
+PC_fuzzy <- cowplot::plot_grid(PCA.fuzz1_2,PCA.fuzz1_3, nrow=1)
+ggsave("_pics/S11_PC_Fuzzy_1-3.png", width=10, height=5, dpi=300, last_plot())
#### 4.0.1 Alternative showing species scores ####
tmp <- as.data.frame(pca.fuzz$CA$v[,1:3]*7) %>%
@@ -760,35 +791,35 @@ PCAfuzz1_2.sp <- basemap0 %+% tmp +
xlab(paste("PC1 (", varexpl[1], "%)", sep="")) +
ylab(paste("PC2 (", varexpl[2], "%)", sep=""))
-PCAfuzz2_3.sp <- basemap0 %+% tmp +
- geom_text(aes(x=PC2, y=PC3, label=labels), size=2, alpha=0.8) +
+PCAfuzz1_3.sp <- basemap0 %+% tmp +
+ geom_text(aes(x=PC1, y=PC3, label=labels), size=2, alpha=0.8) +
geom_segment(data=myvectors %>%
filter(categorical==0),
- aes(x=0, xend=PC2, y=0, yend=PC3, col=mycol),
+ aes(x=0, xend=PC1, y=0, yend=PC3, col=mycol),
arrow = arrow(length = unit(0.08, "inches")), alpha=0.8) +
geom_path(data=dat,aes(x,y), col=gray(0.8), lwd=0.5) + #add correlation circle
# geom_point(data=myvectors %>%
# filter(categorical==1),
# aes(x=PC2, y=PC3, col=mycol), pch=16, size=1, alpha=0.8) +
geom_label_repel(data=myvectors,
- aes(x=PC2, y=PC3, label=mylab, col=mycol, fontface=fontface0), size=2,
+ aes(x=PC1, y=PC3, label=mylab, col=mycol, fontface=fontface0), size=2,
position = position_dodge(1), segment.alpha=0.8, segment.colour=gray(0.7), segment.size = 0.5) +
- xlab(paste("PC2 (", varexpl[2], "%)", sep="")) +
+ xlab(paste("PC1 (", varexpl[2], "%)", sep="")) +
ylab(paste("PC3 (", varexpl[3], "%)", sep=""))
-ggsave("_pics/FigS991a_PCA_Fuzzy_1-2_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz1_2.sp)
-ggsave("_pics/FigS991b_PCA_Fuzzy_2-3_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz2_3.sp)
-
+ggsave("_pics/S11a_PCA_Fuzzy_1-2_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz1_2.sp)
+ggsave("_pics/S11b_PCA_Fuzzy_1-3_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz1_3.sp)
+###### _ ######
#### 4.1 PCA of Y (Bealls) matrix + CWM ####
W.beals <- as.data.frame(beals(species.cov %>%
column_to_rownames("RELEVE_NR"),
include=T, type=2))
-write.table(W.beals, sep="\t", file="_derived/Mesobromion/MatrixY_Beals.csv")
+#write.table(W.beals, sep="\t", file="_derived/Mesobromion/MatrixY_Beals.csv")
pca.out <- rda(W.beals, scale=F)
varexpl <- round((pca.out$CA$eig)/sum(pca.out$CA$eig)*100,1)
scores.pca <- pca.out$CA$u[,1:4]
@@ -813,7 +844,7 @@ myvectors <- as.data.frame(env.cor) %>%
mutate(category="Trait")) %>%
bind_rows(as.data.frame(fuzz.cor) %>%
rownames_to_column("mylab") %>%
- mutate(mylab=paste0("FWM-", mylab)) %>%
+ mutate(mylab=paste0("FW-", mylab)) %>%
mutate(category="Fuzzy-Weighted")) %>%
mutate(fontface0=ifelse(mylab %in% best.traits.cov, "bold", "italic")) %>%
mutate(category=as.factor(category)) %>%
@@ -922,7 +953,7 @@ ggsave("_pics/Fig6a_PCA_Beals_1-2_wSpecies.png", width=8, height=8, dpi=300, PCA
ggsave("_pics/Fig6b_PCA_Beals_3-4_wSpecies.png", width=8, height=8, dpi=300, PCA3_4.sp)
-
+###### _ ######
#### 4.2 RDA of Beals ~ FWMs ####
RDA.beals <- rda(W.beals ~ scores(pca.fuzz, choices=1:3)$sites, scale=F)
# var explained by CONSTRAINED axes
@@ -999,7 +1030,6 @@ panel.RDA_beals <- cowplot::plot_grid(RDA1_2,RDA1_3, nrow=1)
ggsave("_pics/Fig6v2_RDA_Beals_1-2-3.png", width=10, height=5, dpi=300, panel.RDA_beals)
-
#### 4.2.1 Alternative showing species scores ####
tmp <- as.data.frame(scores(RDA.beals, choices = 1:3)$species*4) %>%
#tmp <- as.data.frame(RDA.beals$CCA$v*4) %>%
@@ -1049,7 +1079,7 @@ ggsave("_pics/Fig6v2b_RDA_Beals_1-3_wSpecies.png", width=8, height=8, dpi=300, R
-
+###### _ ######
#### 4.3 PCA 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))
@@ -1203,7 +1233,7 @@ write_delim(tmp %>%
-
+###### ___________________ ######
#### 5 Map of Plot distribution #####
library(rgdal)
library(sp)
@@ -1253,7 +1283,7 @@ ggsave(filename="_pics/S2_PlotDistribution.png", height
-
+###### ___________________ ######
#### 6 Correlation matrices ####
#### 6.1 Trait level #####
library(corrplot)
--
GitLab