From a86030d9ccd385e700fcc90f83e95a67781d330c Mon Sep 17 00:00:00 2001
From: Francesco Sabatini <francesco.sabatini@idiv.de>
Date: Fri, 28 Aug 2020 18:53:18 +0200
Subject: [PATCH] Corrected matrix.x function to account for factors

---
 02_Mesobromion_ExamineOutput.R | 42 +++++++++++++++++++---------------
 99_HIDDEN_functions.R          | 15 ++++++++----
 2 files changed, 34 insertions(+), 23 deletions(-)

diff --git a/02_Mesobromion_ExamineOutput.R b/02_Mesobromion_ExamineOutput.R
index da518f6..a2e1fd3 100644
--- a/02_Mesobromion_ExamineOutput.R
+++ b/02_Mesobromion_ExamineOutput.R
@@ -132,23 +132,24 @@ 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.tbl() %>% 
-#   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$Short_english_name)) %>% 
-#   arrange(ntraits, desc(Coef.obs)) %>% 
-#   dplyr::select(Trait.comb, Test, n, ntraits, everything()) %>% 
-#   mutate(Trait.comb=NA) %>% 
-#   mutate(run="alone")
+ 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
@@ -180,7 +181,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.ci %>% 
+traits.sign.alone.cov <- corXY.all.ci %>% 
   #filter(run=="alone") %>% 
   filter(ntraits==1) %>% 
   filter(sign_plus==T) %>% 
@@ -195,6 +196,11 @@ mydata.byone <- corXY.ci %>%
   arrange(Coef.obs) %>% 
   mutate(seq=1:n())
 
+write_csv(mydata.byone %>% 
+            dplyr::select(Trait.comb:sign_plus) %>% 
+            arrange(desc(sign_plus), desc(Coef.obs)), 
+          path = "_output/S992_SolutionsOneByOIne.cov.csv")
+
 (top.first <- ggplot(data=mydata.byone %>% 
                        mutate(sign_plus=fct_rev(as.factor(sign_plus)))) + 
   geom_segment(aes(x=q025, xend=q975, y=seq, yend=seq, 
diff --git a/99_HIDDEN_functions.R b/99_HIDDEN_functions.R
index d66bbd1..77c43e9 100644
--- a/99_HIDDEN_functions.R
+++ b/99_HIDDEN_functions.R
@@ -39,6 +39,9 @@ get.corXY.bootstrap <- function(comm, traits, trait.sel="all", bootstrap=199){
   
   syn.out.tmp <- matrix.x(comm=comm, traits=traits[,ii,drop=F], scale=T)$matrix.X #, asym.bin=binary.traits
   W.beals <- as.data.frame(beals(comm, include=T, type=2))
+  (RD.tmp <- dcor(W.beals, as.data.frame(syn.out.tmp))^2)
+  
+  
   # permute traits 
   ### Create vector of species to resample from, which exclude species with NAs
   traits.ii <- traits[,ii, drop=F]
@@ -57,7 +60,6 @@ get.corXY.bootstrap <- function(comm, traits, trait.sel="all", bootstrap=199){
   syn.out.perm.tmp <- matrix.x(comm=comm, traits=traits[index.traits[[bootstrap+1]],ii,drop=F], 
                                scale=T)$matrix.X #, asym.bin=binary.traits
   
-  (RD.tmp <- dcor(W.beals, as.data.frame(syn.out.tmp))^2)
   (RD.perm.tmp <- dcor(W.beals, as.data.frame(syn.out.perm.tmp))^2)
   corXY <- NULL
   corXY <- data.frame(Trait.comb=lab.tmp, bootstr.n=0, Test="RD", 
@@ -90,14 +92,17 @@ matrix.x <- function (comm, traits, scale = TRUE, ranks = TRUE, ord, notificatio
   if (any(vartype == "n")) {
     stop("\n trait must contain only numeric, binary or ordinal variables \n")
   }
+  if (any(vartype == "f")) {
+    warning("Warning: There is a factor trait - Using gower")
+  }
   if (missing(ord)) {
     for (i in 1:length(vartype)) {
       if (ranks & vartype[i] == "o") {
         traits[, i] <- rank(traits[, i], na.last = "keep")
+        traits[, i] <- as.numeric(traits[, i])
       }
-      traits[, i] <- as.numeric(traits[, i])
     }
-    traits <- as.matrix(traits)
+    #traits <- as.matrix(traits)
   }
   matrix.w <- sweep(comm, 1, rowSums(comm, na.rm = TRUE), "/")
   w.NA <- apply(matrix.w, 2, is.na)
@@ -114,8 +119,8 @@ matrix.x <- function (comm, traits, scale = TRUE, ranks = TRUE, ord, notificatio
       warning("Warning: NA in traits matrix", call. = FALSE)
     }
   }
-  if (scale) {
-    dist.traits <- FD::gowdis(traits, ...)
+  if (scale | any(vartype=="f")) {
+    dist.traits <- FD::gowdis(traits)#, ...)
     dist.traits <- as.matrix(dist.traits)
     dist.traits[which(rowSums(is.na(dist.traits))==(ncol(dist.traits)-1)),] <- NA   ### In case a species has NO trait info, replace the diag with NA to avoid propagation
     similar.traits <- 1 - as.matrix(dist.traits)
-- 
GitLab