diff --git a/2-Nutrient_Flux_Bruch.R b/2-Nutrient_Flux_Bruch.R new file mode 100644 index 0000000000000000000000000000000000000000..61f90a1903d05b443d5ffe1de45452ecef130238 --- /dev/null +++ b/2-Nutrient_Flux_Bruch.R @@ -0,0 +1,400 @@ +rm(list=ls()) +setwd("C:/users/camille.poulet/Work Folders/Documents/Thèse/Script/Nutrients_Fluxes") + +###### Attention, perte de masse = 2* perte de gonade du à energie. +###### Donc on ne peut pas approximer le poids avant et après comme étant la réduction de la gonade, mais comme étant la réduction de la gonade + perte enrgétique du soma ! + +############## Packages ############### + +library(tidyverse) +library(devtools) +library(ggpmisc) + +############# Préparation & Nettoyage du jeu de données ######## + +BDalosesBruch<-read.csv("BDalosesBruch.csv", header= T, sep=";", dec=",") %>% + rename(Season = Année, + River = Rivière.piégeage, + Sex = Sexe, + LT = Lt..cm., + LF = Lf..cm. , + WT = M.tot..g., + WG = M.gonades..g., + Fish_ID = N..poisson) %>% + na.omit(.) %>% + filter (River != "Dordogne/Garonne") # On enlève les sites d'échantillonage indéterminés + +BDalosesBruch$Sex = str_to_upper(BDalosesBruch$Sex) #Homogénéistaion des caractéres sur la variable "sex": tout est noté en MAJ + +#On isole les données pré-reproduction + +BDalosesBruchPreSpawn<- BDalosesBruch %>% + filter (LOT == "Golfech" | LOT=="Tuilières") + +#On les enlèves de la BDD initial pour n'avoir que les données pots-reproduction +BDalosesBruch<- BDalosesBruch %>% + filter (LOT != "Golfech" & LOT != "Tuilières") + + +######## Calcul des moyennes WT, WG, LT, LF ######### + +#============ Données Post reproduction =========== +MeanFeatures <- BDalosesBruch %>% + select("Sex","Season","WT","WG","LF","LT") %>% + group_by(Sex) %>% + #mutate(LF_Quantile =quantile(.$LF, prob = 0.5)) %>% + summarise_each(.,funs(mean)) + +#============ Données Pré-reproduction ============ +MeanFeaturesPreSpawn <- BDalosesBruchPreSpawn %>% + select("Sex","Season","WT","WG","LF","LT") %>% + group_by(Sex) %>% + summarise_each(.,funs(mean)) + +######## LF ######### + +tapply(BDalosesBruch$LF, BDalosesBruch[,c("Sex")],min, na.rm = TRUE) +tapply(BDalosesBruch$LF, BDalosesBruch[,c("Sex")],max, na.rm = TRUE) + + +tapply(BDalosesBruch$LF, BDalosesBruch[,c("Season","Sex")],quantile, na.rm = TRUE, probs= 0.05) +tapply(BDalosesBruch$LF, BDalosesBruch[,c("Season","Sex")],quantile, na.rm = TRUE, probs= 0.5) +tapply(BDalosesBruch$LF, BDalosesBruch[,c("Season","Sex")],quantile, na.rm = TRUE, probs= 0.95) + +######## WG ######### + +#=========== Post reproduction ================ + +tapply(BDalosesBruch$WG, BDalosesBruch[,c("Season","Sex")],min, na.rm = TRUE) +tapply(BDalosesBruch$WG, BDalosesBruch[,c("Season","Sex")],max, na.rm = TRUE) + +#=========== Pré reproduction ================ + +tapply(BDalosesBruchPreSpawn$WG,BDalosesBruchPreSpawn [,c("Season","Sex")],min, na.rm = TRUE) +tapply(BDalosesBruchPreSpawn$WG, BDalosesBruchPreSpawn[,c("Season","Sex")],max, na.rm = TRUE) + +########### Relation taille-poids ############### + +#========== Vérification des Wpre avec les données de Taverny + +#Taverny (p44 livre Baglinière & Elie), 1991 + +#WT = aLT^b + +#Female +#a=1.2102e-6 +#b=3.3429 + +#Male +#a= 2.4386e-6 +#b= 3.2252 + +WeigtSizeRelationshipTaverny = function (a,Lt,b){ + + a * Lt^b +} + +BDalosesBruchPreSpawn<-BDalosesBruchPreSpawn %>% + select("Sex","Season","WT","WG","LT") %>% + group_by(Sex,Season) %>% + mutate(LTmm= LT*10) %>% + mutate(WTTaverny= case_when( + Sex=="F" ~ WeigtSizeRelationshipTaverny(1.2102e-6,LTmm,3.3429), + Sex=="M" ~ WeigtSizeRelationshipTaverny(2.4386e-6,LTmm,3.2252))) %>% + mutate(DiffWeight= abs(WT-WTTaverny)) + +#======= Graphique relation taille-poids ======== + +#TabBdMix<-gather(data = BDalosesBruchPreSpawn,key = Weight,value= Value,WTTaverny,WT,factor_key = TRUE) + +WT.Taverny<-BDalosesBruchPreSpawn %>% + select("WTTaverny") +ggplot(TabBdMix, aes (LT,Value))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=paste(stat(eq.label),stat(adj.rr.label), sep="--")), formula = y~x, coef.digits = 5, rr.digits =2)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + facet_grid(Sex~.)+ + ggtitle("Pre spawning") + + +#Pre reproduction +ggplot(BDalosesBruchPreSpawn, aes(log(LTmm),log(WTTaverny)))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=paste(stat(eq.label),stat(adj.rr.label), sep="--")), formula = y~x, coef.digits = 5, rr.digits =2)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + facet_grid(Sex~.)+ + ggtitle("Pre spawning") + +#Post Reproduction +ggplot(BDalosesBruch, aes(LF,WT))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + #stat_poly_eq(parse=TRUE, aes(label=..eq.label..), formula = y~x, coef.digits = 5)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + facet_grid(Sex~.)+ + ggtitle("Post spawning") + + +#========= Vérification de la relation entre LT et LF ========# + +# LT = a+bFL (from Gaygusuz et al 2006) + + +ggplot(BDalosesBruchPreSpawn, aes(LF,LT))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=..eq.label..), formula = y~x, coef.digits = 5)+ + xlab("LF (cm)")+ + ylab("LT(cm)")+ + theme_bw() + + + +SSRelationship = function (a,b,Lf){ + + a + b*Lf +} + +BDalosesBruch<-BDalosesBruch %>% + select("Sex","Season","WT","WG","LF","LT") %>% + group_by(Sex,Season) %>% + mutate(LT.SSRelationship= SSRelationship(6.9324,0.9611, LF)) %>% + mutate(DiffWeight=(LT-LT.SSRelationship)) + +BDalosesBruch<-BDalosesBruch %>% + select("Sex","Season","WT","WG","LF","LT") %>% + group_by(Sex,Season) %>% + mutate(LT.SSRelationship= case_when( + Sex=="F" ~ SSRelationship(exp(0.64858),0.85885, LF), + Sex=="M" ~ SSRelationship(exp(0.59854), 0.87097,LF))) %>% + mutate(DiffWeight=(LT-LT.SSRelationship)) + + +#Vérification sur les données de Bruch + + + + + + + +########## Calcul des GSI ################ + +# GSI = (WG / (WT-WG))*100 + +#========== Post reproduction ============== + +DataGSI <- BDalosesBruch %>% + select("Sex","WT","WG") %>% + group_by("Sex") %>% + #summarise_each(.,funs(mean)) %>% + mutate(GSI = WG/(WT)) + +######## Calcul des GSI based on Taverny WSR ##### + +BDalosesBruchPreSpawn %>% + select("Sex","WG","DiffWeight") %>% + group_by(Sex) %>% + mutate(GSI = WG/(DiffWeight)) %>% + summarise(., avg=mean(GSI)) + + +#########" TEST trouver l'erreur ######### + +hist(DataGSI$GSI[DataGSI$Sex=="F"], breaks=200, col="red") +abline(v=0.15) +abline(v=0.12, col="green") + +DataGSI <- BDalosesBruch %>% + select("Sex","WT","WG", "LT", "LF") %>% + group_by("Sex") %>% + filter(Sex=="F") %>% + #summarise_each(.,funs(mean)) %>% + mutate(GSI = WG/(WT)) %>% + filter(GSI<0.12) + +ggplot(DataGSI, aes(log(LF),log(WT)))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=paste(stat(eq.label),stat(adj.rr.label), sep="--")), formula = y~x, coef.digits = 5, rr.digits =2)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + ggtitle("Post spawning") + + +#========== Pré reproduction ============== + +DataGSIPreSpawn <- BDalosesBruchPreSpawn %>% + select("Sex","WT","WG") %>% + group_by(Sex) %>% + summarise_each(.,funs(mean)) %>% + mutate(GSI = WG/(WT)) + +#=========== Calcul d'une perte de poids total ========== + +WTpre<-tapply(BDalosesBruchPreSpawn$WT, BDalosesBruchPreSpawn[,"Sex"],mean, na.rm = TRUE) +WTpost<-tapply(BDalosesBruch$WT, BDalosesBruch[,"Sex"],mean, na.rm = TRUE) +WTloss <- (WTpre - WTpost) + +#=========== Calcul d'une perte de poids gonade ========== + +#Barber et al, 2018: + +#Wgpré-Wgpost = mass of spawned egg/sperm + +WGpre<-tapply(BDalosesBruchPreSpawn$WG, BDalosesBruchPreSpawn[,"Sex"],mean, na.rm = TRUE) +WGpost<-tapply(BDalosesBruch$WG, BDalosesBruch[,"Sex"],mean, na.rm = TRUE) +WGloss <- (WGpre - WGpost) + +#========== Calcul de WT post-reproduction par la relation allométrie ======== + +# WT = a * Lt ^b /( GSI +1) +# WT = Wpre/(GSI +1) + +#==== Based on Taverny weight-size relationship + +ComputeWeigtPostSpawning = function (a, Lt, b, GSI){ + + Wpre <- WeigtSizeRelationshipTaverny (a,Lt,b) + + return (Wpost = Wpre/(GSI+1)) + +} + +BDalosesBruch<-BDalosesBruch%>% + select("Sex","Season","WT","WG","LF","LT") %>% + group_by(Sex,Season) %>% + mutate(LT.mm= LT*10) %>% + mutate(WT.CPL= case_when( + Sex=="F" ~ ComputeWeigtPostSpawning(1.2102e-6,LT.mm,3.3429,0.10389710), + Sex=="M" ~ ComputeWeigtPostSpawning(2.4386e-6,LT.mm,3.2252,0.05026625))) %>% + mutate(DiffWeight=(WT-WT.CPL)) + +ggplot(BDalosesBruch, aes(WT.CPL,WT))+ + geom_point(stat="identity")+ + xlab("Weight size relationship (Taverny)")+ + ylab("Weight (g) ")+ + theme_bw()+ + facet_grid(Sex~.) + + +BDalosesBruch %>% + filter (LT > 65) + +########### Calcul d'une relation d'allométrie WG et LT ##### + +#Pre spawning +ggplot(BDalosesBruchPreSpawn, aes(log(WT),log(WG)))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=paste(stat(eq.label),stat(adj.rr.label), sep="--")), formula = y~x, coef.digits = 5, rr.digits =2)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + facet_grid(Sex~.)+ + ggtitle("Pre spawning") + +#Post spawning +ggplot(BDalosesBruch, aes(log(WT),log(WG)))+ + geom_point(stat="identity")+ + geom_smooth(method="lm")+ + stat_poly_eq(parse=TRUE, aes(label=paste(stat(eq.label),stat(adj.rr.label), sep="--")), formula = y~x, coef.digits = 5, rr.digits =2)+ + xlab("Size (cm)")+ + ylab("Weight (g)")+ + theme_bw()+ + facet_grid(Sex~.)+ + ggtitle("Post spawning") + +BDalosesBruch<-BDalosesBruch %>% + filter (log(LF)<4.2) + + + +#========= Calcul des Wpost et comparaison avec la BDD ====== + +#totalWeightPost = totalWeightPre * (1-GSIfemalePost)+ totalWeightPost * GSIfemalePost * CoeffLossWeight + +GSI <- DataGSI$GSI[,DataGSI$Sex=="F"] +WTpost = mean (BDalosesBruchPreSpawn$WT)*(1-DataGSI$GSI[,sex=="F"])+ mean(BD) + +############ Boxplot exploratoires ########## + +par(mfrow=c(1,3)) +boxplot(BDalosesBruch$WG ~ BDalosesBruch$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean Gonad Mass (g)",main="Gonad Mass") +boxplot(BDalosesBruch$WT ~ BDalosesBruch$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean Total Mass (g)", main = "Total Mass") +boxplot(BDalosesBruch$LF ~ BDalosesBruch$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean LF (cm)", main = "LF") + + +par(mfrow=c(1,3)) +boxplot(MeanFeatures$WG ~ MeanFeatures$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean Gonad Mass (g)",main="Gonad Mass") +boxplot(MeanFeatures$WT ~ MeanFeatures$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean Total Mass (g)", main = "Total Mass") +boxplot(MeanFeatures$LF ~ MeanFeatures$Sex, names=c("F","M"), xlab = "Sex", ylab= "Mean LF (cm)", main = "LF") + + +############## Mean Gonad Mass ############# +ggplot(MeanFeatures, aes(Season,(WG)))+ + geom_bar(stat="identity")+ + xlab("Season")+ + ylab("Gonad Mass (g)")+ + theme_bw()+ + facet_grid(Sex~.) + +############## Mean Total Mass ############# +ggplot(MeanFeatures, aes(Season,(WT)))+ + geom_bar(stat="identity")+ + xlab("Season")+ + ylab("Total Mass (g)")+ + theme_bw()+ + facet_grid(Sex~.) + +############## Mean Total Lenght ############# +ggplot(MeanFeatures, aes(Season,(LF)))+ + geom_bar(stat="identity")+ + xlab("Season")+ + ylab("Total length (cm)")+ + theme_bw()+ + facet_grid(Sex~.) + +############## GSI ############# + +ggplot(BDalosesBruch, aes(Season, (WT-WG)/(WT))) +geom_bar(stat="identity")+ + xlab("Season")+ + ylab("Total length (cm)")+ + theme_bw()+ + facet_grid(Sex~.) + + +#TODO + +#Faire un graphique avec les GSI pour regarder si on obtient des données similaires à celles trouvées dans la littérature +#Regarder les relation taille-poids pour définir des coefficient a et b +#Faire des quantiles 5% sur la taille et le poids pour enlever les données abérentes +#Comparer les valeurs de poids moyen des gonades avec ceux de la littérature + +################ StringR essay ################# + +BDalosesBruch$Fish_ID<-as.character(BDalosesBruch$Fish_ID) + +## Detect function +Pre_Spawn<-str_detect(BDalosesBruch$Fish_ID, regex("Golfech.*$", ignore_case = TRUE)) +Pre_Spawn<-str_detect(BDalosesBruch$Fish_ID, regex("Tuilières.*$", ignore_case = TRUE)) + +str_count(Pre_Spawn, "TRUE") + +## Subset function +Pre_Spawn<-str_subset(BDalosesBruch$Fish_ID, "Golfech.*$") +Pre_Spawn<-str_subset(BDalosesBruch$Fish_ID, regex ("Tuilières.*$", ignore_case = TRUE)) + +str_count(Pre_Spawn, "TRUE") + +