Skip to content
Snippets Groups Projects
Commit 021f1ed4 authored by Francesco Sabatini's avatar Francesco Sabatini
Browse files

Aligned Figure numbers - Output for all combinations

parent 6db0af12
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
02_Mesobromion_ExamineOutput.R
+ 194
− 12
View file @ 021f1ed4
......@@ -348,7 +348,7 @@ mydata.best <- mydata %>%
write_csv(mydata.best %>%
dplyr::select(Trait.comb:sign_plus),
path = "_output/S5_BestSolutionTiers.cov.csv")
path = "_output/S9_BestSolutionTiers.cov.csv")
### Graph of all the best combinations with text legend
(top.all <- ggplot(data=(mydata.best %>%
......@@ -397,6 +397,188 @@ ggsave(filename = "_pics/Fig5_Best_AllCombinations_CI_cov.png", dpi=400,
width=6, height=2, topall.leg)
###### ________R1___________ ######
###### R1.2. Import output from Cluster ####
##### R1.2.0. Trait labs for significant traits
traits.sign.cov <- read_delim(file="_data/Mesobromion/traits.v2.10perc.cov.sign.txt", delim="\t")
traits.sign.cov <- traits.sign.cov %>%
as.data.frame() %>%
mutate_if(~is.character(.), .funs=~as.factor(.)) %>%
column_to_rownames("species0")
## adapt trait labs to sign traits only
trait.labs.sign.cov <- trait.labs %>%
filter(trait.name %in% colnames(traits.sign.cov)) %>%
arrange(match(trait.name, colnames(traits.sign.cov))) %>%
rename(Trait.comb.new=Trait.comb) %>%
mutate(Trait.comb=1:n()) %>%
dplyr::select(Trait.comb, everything(), -Trait.comb.new)
##### R1.2.1 Cover values ######
### sequential trait combo
myfilelist1 <- list.files(path="_derived/Mesobromion/Cover/R1_all/", pattern="HIDDENcov-nona2_[0-9]+_.RData", full.names = T)
dataFiles1 = purrr::map(myfilelist1, function(x){get(load(x))})
#load("_derived/Mesobromion/PresAbs/HIDDEN_round_11.RData")
corXY.all = bind_rows(dataFiles1) %>%
as_tibble() %>%
distinct()
corXY.all.ci <- get.ci(corXY.all)
corXY.all.ci <- corXY.all.ci %>%
mutate(Trait.comb2=Trait.comb) %>%
separate(Trait.comb2, into=paste0("trait", 1:7)) %>%
mutate_at(.vars=vars(trait1:trait7),
.funs=~factor(.,
levels=trait.labs.sign.cov$Trait.comb,
labels=trait.labs.sign.cov$trait.name)) %>%
arrange(ntraits, desc(Coef.obs)) %>%
#filter(ntraits>1) %>%
dplyr::select(Trait.comb, Test, n, ntraits, everything()) %>%
mutate(run="seq")
rm( dataFiles1) #dataFiles0,
### merge together
corXY.ci <- corXY.all.ci # %>%
mydata <- corXY.ci
######## R1.2.1.4 Best - Graph of r(XY) using best combination of traits at each level of interaction N ########
### extract best combinations of traits
top.one.by.one <- get.best(mydata, N=1, labs=trait.labs.sign.cov)
## Routine to extract the best combination at each level of interaction (up to max traits)
maxtraits <- 7
for(nn in 1:maxtraits){
if(nn==1) {
best.at.1 <- get.best(mydata, N=nn, labs=trait.labs.sign.cov)
newdata <- mydata %>%
filter_at(.vars=vars(trait1:trait7),
.vars_predicate = any_vars(. %in% best.at.1$trait.name | is.na(.)))
new.best.row <- newdata %>%
filter(Trait.comb==best.at.1$Trait.comb)
upper <- new.best.row$q975
lower <- new.best.row$q025
print(paste("new best at nn", nn, best.at.1$trait.name))
best.progr <- best.at.1$Trait.comb
}
if(nn>1){
better <- list()
better$Trait.comb <- newdata %>%
filter(ntraits==nn) %>%
filter(q025>upper) %>%
arrange(desc(Coef.obs)) %>%
slice(1) %>%
pull(Trait.comb)
if(length(better$Trait.comb>0)){
better$trait.name <- trait.labs.sign.cov %>%
filter(Trait.comb %in% strsplit(better$Trait.comb, split = "_")[[1]]) %>%
pull(trait.name)
newdata <- newdata %>%
rowwise() %>%
mutate(nmatching= sum(unlist(strsplit(Trait.comb, "_")) %in%
unlist(strsplit(better$Trait.comb, "_")),
na.rm=T)) %>%
ungroup() %>%
filter(nmatching==nn)
new.best.row <- newdata %>%
filter(Trait.comb==better$Trait.comb)
upper <- new.best.row$q975
lower <- new.best.row$q025
print(paste("new best at nn", nn, paste(better$trait.name, collapse=" ")))
best <- better
best.progr <- c(best.progr, better$Trait.comb)
}
}
}
best.traits.cov <- corXY.ci %>%
filter(as.character(Trait.comb)==best.progr[length(best.progr)]) %>%
dplyr::select(trait1:trait7) %>%
mutate_all(~as.character(.)) %>%
dplyr::select(colnames(.)[which(colSums(is.na(.))==0)])
best.traits.cov <- as.character(best.traits.cov[1,])
#"Leaf_Scleroph" "FP_Dur" "VP_Fragm" "Height" "SLA"
### Create dataset with best combinations + all the one-way combinations
mydata.best <- mydata %>%
#filter_at(.vars=vars(trait1:trait7),
# # .vars_predicate = all_vars(. %in% best$trait.name | is.na(.))) %>%
# .vars_predicate = all_vars(. %in% traits.sign.alone.cov | is.na(.))) %>%
filter(ntraits>1) %>%
filter(sign_plus==T) %>%
arrange(ntraits, Coef.obs) %>%
group_by(ntraits) %>%
slice(n()) %>%
ungroup() %>%
bind_rows(corXY.ci %>%
filter(run=="seq") %>%
filter(ntraits==1)) %>%
#filter(trait1 %in% traits.sign.alone.cov)) %>%
arrange(ntraits, Coef.obs) %>%
mutate(seq=1:n()) %>%
mutate(sign_plus=factor(Trait.comb %in% best.progr))
write_csv(mydata.best %>%
dplyr::select(Trait.comb:sign_plus),
path = "_output/R1.S5_BestSolutionTiers.cov_allcombos.csv")
### Graph of all the best combinations with text legend
(top.all <- ggplot(data=(mydata.best %>%
mutate(size0=.6+(as.numeric(sign_plus)-1)*.6))) +
geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq, col="a",
lwd=size0)) +
geom_point(aes(x=Coef.obs, y=seq), pch=15) +
scale_y_continuous(breaks=mydata.best$seq,
labels=mydata.best$Trait.comb, name=NULL) +
scale_x_continuous(name="RD correlation") +
scale_size_identity() +
theme_bw() +
theme(panel.grid.minor = element_blank(),
axis.text = element_text(size=7),
legend.position = "none"))
# create legend of names
tt2 <- ttheme_minimal(
core = list(fg_params=list(cex = .7),
padding=unit(c(1, 1), "mm")),
colhead = list(fg_params=list(cex = .7)),
rowhead = list(fg_params=list(cex = .7)))
ttlabs <- trait.labs.sign.cov %>%
mutate(Code=1:n())
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)
## Make significant traits bold
for (i in (11 + tobold)) {
tg$grobs[[i]] <- editGrob(tg$grobs[[i]], gp=gpar(fontface="bold"))
}
#arrange into a panel
(topall.leg <- cowplot::plot_grid(top.all, tg,
nrow=1, rel_widths=c(0.60, 0.4)))
ggsave(filename = "_pics/R1/Fig5_R1_Best_AllCombinations_CI_cov.png", dpi=400,
width=6, height=2, topall.leg)
break()
###### ___________________ ######
......@@ -755,7 +937,7 @@ PCA.fuzz1_3 <- basemap0 +
ylab(paste("PC3 (", varexpl[3], "%)", sep=""))
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())
ggsave("_pics/S13_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) %>%
......@@ -812,8 +994,8 @@ PCAfuzz1_3.sp <- basemap0 %+% tmp +
ylab(paste("PC3 (", varexpl[3], "%)", sep=""))
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)
ggsave("_pics/S13a_PCA_Fuzzy_1-2_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz1_2.sp)
ggsave("_pics/S13b_PCA_Fuzzy_1-3_wSpecies.png", width=8, height=8, dpi=300, PCAfuzz1_3.sp)
......@@ -833,7 +1015,7 @@ env.cor <- cor(env %>%
scores.pca, use = "pairwise.complete.obs") #double check
fuzz.cor <- cor(pca.fuzz$CA$u[,1:3], scores.pca)
sink("_output/S9_EnvFit_CWMs_env.txt")
sink("_output/S14_EnvFit_CWMs_env.txt")
cwms.cor
env.cor
fuzz.cor
......@@ -953,8 +1135,8 @@ PCA3_4.sp <- basemap0 %+% tmp +
ylab(paste("PC4 (", varexpl[4], "%)", sep=""))
ggsave("_pics/S10a_PCA_Beals_1-2_wSpecies.png", width=8, height=8, dpi=300, PCA1_2.sp)
ggsave("_pics/S10b_PCA_Beals_3-4_wSpecies.png", width=8, height=8, dpi=300, PCA3_4.sp)
ggsave("_pics/S14a_PCA_Beals_1-2_wSpecies.png", width=8, height=8, dpi=300, PCA1_2.sp)
ggsave("_pics/S14b_PCA_Beals_3-4_wSpecies.png", width=8, height=8, dpi=300, PCA3_4.sp)
###### _ ######
......@@ -1155,7 +1337,7 @@ PCA.t2 <- baseplot +
PC_traits <- cowplot::plot_grid(PCA.t1, PCA.t2, nrow=1)
ggsave("_pics/S6_PCA_Traits_1-4_only7.png", width=10, height=5, dpi=300, PC_traits)
ggsave("_pics/S10c_PCA_Traits_1-4_only7.png", width=10, height=5, dpi=300, PC_traits)
##### 4.3b Alternative version of figS6, showing the species ####
tmp <- as.data.frame(pca.scores[,1:4]*.2) %>%
......@@ -1200,8 +1382,8 @@ PCA.t2.sp <- baseplot %+% tmp +
xlab(paste("PC3 (", varexpl[3], "%)", sep="")) +
ylab(paste("PC4 (", varexpl[4], "%)", sep=""))
ggsave("_pics/S6a_PCA_Traits_1-2_wSpecies.png", width=8, height=8, dpi=300, PCA.t1.sp)
ggsave("_pics/S6b_PCA_Traits_3-4_wSpecies.png", width=8, height=8, dpi=300, PCA.t2.sp)
ggsave("_pics/S10_PCA_Traits_1-2_wSpecies.png", width=8, height=8, dpi=300, PCA.t1.sp)
ggsave("_pics/S10b_PCA_Traits_3-4_wSpecies.png", width=8, height=8, dpi=300, PCA.t2.sp)
#traits.dummy %>% filter(species0 %in% (tmp %>% filter(labels %in% c("Fes_pal", "Ses_alb", "Car_hum")) %>% pull(species))) %>% dplyr::select(any_of(starts_with(as.character(traits.sign.alone.cov))))
......@@ -1304,7 +1486,7 @@ traits7 <- traits %>%
# relocate(all_of(starts_with(as.character(best.traits.cov))), everything())
res <- cor(traits7, use = "pairwise.complete.obs")
png(file="_pics/S7_Correlations_Trait.png", width=8, height=6.5, units = "in", res=300)
png(file="_pics/S11_Correlations_Trait.png", width=8, height=6.5, units = "in", res=300)
corrplot(res, type = "upper",
tl.col = "black", tl.srt = 45, number.cex=0.6, addCoef.col = "black", diag=F)
dev.off()
......@@ -1315,7 +1497,7 @@ res2 <- cor(CWM.wide %>%
dplyr::select(any_of(traits.sign.alone.cov)) %>% ## caution selecting only numerical variables
dplyr::select(sort(tidyselect::peek_vars())) %>%
relocate(any_of(best.traits.cov), everything()))
png(file="_pics/S8_Correlations_CWMs.png", width=8, height=6.5, units = "in", res=300)
png(file="_pics/S12_Correlations_CWMs.png", width=8, height=6.5, units = "in", res=300)
corrplot(res2, type = "upper",
tl.col = "black", tl.srt = 45, number.cex=0.6, addCoef.col = "black", diag=F)
dev.off()
......
03_Figures_Simulations.R
+ 217
− 93
View file @ 021f1ed4
This diff is collapsed.
04_Additional_Figs.R
+ 51
− 37
View file @ 021f1ed4
......@@ -122,6 +122,7 @@ do.the.parse <- function(toparse) {
}
### Figure R1.S8 - Comparison beta vs cor ####
## Function to create Figure SXXV
create.panel <- function(x){
gg.betaW <- ggplot(data=x %>%
......@@ -132,7 +133,6 @@ create.panel <- function(x){
geom_density(aes(value)) +
xlab("Proportional Beta Diversity (W)") +
xlim(c(0,1.1)) +
ylim(c(0,60)) +
theme_classic()
gg.betaY <- gg.betaW %+%
......@@ -149,7 +149,8 @@ create.panel <- function(x){
group_by(dataset, metric, matrix) %>%
summarize(value=max(abs(value))))+
xlab("Cor(WE)") +
ylab(NULL)
ylab(NULL) +
ylim(c(0,13.5))
gg.corY <- gg.betaW %+%
(x %>%
......@@ -159,56 +160,62 @@ create.panel <- function(x){
group_by(dataset, metric, matrix) %>%
summarize(value=max(abs(value)))) +
xlab("Cor(YE)") +
ylab(NULL)
ylab(NULL)+
ylim(c(0,13.5))
gg.panel <- cowplot::plot_grid(gg.betaW, gg.corW,
gg.betaY, gg.corY,
nrow=2, rel_widths = c(1,1.06))
nrow=2, rel_widths = c(1,1.06), align = "hv")
return(gg.panel)
}
### Figure SXXV - comparison beta vs cor ####
# Set path of files to import and parse
mypath <- "_data/Experiment_30Oct2020_FactorInteraction&TraitCorr_XY_SampleSize_Main=040_Inter=00_Corr=00_v21169"
myfiles <- list.files(path=mypath, pattern = "FinalSimulatedData.txt", recursive = T, full.names = T)
myfiles <- myfiles[grepl("_new", x=myfiles)]
# Loop over files to import and parse. Create graphs
for(i in 1:length(myfiles)){
toparse <- myfiles[i]
sampleN <- regmatches(toparse, gregexpr("N=[[:digit:]]+", toparse))[[1]]
div.summary <- do.the.parse(toparse)
gg.out <- create.panel(div.summary)
ggsave(filename = paste0("_pics/R1/FigSXXZ_Panel_BetaCor_", sampleN, ".png"),
ggsave(filename = paste0("_pics/R1/FigS8_Panel_BetaCor_", sampleN, ".png"),
width=6, height=5, device="png", dpi = 300, plot = gg.out)
}
#### Figure SXXK - Comparison abg across combinatios Inter X Corr (Main=0.3) ####
#### Figure R1.S1d - Comparison abg across combinatios Inter X Corr (Main=0.3) ####
create.panel2 <- function(xx, i, tit){
require(ggpubr)
# alpha, min, mean, max
gg.alpha <- ggplot(data=xx %>%
filter(matrix=="W") %>%
filter(metric %in% c("richness", "alpha")) %>%
mutate(with=ifelse(!is.na(with), paste0("OTU Rich (", with, ")"), with)) %>%
mutate(with=ifelse(is.na(with), "Eq. OTU (mean)", with))) +
geom_density(aes(value, group=with, col=with), alpha=0.7, show.legend=FALSE)+
stat_density(aes(x=value, colour=with),
geom="line",position="identity") +
scale_color_brewer(palette="Dark2", name=NULL) +
theme_classic() +
xlab("Alpha diversity (OTU Richness)") +
theme(legend.position = c(0.65, 0.9)#,
#legend.text = element_text(size=7)
) +
gg.rich <- ggplot(data=xx %>%
filter(matrix=="W") %>%
filter(metric %in% c("richness"))) +
geom_density(aes(value, group=with, col=with), alpha=0.7, show.legend=FALSE)+
stat_density(aes(x=value, colour=with),
geom="line",position="identity") +
scale_color_brewer(palette="Dark2", name=NULL) +
theme_classic() +
xlab("Species richness") +
theme(legend.position = c(0.75, 0.9)) +
xlim(c(-.1,100)) +
ylim(c(0,0.25))
ylim(c(0,0.2))
# alpha, min, mean, max
gg.alpha <- gg.rich %+%
(xx %>%
filter(matrix == "W") %>%
filter(metric == "alpha")) +
geom_density(aes(value), show.legend = FALSE) +
xlab("Mean alpha diversity") +
ylim(c(0, 0.07)) +
ylab(NULL)
# beta
gg.beta <- ggplot(data=xx %>%
......@@ -225,7 +232,7 @@ create.panel2 <- function(xx, i, tit){
(xx %>%
filter(matrix=="W") %>%
filter(metric %in% "propbeta")) +
xlab("Proportional Beta Diversity") +
xlab("Proportional beta Diversity") +
xlim(c(0,0.7)) +
ylim(c(0,4))
#gamma
......@@ -243,12 +250,14 @@ create.panel2 <- function(xx, i, tit){
))
if(i!=1){
gg.alpha <- gg.alpha +
gg.rich <- gg.rich +
theme(legend.position="none")
}
}
if(i!=3){
gg.rich <- gg.rich +
xlab(NULL) +
theme(axis.text.x = element_blank())
gg.alpha <- gg.alpha +
xlab(NULL) +
theme(axis.text.x = element_blank())
......@@ -263,27 +272,32 @@ create.panel2 <- function(xx, i, tit){
theme(axis.text.x = element_blank())
}
gg.panel <- cowplot::plot_grid(gg.title, gg.alpha, gg.beta, gg.propbeta, gg.gamma,
nrow=1, rel_widths = c(0.08, 1,.94, .94, .94))
gg.panel <- cowplot::plot_grid(gg.title, gg.rich, gg.alpha, gg.beta, gg.propbeta, gg.gamma,
nrow=1, rel_widths = c(0.08, 1, .94, .94, .94, .94))
return(gg.panel)
}
## set path of files to import and parse
mypath <- "_data/Experiment_30Oct2020_FactorInteraction&TraitCorr_XY_DataExamples/"
myfiles <- list.files(path=mypath, pattern = "FinalSimulatedData.txt", recursive = T, full.names = T)
#loop of files, parse and create graphs
panel.list <- list()
for(i in 1:length(myfiles)){
toparse <- myfiles[i]
Inter <- regmatches(toparse, gregexpr("Inter=[[:digit:]]+", toparse))[[1]]
Inter <- gsub(pattern="=0", replacement=" = 0.", x = Inter)
Inter <- gsub(pattern="0.0", replacement="0", x = Inter)
Corr <- regmatches(toparse, gregexpr("Corr=[[:digit:]]+", toparse))[[1]]
Corr <- gsub(pattern="=0", replacement=" = 0.", x = Corr)
Corr <- gsub(pattern="0.0", replacement="0", x = Corr)
div.summary <- do.the.parse(toparse)
panel.list[[i]] <- create.panel2(div.summary, i=i,tit=paste(Inter, Corr))
}
gg.SXXK <- cowplot::plot_grid(plotlist=panel.list, nrow=3)
ggsave(filename = "_pics/R1/FigSXXK_Panel_abg.png",
width=10, height=6, device="png", dpi = 300, plot = gg.SXXK)
gg.S1d <- cowplot::plot_grid(plotlist=panel.list, nrow=3)
ggsave(filename = "_pics/R1/FigS1d_Panel_abg.png",
width=10, height=6, device="png", dpi = 300, plot = gg.S1d)
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment