diff --git a/01_Mesobromion.R b/01_Mesobromion.R
index 9d6bf71f46292ea21ebf34a52e70ee765adebeec..09c41024666c06b6ce648e95ab119259ea702e07 100644
--- a/01_Mesobromion.R
+++ b/01_Mesobromion.R
@@ -169,127 +169,69 @@ traits <- traits %>%
#### ## Import output ####
-
myfilelist <- list.files(path="_derived/Mesobromion/", pattern="HIDDEN(\\.[A-Za-z]+|_[A-Za-z]+)_[0-9]+.RData", full.names = T)
-corXY <- NULL
-#corXY.perm <- NULL
-#corXY.bootstr <- NULL
-for(ff in myfilelist){
- index <- as.numeric(regmatches(ff, gregexpr("[[:digit:]]+", ff))[[1]])
- load(ff)
- corXY <- bind_rows(corXY, cor.obs)
- #corXY.perm <- rbind(corXY.perm, cor.perm)
- #corXY.bootstr <- rbind(corXY.bootstr, cor.bootstr)
- print(index)
-}
-
-
-dataFiles = map(myfilelist, load)
-dat = bind_rows(dataFiles)
-
-
+dataFiles = purrr::map(myfilelist, function(x){get(load(x))})
+corXY = bind_rows(dataFiles) %>%
+ as.tbl()
+rm(corXY.output, dataFiles)
-diffs <- corXY.bootstr - corXY.perm
-sign <- apply(diffs, MARGIN=1, function(x){sum(x>0)>190})
-sign <- data.frame(Trait.comb=names(sign), sign=sign)
-
-corXY %>%
- arrange(Test, desc(Coef))
-
-aa <- data.frame(Trait.comb=paste0("t", 1:ncol(traits)),
+trait.labs <- data.frame(Trait.comb=1:ncol(traits),
trait.name=colnames(traits))#[-which(colnames(traits) %in% c("species", "species0"))])
-corXY.ci <- data.frame(Trait.comb=rownames(corXY.bootstr),
- obs=corXY.bootstr[,1],
- q025=apply(corXY.bootstr[,-1], MARGIN=1, quantile, probs=0.025),
- q975=apply(corXY.bootstr[,-1], MARGIN=1, quantile, probs=0.975)) %>%
- left_join(sign, by="Trait.comb") %>%
- as.tbl() %>%
+## calculate confidence intervals for bootstrapped correlations
+corXY.ci <- corXY %>%
+ filter(bootstr.n==0) %>%
+ dplyr::select(-bootstr.n) %>%
+ group_by(Trait.comb) %>%
+ slice(1) %>%
+ left_join(corXY %>%
+ filter(bootstr.n!=0) %>%
+ group_by(Trait.comb) %>%
+ summarize(q025=quantile(Coef.obs, 0.025),
+ q975=quantile(Coef.obs, 0.975),
+ greater.than.perm=mean(Coef.obs>Coef.perm),
+ n=n())) %>%
+ #calculate significance using permuted correlations
+ mutate(sign_plus=greater.than.perm>0.975) %>%
+ mutate(sign_minus=greater.than.perm<0.025) %>%
mutate(Trait.comb2=Trait.comb) %>%
- separate(Trait.comb2, c("trait1", "trait2")) %>% #, "trait3")) %>%
- mutate(trait2=paste0("t", trait2)) %>%
- #mutate(trait3=paste0("t", trait3)) %>%
- left_join(aa %>%
- rename(trait1=Trait.comb), by="trait1") %>%
- dplyr::select(-trait1) %>%
- rename(trait1=trait.name) %>%
- left_join(aa %>%
- rename(trait2=Trait.comb), by="trait2") %>%
- dplyr::select(-trait2) %>%
- rename(trait2=trait.name) %>%
- #left_join(aa %>%
- # rename(trait3=Trait.comb), by="trait3") %>%
- #dplyr::select(-trait3) %>%
- #rename(trait3=trait.name) %>%
- mutate(ntraits= (!is.na(trait1)) + (!is.na(trait2))) %>%
- arrange(ntraits, obs)
+ #split strings of trait combinations and add labels
+ separate(Trait.comb2, into=paste0("trait", 1:7)) %>%
+ mutate_at(.vars=vars(trait1:trait7),
+ .funs=~factor(., levels=trait.labs$Trait.comb, labels=trait.labs$trait.name)) %>%
+ rowwise() %>%
+ mutate(ntraits= length(unlist(strsplit(Trait.comb, "_")))) %>%
+ ungroup() %>%
+ arrange(ntraits, Coef.obs) %>%
+ dplyr::select(Trait.comb, Test, n, ntraits, everything())
-### spider plot
-# faceted spider plot
-# mydata <- corXY.ci %>%
-# filter(grepl(pattern = "t2", Trait.comb)) %>%
-# #filter(grepl(pattern = "_6", Trait.comb)) %>%
-# #filter(sign==T) %>%
-# filter(ntraits<=2) %>%
-# mutate(Trait.comb=factor(Trait.comb)) %>%
-# mutate(seq=1:nrow(.))
-
-
-
-corXY.SES <- get.SES(obs.df = corXY, perm.df = corXY.perm, stat="RV") %>%
- mutate(Trait.comb0=Trait.comb) %>%
- separate(Trait.comb0, c("trait1", "trait2")) %>% #, "trait3")) %>%
- mutate(trait2=paste0("t", trait2)) %>%
- #mutate(trait3=paste0("t", trait3)) %>%
- left_join(aa %>%
- rename(trait1=Trait.comb), by="trait1") %>%
- dplyr::select(-trait1) %>%
- rename(trait1=trait.name) %>%
- left_join(aa %>%
- rename(trait2=Trait.comb), by="trait2") %>%
- dplyr::select(-trait2) %>%
- rename(trait2=trait.name) %>%
-# left_join(aa %>%
-# rename(trait3=Trait.comb), by="trait3") %>%
-# dplyr::select(-trait3) %>%
-# rename(trait3=trait.name) %>%
- dplyr::select(Trait.comb, trait1, trait2, q025:conf.p) %>% #trait3,
- arrange(desc(SES.np))
-print(corXY.SES, n=20)
-
-
-mydata.SES <- corXY.ci %>%
- left_join(corXY.SES %>%
- mutate(sign.SES=ifelse(smaller.than.obs>0.975, T, F)) %>% # smaller.than.obs>0.995 |
- dplyr::select(Trait.comb, sign.SES),
- by="Trait.comb")
-
-
-traits.sign.alone <- mydata.SES %>%
+traits.sign.alone <- corXY.ci %>%
filter(ntraits==1) %>%
- filter(sign.SES==T) %>%
+ filter(sign_plus==T) %>%
pull(trait1)
-mydata <- mydata.SES %>%
+mydata <- corXY.ci %>%
filter(ntraits==1)%>%
- mutate(seq=1:nrow(.))
+ mutate(seq=1:n())
ggplot(data=mydata) +
- geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq, col=sign.SES)) +
- geom_point(aes(x=obs, y=seq), pch=15) +
+ geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq,
+ col=sign_plus)) +
+ geom_point(aes(x=Coef.obs, y=seq), pch=15) +
scale_y_continuous(breaks=mydata$seq,
labels=mydata$Trait.comb)
-mydata <- mydata.SES %>%
+mydata <- mydata %>%
#filter(ntraits==2)%>%
- filter(sign.SES==T) %>%
+ filter(sign_plus==T) %>%
filter(trait1 %in% traits.sign.alone & (is.na(trait2) | trait2 %in% traits.sign.alone)) %>%
mutate(seq=1:nrow(.)) # %>%
# arrange(desc(obs)) %>%
# slice(1:20)
ggplot(data=mydata) +
- geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq, col=sign.SES)) +
- geom_point(aes(x=obs, y=seq), pch=15) +
+ geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq,
+ col=sign_plus)) +
+ geom_point(aes(x=Coef.obs, y=seq), pch=15) +
scale_y_continuous(breaks=mydata$seq,
labels=mydata$Trait.comb,
name="RV correlation")
@@ -297,9 +239,9 @@ ggplot(data=mydata) +
-traits.sign <- traits %>%
- dplyr::select(species0, traits.sign.alone)
-write_delim(traits.sign, path="_data/Mesobromion/traits.out.10perc.sign.txt", delim="\t")
+#traits.sign <- traits %>%
+# dplyr::select(species0, traits.sign.alone)
+#write_delim(traits.sign, path="_data/Mesobromion/traits.out.10perc.sign.txt", delim="\t")
break()