# HG changeset patch
# User ecology
# Date 1672838513 0
# Node ID ffaa9a71b2d277145154c5f4ef7f1e24d946ca15
# Parent  7e6cc3da11894ebba1f7a2fc28ea8cba8799353e
planemo upload for repository https://github.com/Marie59/champ_blocs commit 0d86db7d42b608c386a54500064f5f9c9d7019a4

diff -r 7e6cc3da1189 -r ffaa9a71b2d2 cb_dissimilarity.xml
--- a/cb_dissimilarity.xml	Fri Dec 02 16:13:18 2022 +0000
+++ b/cb_dissimilarity.xml	Wed Jan 04 13:21:53 2023 +0000
@@ -1,14 +1,14 @@
-<tool id="cb_dissim" name="Dissimilarity" version="@VERSION@" profile = "20.01">
+<tool id="cb_dissim" name="QECB-Dissimilarity" version="@VERSION@" profile = "20.01">
     <description>indicator</description>
     <macros>
         <import>macro.xml</import>
     </macros>
     <expand macro="cb_requirements">
-        <requirement type="package" version="2.16">r-rmarkdown</requirement>
+        <requirement type="package" version="2.19">r-rmarkdown</requirement>
         <requirement type="package" version="0.2.4">r-officedown</requirement>
-        <requirement type="package" version="2.6_2">r-vegan</requirement>
+        <requirement type="package" version="2.6_4">r-vegan</requirement>
         <requirement type="package" version="1.3.2">r-tidyverse</requirement>
-        <requirement type="package" version="3.3.6">r-ggplot2</requirement>
+        <requirement type="package" version="3.4">r-ggplot2</requirement>
     </expand>
     <required_files>
         <include type="literal" path="test-data/champbloc_qecb.csv"/>
@@ -40,9 +40,9 @@
         <param name="choice_date" type="integer" value="2021" label="Until when do you have data (write only the YEAR) ?" help="Example write : 2021 "/>
     </inputs>
     <outputs>
-        <data name="output_qecb" from_work_dir="qecbnato0.RDS" format="Rdata" label="Clean Rdata"/>
-        <data name="output_eco" from_work_dir="bret_egmp_basq_qecb.RDS" format="Rdata" label="Bret Rdata"/>
-        <data name="output_matri" from_work_dir="matri_full.RDS" format="Rdata" label="Matrice full Rdata"/>
+        <data name="output_qecb" from_work_dir="qecbnato0.RDS" format="rdata" label="Clean Rdata"/>
+        <data name="output_eco" from_work_dir="bret_egmp_basq_qecb.RDS" format="rdata" label="Bret Rdata"/>
+        <data name="output_matri" from_work_dir="matri_full.RDS" format="rdata" label="Matrice full Rdata"/>
         <collection type="list" name="plots_bf" label="Bloc fixe et bloc mobile face sup plot">
             <discover_datasets pattern="bm(?P&lt;designation&gt;.+)\.png" visible="false" format="png"/>
         </collection>
