Repository 'ecology_stat_presence_abs'
hg clone https://toolshed.g2.bx.psu.edu/repos/ecology/ecology_stat_presence_abs

Changeset 1:9a2e0195bb43 (2021-08-13)
Previous changeset 0:e5552099d0e3 (2021-07-27)
Commit message:
"planemo upload for repository https://github.com/Marie59/Data_explo_tools commit 60627aba07951226c8fd6bb3115be4bd118edd4e"
modified:
graph_stat_presence_abs.r
stat_presence_abs.xml
added:
funct_anomy.r
functions.r
graph_lcbd.r
test-data/SCBD.txt
b
diff -r e5552099d0e3 -r 9a2e0195bb43 funct_anomy.r
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/funct_anomy.r Fri Aug 13 18:18:21 2021 +0000
[
@@ -0,0 +1,39 @@
+#Rscript
+
+###########################
+##     Anonymization     ##
+###########################
+
+#####Packages : tangles
+
+#Load arguments
+
+args <- commandArgs(trailingOnly = TRUE)
+
+if (length(args) == 0) {
+    stop("This tool needs at least one argument")
+}else{
+    table <- args[1]
+    hr <- args[2]
+    latitude <- as.numeric(args[3])
+    longitude <- as.numeric(args[4])
+}
+
+if (hr == "false") {
+  hr <- FALSE
+}else{
+  hr <- TRUE
+}
+
+#####Import data
+data <- read.table(table, sep = "\t", dec = ".", header = hr, fill = TRUE, encoding = "UTF-8")
+
+randomized_data <- tangles::tangles(data = as.matrix(data[, c(latitude, longitude)]), depth = 3, rasterdata = FALSE, raster_object = FALSE, saveTangles = FALSE, path = NULL)

+data[, c(latitude, longitude)] <- NULL
+
+tab_anon <- data.frame(longitude = randomized_data[[1]]$X, latitude = randomized_data[[1]]$Y)
+
+tab_anon <- cbind(data, tab_anon)
+
+write.table(tab_anon, "anonym_data.tabular", row.names = FALSE, quote = FALSE, sep = "\t", dec = ".", fileEncoding = "UTF-8")
b
diff -r e5552099d0e3 -r 9a2e0195bb43 functions.r
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/functions.r Fri Aug 13 18:18:21 2021 +0000
[
@@ -0,0 +1,22 @@
+#Rscript
+
+#########################################################################################
+####################### Exploration data tools function #################################
+#########################################################################################
+#### Based on Romain Lorrillière R script
+#### Modified by Alan Amosse, Benjamin Yguel and Marie Jossé for integrating within Galaxy-E
+
+######################################### start of the function makeTableAnalyse
+##Species are placed in separated columns and addition of zero on plots where at least one selected species is present
+make_table_analyse <- function(data, var, spe, var2, var3) {
+    tab <- reshape(data
+                  , v.names = var
+                  , idvar = c(var2, var3)
+                  , timevar = spe
+                  , direction = "wide")
+    tab[is.na(tab)] <- 0 ###### remplace les na par des 0 / replace NAs by 0
+
+    colnames(tab) <- sub(paste0(var, "."), "", colnames(tab))### remplace le premier pattern "abond." par le second "" / replace the column names "abond." by ""
+    return(tab)
+}
+######################################### end of the function makeTableAnalyse
b
diff -r e5552099d0e3 -r 9a2e0195bb43 graph_lcbd.r
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/graph_lcbd.r Fri Aug 13 18:18:21 2021 +0000
[
b'@@ -0,0 +1,225 @@\n+#Rscript\n+\n+#########################\n+##    Beta diversity   ##\n+#########################\n+\n+#####Packages : ggplot2\n+#               vegan\n+#               adespatial\n+#               dplyr\n+#               tibble\n+#               tdyr\n+\n+#####Load arguments\n+\n+args <- commandArgs(trailingOnly = TRUE)\n+\n+if (length(args) < 2) {\n+    stop("This tool needs at least 2 arguments")\n+}else{\n+    table <- args[1]\n+    hr <- args[2]\n+    abund <- as.numeric(args[3])\n+    loc <- as.numeric(args[4])\n+    spe <- as.numeric(args[5])\n+    date <- as.numeric(args[6])\n+    map <- as.logical(args[7])\n+    sepa <- as.logical(args[8])\n+    not <- as.logical(args[9])\n+    lat <- as.numeric(args[10])\n+    long <- as.numeric(args[11])\n+    var <- as.numeric(args[12])\n+    source(args[13])\n+}\n+\n+if (hr == "false") {\n+  hr <- FALSE\n+}else{\n+  hr <- TRUE\n+}\n+\n+#####Import data\n+data <- read.table(table, sep = "\\t", dec = ".", header = hr, fill = TRUE, encoding = "UTF-8")\n+colabund <- colnames(data)[abund]\n+colloc <- colnames(data)[loc]\n+if (map) {\n+  collat <- colnames(data)[lat]\n+  collong <- colnames(data)[long]\n+}\n+colspe <- colnames(data)[spe]\n+coldate <- colnames(data)[date]\n+data[, coldate] <- as.factor(data[, coldate])\n+\n+data <- data[grep("^$", data[, spe], invert = TRUE), ]\n+\n+if (sepa) {\n+colvar <- colnames(data)[var]\n+}\n+\n+# Data for species\n+data_num <- make_table_analyse(data, colabund, colspe, colloc, coldate)\n+nb_spe <- length(unique(data[, spe]))\n+nb_col <- ncol(data_num) - nb_spe + 1\n+\n+#Data with coordinates and environmental\n+if (map) {\n+  data_xy <- data_num[, c(collat, collong)]\n+  colnames(data_xy) <- c("latitude", "longitude")\n+  # Data for environment\n+  data_env <- data_num[, c(colloc, collat, collong)]\n+  colnames(data_env) <- c("site", "latitude", "longitude")\n+}\n+\n+# Data with only species and their abundance\n+data_spe <- data_num[, nb_col:ncol(data_num)]\n+rownames(data_spe) <- paste0(data_num[, colloc], " - ", data_num[, coldate])\n+\n+#####Your analysis\n+\n+# Computation beta.div {adespatial}\n+# Beta.div on Hellinger-transformed species data\n+data_beta <- adespatial::beta.div(data_spe, method = "hellinger", nperm = 9999)\n+\n+save(data_beta, file = "beta_diversity.Rdata")\n+cat("##############################################################################",\n+    "\\n########################### Beta Diversity Summary ###########################",\n+    "\\n##############################################################################",\n+    "\\n\\n### All data ###",\n+    "\\nBeta diversity: ", data_beta$beta[[2]],\n+    "\\nSum of Squares: ", data_beta$beta[[1]],\n+    "\\n\\n### Vector of Local Contributions to Beta Diversity (LCBD) for the sites each date ###",\n+    "\\n", capture.output(data_beta$LCBD),\n+    "\\n\\n### Vector of P-values associated with the LCBD indices ###",\n+    "\\n", capture.output(data_beta$p.LCBD),\n+    "\\n\\n### Vector of Corrected P-values for the LCBD indices, Holm correction ###",\n+    "\\n", capture.output(data_beta$p.adj),\n+    "\\n\\n### Vector of Species contributions to beta diversity (SCBD) ###",\n+    "\\n", capture.output(data_beta$SCBD), file = "LCBD.txt", fill = 1, append = TRUE)\n+\n+# Which species have a SCBD larger than the mean SCBD?\n+scbd <- capture.output(data_beta$SCBD[data_beta$SCBD >= mean(data_beta$SCBD)])\n+write(scbd, "SCBD.txt")\n+\n+##1st fonction\n+beta_div_ext <- function(data_beta, data_xy, data_env) {\n+   data_beta_ext <- data.frame(data_xy, data_env, LCBD = data_beta$LCBD * 100, p.LCBD = data_beta$p.LCBD, signif = data_beta$p.LCBD < 0.05)\n+\n+  graph_beta_ext <- ggplot2::ggplot(data = data_beta_ext, ggplot2::aes(x = latitude, y = longitude, size = LCBD, col = signif)) +\n+  ggplot2::geom_point() +\n+  ggplot2::scale_colour_manual(values = c("#57bce0", "#ce0b0b"), labels = c("Non significant", "Significant"), name = "Significance at 0.05") +\n+  ggplot2::xlab("Longitude") + ggplot2::ylab("Latitude")\n+\n+  ggplot2::ggsave("Beta_diversity_through_space.png", graph_beta_ext)\n+}\n+\n'..b'###\n+\n+####LCBD####\n+lcbd_site <- adespatial::beta.div(data_spe, "hellinger", nperm = 999)\n+\n+compute_lcbd <- function(data_beta, data_spe, data_num) {\n+\n+#############\n+  mat_lcbd_site <- data.frame(data_spe, LCBD = data_beta$LCBD * 100, p.LCBD = data_beta$p.LCBD, signif = data_beta$p.LCBD < 0.05, site = data_num[, colloc], date = data_num[, coldate])\n+\n+## Map spatio-temp\n+##################\n+  p1 <- ggplot2::qplot(date, site, size = LCBD, col = signif, data = mat_lcbd_site)\n+  p1 <- p1 + ggplot2::scale_colour_manual(values = c("#57bce0", "#ce0b0b"), labels = c("Non significant", "Significant"), name = "Significance at 0.05")\n+  p1 <- p1 + ggplot2::theme_bw() + ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 90)) + ggplot2::xlab("Date") + ggplot2::ylab("Site")\n+\n+  ggplot2::ggsave("LCBD_sites_time.png", p1)\n+\n+\n+## Par ann\xc3\xa9es\n+#############\n+  mean_time <- tapply(mat_lcbd_site$LCBD, mat_lcbd_site$date, mean)\n+  sd_time <- tapply(mat_lcbd_site$LCBD, mat_lcbd_site$date, sd)\n+  date <- unique(mat_lcbd_site$date)\n+\n+  data <- data.frame(date, mean_time, sd_time)\n+\n+  time <- ggplot2::ggplot() + ggplot2::geom_pointrange(ggplot2::aes(x = date, y = mean_time, ymin = mean_time - sd_time, ymax = mean_time + sd_time), data = data)\n+  time <- time + ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 90), axis.line.y = ggplot2::element_line(size = 0.5)) + ggplot2::ylab("mean LCBD")\n+\n+  ggplot2::ggsave("Mean_LCBD_through_time.png", time)\n+}\n+\n+## Choose another graph\n+#######################\n+compute_lcbd2 <- function(data_beta, data_spe, data_num) {\n+\n+#############\n+  mat_lcbd_site <- data.frame(data_spe, LCBD = data_beta$LCBD * 100, p.LCBD = data_beta$p.LCBD, signif = data_beta$p.LCBD < 0.05, site = data_num[, colloc], date = data_num[, coldate], variable = data_num[, colvar])\n+\n+  p1 <- ggplot2::qplot(date, variable, size = LCBD, col = signif, data = mat_lcbd_site)\n+  p1 <- p1 + ggplot2::scale_colour_manual(values = c("#57bce0", "#ce0b0b"), labels = c("Non significant", "Significant"), name = "Significance at 0.05")\n+  p1 <- p1 + ggplot2::theme_bw() + ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 90)) + ggplot2::xlab("Date") + ggplot2::ylab(colvar)\n+\n+  ggplot2::ggsave(paste0("LCBD_per_", colvar, "_through_time.png"), p1)\n+}\n+\n+####SCBD###\n+# Function to compute SCBD\n+library(dplyr)\n+make_scbd_uvc <- function(data_spe, z, data_beta) {\n+  # Computation using beta.div {adespatial} on\n+  # Hellinger-transformed species data\n+\n+  # Which species have a SCBD larger than the mean SCBD?\n+  spe_scbd <-  data_beta$SCBD[data_beta$SCBD >= mean(data_beta$SCBD)] %>%\n+    as.data.frame() %>%\n+    tibble::rownames_to_column(var = "Taxon") %>%\n+    dplyr::mutate("Methode" = z)\n+\n+  return(spe_scbd)\n+}\n+\n+# Function to make a radar plot\n+\n+coord_radar <- function(theta = "x", start = 0, direction = 1) {\n+  theta <- match.arg(theta, c("x", "y"))\n+  r <- if (theta == "x") "y" else "x"\n+  ggplot2::ggproto("CordRadar", ggplot2::coord_polar(theta = theta, start = start,\n+          direction = sign(direction)),\n+          is_linear = function(coord) TRUE)\n+}\n+\n+# Make the radar plot\n+radar_plot <- function(scbd_uvc_tc) {\n+  uvc_rd_plot_data <- scbd_uvc_tc %>%\n+    rename(scbd_score = ".")\n+\n+  rad_uvc <- ggplot2::ggplot(uvc_rd_plot_data, ggplot2::aes(x = Taxon, y = scbd_score, group = Methode)) +\n+    ggplot2::geom_line() +\n+    ggplot2::geom_point(size = 3) +\n+    coord_radar() +\n+    ggplot2::theme_bw() +\n+    ggplot2::theme(axis.text.x = ggplot2::element_text(size = 10),\n+        legend.position = "bottom")\n+\n+  ggplot2::ggsave("SCBD_Species_Radar_plot.png", rad_uvc)\n+}\n+\n+## LCBD\n+\n+if (map) {\n+  #Beta diversity\n+  beta_div_ext(data_beta, data_xy, data_env)\n+}\n+\n+#Lcbd per places and time\n+compute_lcbd(data_beta, data_spe, data_num)\n+\n+#Lcbd of your choice\n+if (sepa) {\n+  compute_lcbd2(data_beta, data_spe, data_num)\n+}\n+\n+##SCBD\n+\n+scbd_uvc_tc <- make_scbd_uvc(data_spe, z = "TC", data_beta)\n+\n+radar_plot(scbd_uvc_tc)\n'
b
diff -r e5552099d0e3 -r 9a2e0195bb43 graph_stat_presence_abs.r
--- a/graph_stat_presence_abs.r Tue Jul 27 16:57:02 2021 +0000
+++ b/graph_stat_presence_abs.r Fri Aug 13 18:18:21 2021 +0000
[
@@ -20,6 +20,7 @@
     spe <- as.numeric(args[4])
     loc <- as.numeric(args[5])
     time <- as.numeric(args[6])
+    source(args[7])
 }
 
 if (hr == "false") {
@@ -37,7 +38,6 @@
 coltime <- colnames(data)[time]
 
 data <- data[grep("^$", data[, spe], invert = TRUE), ]
-time <- as.integer(substring(data[, time], first = 1, last = 4))
 
 #####Your analysis
 
@@ -75,17 +75,7 @@
 #### Zero problem in data ####
 
 #Put data in form
-make_table_analyse <- function(data, var, spe, var2, var3) {
-    tab <- reshape(data
-                  , v.names = var
-                  , idvar = c(var2, var3)
-                  , timevar = spe
-                  , direction = "wide")
-    tab[is.na(tab)] <- 0 ###### remplace les na par des 0 / replace NAs by 0
 
-    colnames(tab) <- sub(var, "", colnames(tab))### remplace le premier pattern "abond." par le second "" / replace the column names "abond." by ""
-    return(tab)
-}
 data_num <- make_table_analyse(data, colvar, colspe, colloc, coltime)
 nb_spe <- length(unique(data[, spe]))
 nb_col <- ncol(data_num) - nb_spe + 1
b
diff -r e5552099d0e3 -r 9a2e0195bb43 stat_presence_abs.xml
--- a/stat_presence_abs.xml Tue Jul 27 16:57:02 2021 +0000
+++ b/stat_presence_abs.xml Fri Aug 13 18:18:21 2021 +0000
b
@@ -17,6 +17,7 @@
             '$species'
             '$location'
             '$time'
+            '$__tool_directory__/functions.r'
             '$output_md'
         ]]>
     </command>
b
diff -r e5552099d0e3 -r 9a2e0195bb43 test-data/SCBD.txt
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/SCBD.txt Fri Aug 13 18:18:21 2021 +0000
b
@@ -0,0 +1,4 @@
+    Blenniidae       Gobiidae   Scorpaenidae Tripterygiidae    Plesiopidae 
+    0.13452280     0.12622864     0.05950360     0.15163110     0.07035376 
+    Apogonidae  Nototheniidae 
+    0.06153987     0.14519360