updated output and scripts

RAnalysis/Scripts/F1_Phys/F1_DryWeights_analysis.Rmd

@@ -268,7 +268,7 @@ dev.off()
 * here we do no require a bfactor becuase we are already scalling to the individual!
 
 ### Build poroprtion dataset 
-```{r Poroportion calculate}
+```{r Proportion calculate}
 
 
 Proportion_master <- F1_dryweights_calc.long %>% 
@@ -2109,4 +2109,12 @@ print(ggarrange(Gonad_final_plot, Adductor_final_plot,
 dev.off()
 
 
+pdf("Output/DryWeights/F1/DryWeights_bFactorLengthpCO2.pdf", height=7.5, width =10)
+print(ggarrange(Gonad_final_plot, Adductor_final_plot,
+               Somatic_final_plot, Tissue_final_plot_SUBSET,
+               Shell_final_plot_SUBSET,
+          ncol = 3, nrow = 2
+          ))
+dev.off()
+
 ```

RAnalysis/Scripts/F1_Phys/F1_Lengths_analysis.Rmd

@@ -21,13 +21,14 @@ library(dplyr)
 library(ggplot2)
 library(forcats)
 library(lme4)
-library(lmerTest)
-library(performance)
+# library(lmerTest)
+# library(performance)
+# install.packages('')
 library(car)
 library(ggpubr)
-library(SciViews)
+# library(SciViews)
 library(Rmisc)
-library(hrbrthemes)
+# library(hrbrthemes)lmerTest
 
 # SET WORKING DIRECTORY :::::::::::::::::::::::::::::::::::::::::::::::
 knitr::opts_knit$set(root.dir = "C:/Users/samjg/Documents/Github_repositories/Airradians_multigen_OA/RAnalysis")
@@ -226,8 +227,47 @@ for (i in 1:nrow(ANOVA_Age_DPF)) {
     print(AOVdf_total) # print to monitor progress
 }
 View(AOVdf_total) # view all the anova tests within data 
-View(AOVdf_total %>% dplyr::filter(P_val < 0.05))
+View(AOVdf_total %>% dplyr::filter(P_val < 0.05)) # 243, 275, 318 are significnatly different via one0-way anova
 
 write.csv(AOVdf_total, "C:/Users/samjg/Documents/Github_repositories/Airradians_multigen_OA/RAnalysis/Output/Length/F1/F1_Length_ANOVA_table.csv")
 
+
+# (2)  You can easitly do #1 by piping in the formula here BUT this does not have contingency based on assumptions
+F1_length_MEANS_growout %>%
+  group_by(Age_DPF) %>%
+  t_test(data =., Length_mm ~ pCO2) %>% #243, 275, 318 are significnatly different via one0-way anova
+  adjust_pvalue(method = "bonferroni") %>%
+  add_significance("p.adj") # with bonferroni JUST 273 is significantly different!
+
+
+# (3) GLMM 
+F1_length_MEANS_growout$ID <- paste0(F1_length_MEANS_growout$Replicate,
+                                     '_',
+                                     F1_length_MEANS_growout$pH) # ID variable to run the GLMM withput averaging!
+
+F1_GLMM <- glmer(
+                   Length_mm ~ 
+                    pCO2 * Age_DPF + (1 | ID), 
+                      data = F1_length_MEANS_growout, 
+                     # gaussian(link = "identity")
+                     Gamma(link="identity") 
+                     # Gamma(link="log") 
+                  )
+
+summary(F1_GLMM)
+
+# check the distribution of residuals 
+simulationOutput <- simulateResiduals(fittedModel = F1_GLMM, use.u = T)
+residuals        <- residuals(F1_GLMM, type = "response", retype="normalized")
+
+# pdf(file="Output/1_Survival_Growth/size_growth/Experiment4/Experiment4_Length_GLMM_QQ_Resid.pdf", height = 4, width = 7)
+plot(simulationOutput)
+# dev.off()
+
+# pdf(file="Output/1_Survival_Growth/size_growth/Experiment4/Experiment4_Length_GLMM_Resid.pdf", height = 5, width = 5)
+plot(residuals)
+# dev.off()
+
+# output the anova table 
+car::Anova(F1_GLMM, type=3) # JUST AN EFFECT OF TIME
 ```

RAnalysis/Scripts/F2_Phys/F2_Lengths_analysis.Rmd

@@ -222,4 +222,46 @@ View(AOVdf_total %>% dplyr::filter(P_val < 0.05))
 
 write.csv(AOVdf_total, "C:/Users/samjg/Documents/Github_repositories/Airradians_multigen_OA/RAnalysis/Output/Length/F2/F2_Length_ANOVA_table.csv")
 
+
+
+# (2)  You can easitly do #1 by piping in the formula here BUT this does not have contingency based on assumptions
+F2_length_MEANS_bytank %>%
+  group_by(Age_DPF) %>%
+  t_test(data =., Length_mm ~ pCO2) %>% #243, 275, 318 are significnatly different via one0-way anova
+  adjust_pvalue(method = "bonferroni") %>%
+  add_significance("p.adj") # with bonferroni JUST 273 is significantly different!
+# NO SIGNIFICANT DIFFERENCES WHETHER WITH CONSERVATIVE ADJUSTMENT OR NOT!
+
+# (3) GLMM 
+F2_length_MEANS_bytank$ID <- paste0(F2_length_MEANS_bytank$Replicate,
+                                     '_',
+                                     F2_length_MEANS_bytank$pH) # ID variable to run the GLMM withput averaging!
+
+F2_GLMM <- glmer(
+                   Length_mm ~ 
+                    pCO2 * Age_DPF + (1 | ID), 
+                      data = F2_length_MEANS_bytank, 
+                     # gaussian(link = "identity")
+                     Gamma(link="identity") 
+                     # Gamma(link="log") 
+                  )
+
+summary(F2_GLMM)
+
+# check the distribution of residuals 
+simulationOutput <- simulateResiduals(fittedModel = F2_GLMM, use.u = T)
+residuals        <- residuals(F2_GLMM, type = "response", retype="normalized")
+
+# pdf(file="Output/1_Survival_Growth/size_growth/Experiment4/Experiment4_Length_GLMM_QQ_Resid.pdf", height = 4, width = 7)
+plot(simulationOutput)
+# dev.off()
+
+# pdf(file="Output/1_Survival_Growth/size_growth/Experiment4/Experiment4_Length_GLMM_Resid.pdf", height = 5, width = 5)
+plot(residuals)
+# dev.off()
+
+# output the anova table 
+car::Anova(F2_GLMM, type=3) # JUST AN EFFECT OF TIME
+
+
 ```

RAnalysis/Scripts/Popgen/F1Juveniles_F3Juvneiles_LOWvMOD_workbook.Rmd

@@ -23,7 +23,9 @@ library(pcadapt)
 library(adegenet)
 library(qqman)
 library(ggpubr)
-install.packages("AFLPsim")
+library(radiator) # https://thierrygosselin.github.io/radiator/articles/get_started.html
+# run bayscan in R using radiator https://rdrr.io/github/thierrygosselin/radiator/man/run_bayescan.html
+?run_bayescan
 # library(pegas)
 # library(LDlinkR)
 # library(ldsep) # need a Seqinfo object  Airradians not registered in package
@@ -37,6 +39,8 @@ install.packages("AFLPsim")
 * check out thise use of bayscan for neutral blaancing and outlier loci https://github.com/laurabenestan/Bayescan
 * tutorial to contrast the outliers detected from bayscan and pcadapt https://rpubs.com/lbenestan/outlier
 
+* Narum and Hess 2011 " Comparisons included simulation procedures (FDIST2, ARLEQUIN v.3.5 and BAYESCAN) as well as more conventional tools such as global FST histograms. Of the three simulation methods, FDIST2 and BAYESCAN typically had the lowest type II error, BAYESCAN had the least type I error and Arlequin had highest type I and II error. High error rates in Arlequin with a hierarchical approach were partially because of confounding scenarios where patterns of adaptive variation were contrary to neutral structure; however, Arlequin consistently had highest type I and type II error in all four simulation scenarios tested in this study. Given the results provided here, it is important that outlier loci are interpreted cautiously and error rates of various methods are taken into consideration in studies of adaptive molecular variation, especially when hierarchical structure is included."
+
 ## load files
 
 * vcf.gz files