diff -r 7e6cc3da1189 -r ffaa9a71b2d2 div_index.r
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/div_index.r	Wed Jan 04 13:21:53 2023 +0000
@@ -0,0 +1,187 @@
+# author: "Jonathan Richir"
+# date: "01 October 2022"
+
+
+#Rscript
+
+###############################
+##  ##
+###############################
+
+#####Packages : dplyr
+#               tidyr
+#               readr
+#               writexl
+#               stringr
+#               readxl
+#               tibble
+#               lubridate
+#               cowplot
+#               magrittr
+#               rmarkdown
+library(magrittr)
+library(dplyr)
+#####Load arguments
+
+args <- commandArgs(trailingOnly = TRUE)
+
+#####Import data
+
+if (length(args) < 1) {
+    stop("This tool needs at least 1 argument")
+}else {
+    qecnato0 <- args[1]
+
+}
+
+qecnato0 <- readRDS(qecnato0)
+
+
+# first, create vector (4) for qecb and fishing by region (same as above)
+
+bret_egmp_basq_qecb <- c(
+  "X..algues.brunes",
+  "X..algues.rouges",
+  "X..algues.vertes",
+  "X..Cladophora",
+  "X..Lithophyllum",
+  "Nb.Littorina.obtusata",
+  "Nb.Gibbula.cineraria",
+  "Nb.Gibbula.pennanti",
+  "Nb.Gibbula.umbilicalis",
+  "Nb.Phallusia.mamillata",
+  "Nb.Tethya.aurantium",
+  "Nb.Spirobranchus.lamarckii.total",
+  "Nb.spirorbis.total",
+  "X..Eponges",
+  "X..Ascidies.Coloniales",
+  "X..Ascidies.Solitaires",
+  "X..Bryozoaires.Dresses",
+  "X..Balanes.Vivantes"
+  #, "X..Recouvrement.Sediment"
+  #, "X..Roche.Nue"
+  #, "X..Surface.Accolement"
+  )
+
+egmp_basq_qecb <- c("Nb.Crassostrea.gigas", "Nb.Ostrea.edulis", "X..Mytilus.sp.", "X..Hermelles", "X..Hydraires")
+
+bret_egmp_basq_fishing <- c("Nb.Cancer.pagurus..Tourteau.",
+                            "Nb.Necora.puber..Etrille.",
+                            "Nb.Carcinus.maenas..Crabe.vert.",
+                            "Nb.Nucella.lapilus..Pourpre.",
+                            "Nb.Galathea..Galathées.",
+                            "Nb.Lophozozymus.incisus..ancien.Xantho.incisus.",
+                            "Nb.Palaemon.sp..Crevette.bouquet.ou.crevette.rose.",
+                            "Nb.Haliotis.tuberculata..Ormeau.",
+                            "Nb.Littorina.littorea..Bigorneau.",
+                            "Nb.Xantho.pilipes..Xanthe.poilu.",
+                            "Nb.Mimachlamys.varia..Pétoncle.noir.")
+
+egmp_basq_fishing <- c("Nb.Eriphia.verrucosa..Crabe.verruqueux.", "Nb.Octopus.vulgaris..Poulpe.", "Nb.Paracentrotus.lividus..Oursin.", "Nb.Stramonita.haemastoma..Pourpre.bouche.de.sang.")
+
+# here I can choose to either replace spirorbis and/or spirobranchus by their log10 transformation in bret_egmp_basq_qecb vector
+bret_egmp_basq_qecb <- replace(bret_egmp_basq_qecb, bret_egmp_basq_qecb == "Nb.spirorbis.total", "log10.Nb.spirorbis.total")
+
+## Diversity index
+
+
+# adiv contains two main functions for species diversity indices: speciesdiv, which includes widely used indices such as species richness and the Shannon index, and divparam, which includes indices that have a parameter to control the importance given to rare versus abundant species in diversity measurements (Pavoine (2020) - adiv: An r package to analyse biodiversity in ecology).
+
+# NB: just like for dissimilarity distance matrices, no sense to use the "fishing" variable lists, because either they are present for the bloc mobile and not for the bloc fixe (therefore false higher diversity for bloc mobile), either they are repeated between face supérieure and face inférieure of bloc mobile.
+
+# function in a loop
+
+row.names(qecnato0) <- c(paste0(qecnato0$region.site_year_month_day, "_", qecnato0$Quadrat.bis, "_", qecnato0$Type.Bloc, "_", qecnato0$Numéro.Bloc.échantillon, "_", qecnato0$Face))
+
+# later on I can copy-paste above code to recreate variable names vector
+#bret_egmp_basq_qecb
+#egmp_basq_qecb
+#Bret_EGMP.BASQ_fishing
+#EGMP.BASQ_fishing
+
+# remove boulder variables
+bret_egmp_basq_qecb <- bret_egmp_basq_qecb[! bret_egmp_basq_qecb %in% c("X..Recouvrement.Sediment", "X..Roche.Nue", "X..Surface.Accolement")]
+
+qecnato0$period <- as.character(qecnato0$period)
+qecnato0$Face <- as.character(qecnato0$Face)
+
+div_list <- vector("list", length(unique(qecnato0$site_year_month_day)))
+
+for (i in c(1:nrow(qecnato0))) {
+  div_i <- dplyr::filter(qecnato0, site_year_month_day == qecnato0$site_year_month_day[i])
+
+  ifelse(unique(div_i$region) == "Bretagne", var. <- c(bret_egmp_basq_qecb), var. <- c(bret_egmp_basq_qecb, egmp_basq_qecb)) # Qu. : Why can't R's ifelse statements return vectors? => you can circumvent the problem if you assign the result inside the ifelse.
+
+  #8 remove empty row cfr: In speciesdiv(div_i[, var.]) & divparam(div_i[, var.]) : empty communities should be discarded
+  div_i <- dplyr::filter(div_i, rowSums(div_i[, var.]) > 0)
+
+  div_i_speciesdiv <- adiv::speciesdiv(div_i[, var.])
+  adiv_i_df <- data.frame(div_i_speciesdiv)
+  
+  div_i_divparam <- adiv::divparam(div_i[, var.], q = c(0, 0.25, 0.5, 1, 2, 4, 8)) # When q increases, abundant species are overweighted compared to rare species, we thus expect that the evenness in species weights decreases.
+
+
+  par(mfrow = c (1, 1))
+  plot(adiv::divparam(div_i[, var.], q = 0), main = unique(div_i$site_year_month_day))
+  plot(adiv::divparam(div_i[, var.], q = 0:10), legend = FALSE, main = unique(div_i$site_year_month_day))
+  
+  adiv_i_df$x <- div_i_divparam$div$`1`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[1], " (equi. richness)")
+  adiv_i_df$x <- div_i_divparam$div$`2`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[2])
+  adiv_i_df$x <- div_i_divparam$div$`3`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[3])
+  adiv_i_df$x <- div_i_divparam$div$`4`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[4])
+  adiv_i_df$x <- div_i_divparam$div$`5`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[5], " (equi. Simpson)")
+  adiv_i_df$x <- div_i_divparam$div$`6`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[6])
+  adiv_i_df$x <- div_i_divparam$div$`7`
+  colnames(adiv_i_df)[which(colnames(adiv_i_df) == "x")] <- paste0("Para. ISD, q = ", div_i_divparam$q[7])
+  
+  # plot
+  par(mfrow = c(3, 2))
+  sapply(names(adiv_i_df[, c(1, 8, 2, 3, 12, 4:7, 9:11, 13:ncol(adiv_i_df))]), 
+         function(cname) {
+           png(paste0(cname, "_histo.png"))
+           hist(adiv_i_df[, c(1, 8, 2, 3, 12, 4:7, 9:11, 13:ncol(adiv_i_df))][[cname]], main = "", xlab = cname, breaks = length(unique(adiv_i_df[, c(1, 8, 2, 3, 12, 4:7, 9:11, 13:ncol(adiv_i_df))][[cname]])))
+dev.off()
+         }
+
+)
+  par(mfrow = c(1,1))
+  
+  div_list[[i]] <- adiv_i_df
+  
+  rm(div_i, adiv_i_df, div_i_speciesdiv, div_i_divparam)
+  
+}
+
+# for the error message due to richness NA data => 35 observations in 21 surveys; no reason to remove these data "...remo", was checked in the complete script.
+
+
+div_df <- do.call("rbind", div_list)
+
+
+# There is an issue with region.terri that are merged with no "_" ...
+
+div_df <- tibble::add_column(div_df, rownames. = rownames(div_df), .before = "richness")
+div_df <- tidyr::separate(div_df, rownames., into = c("region.terri.", "site_year_month_day", "Quadrat.bis", "Type.Bloc", "Numéro.Bloc.échantillon", "Face"), sep = "_")
+
+# I therefore add these lines to solve that issue
+
+div_df <- tibble::add_column(div_df, terri. = substring(div_df$region.terri., nchar(div_df$region.terri.)-3), .after = "region.terri.")
+
+div_df$region.terri. <- substring(div_df$region.terri., 1, nchar(div_df$region.terri)-4)
+div_df <- dplyr::rename(div_df, region = region.terri.)
+
+div_df$site_year_month_day <- paste0(div_df$terri., "_", div_df$site_year_month_day)
+div_df <- subset(div_df, select = -c(terri.))
+
+div_df$Type.Bloc <- as.factor(div_df$Type.Bloc)
+div_df$Face <- as.factor(div_df$Face)
+div_df$Numéro.Bloc.échantillon <- as.integer(div_df$Numéro.Bloc.échantillon)
+
+saveRDS(div_df, "div_df.RDS")
+write.table(div_df, "Valeurs_stat.tabular", row.names = FALSE, quote = FALSE, sep = "\t", dec = ".", fileEncoding = "UTF-8")
diff -r 7e6cc3da1189 -r ffaa9a71b2d2 macro.xml
--- a/macro.xml	Fri Dec 02 16:13:18 2022 +0000
+++ b/macro.xml	Wed Jan 04 13:21:53 2023 +0000
@@ -2,11 +2,12 @@
     <token name="@VERSION@">0.0.0</token>
     <xml name="cb_requirements">
         <requirements>
+            <requirement type="package" version="4.2.2">r-base</requirement>
             <requirement type="package" version="1.0.10">r-dplyr</requirement>
             <requirement type="package" version="1.2.1">r-tidyr</requirement>
-            <requirement type="package" version="1.4.1">r-stringr</requirement>
+            <requirement type="package" version="1.5.0">r-stringr</requirement>
             <requirement type="package" version="3.1.8">r-tibble</requirement>
-            <requirement type="package" version="1.40">r-knitr</requirement>
+            <requirement type="package" version="1.41">r-knitr</requirement>
             <requirement type="package" version="2.0.3">r-magrittr</requirement>
            <yield/>
         </requirements>
diff -r 7e6cc3da1189 -r ffaa9a71b2d2 test-data/Clean_Rdata.rdata
Binary file test-data/Clean_Rdata.rdata has changed