| Previous changeset 0:c7dd4706f982 (2021-07-27) |
|
Commit message:
"planemo upload for repository https://github.com/Marie59/Data_explo_tools commit 60627aba07951226c8fd6bb3115be4bd118edd4e" |
|
modified:
graph_stat_presence_abs.r |
|
added:
funct_anomy.r functions.r graph_lcbd.r test-data/SCBD.txt |
| b |
| diff -r c7dd4706f982 -r 8e8867bf491a funct_anomy.r --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/funct_anomy.r Fri Aug 13 18:18:00 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 c7dd4706f982 -r 8e8867bf491a functions.r --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/functions.r Fri Aug 13 18:18:00 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 c7dd4706f982 -r 8e8867bf491a graph_lcbd.r --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/graph_lcbd.r Fri Aug 13 18:18:00 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 c7dd4706f982 -r 8e8867bf491a graph_stat_presence_abs.r --- a/graph_stat_presence_abs.r Tue Jul 27 16:55:49 2021 +0000 +++ b/graph_stat_presence_abs.r Fri Aug 13 18:18:00 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 c7dd4706f982 -r 8e8867bf491a test-data/SCBD.txt --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/SCBD.txt Fri Aug 13 18:18:00 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 |