Early evidence suggests that DNA methylation can mediate phenotypic responses of marine calcifying species to ocean acidification (OA). Few studies, however, have explicitly studied DNA methylation in calcifying tissues through time. Here, we examined the phenotypic and molecular responses in the extrapallial fluid and mantle (fluid and tissue at the calcification site) in adult eastern oyster (Crassostrea virginica) exposed to experimental OA over 80 days. Oysters were reared under three experimental pCO2 treatments (‘control’, 580 μatm; ‘moderate OA’, 1000 μatm; ‘high OA’, 2800 μatm) and sampled at 6 time points (24 hours - 80 days). We found that high OA initially induced an increase in the pH of the extrapallial fluid (pHEPF) relative to the external seawater that peaked at day 9, but then diminished over time. Calcification rates were significantly lower in the high OA treatment compared to the other treatments. To explore how oysters regulate their extrapallial fluid, gene expression and DNA methylation were examined in the mantle-edge tissue of oysters from days 9 and 80 in the control and high OA treatments. Mantle tissue mounted a significant global molecular response (both in the transcriptome and methylome) to OA that shifted through time. Although we did not find individual genes that were significantly differentially expressed under OA, the pHEPF was significantly correlated with the eigengene expression of several co-expressed gene clusters. A small number of OA-induced differentially methylated loci were discovered, which corresponded with a weak association between OA-induced changes in genome-wide gene body DNA methylation and gene expression. Gene body methylation, however, was not significantly correlated with the eigengene expression of pHEPF-correlated gene clusters. These results suggest that OA induces a subtle response in a large number of genes in C. virginica, but also indicate that plasticity at the molecular level may be limited. Our study highlights the need to reassess our understanding of tissue-specific molecular responses in marine calcifiers, as well as the role of DNA methylation and gene expression in mediating physiological and biomineralization responses to OA.

Collection and Experimental Conditions

Adult C. virginica were collected from three intertidal sites within Plum Island Sound, Massachusetts, USA (Site 1, 42.751636, -70.837023; Site 2, 42.725186, -70.855022; Site 3, 42.681764, -70.813498) in late April 2017. These sites are all within 8 km of each other and the oysters from these sites are not distinct genetic populations. 

In the lab, oysters were cleaned of epibionts and shell ports were installed over eight days while being maintained in 50-L flow-through tanks. To measure the carbonate chemistry of the EPF, a 2 mm hole was drilled into the right valve approximately 2 cm from the hinge to expose the EPF cavity, without damaging the underlying mantle tissue. The drilled hole was gently rinsed with filtered seawater and patted dry. The barbed end of a nylon luer-lock coupling (McMaster-Carr 51525K123) was trimmed to the length of the shell thickness, inserted into the hole and sealed in place using marine-safe cyanoacrylate (Starbond EM-2000 CA USA). The coupling was then plugged with a nylon cap socket (McMaster-Carr 51525K315) to prevent exchange between seawater and the EPF via the port and limit the potential for increased dissolution near the site of the implant. Oysters were allowed to recover for 4 days, after which they were transferred to the experimental tank array (described below) for acclimation and exposure. During the transfer, oysters were assigned one of three pCO2 treatments: ‘control’ (ca. 580 μatm) – corresponding to near-present-day conditions in the estuary where they were collected; ‘moderate OA’ (ca. 1000 μatm) – corresponding to end-of-century predictions (IPCC 2007); and ‘high OA’ (ca. 2800 μatm) – corresponding to seawater conditions that are undersaturated with respect to calcite (Ωcalcite < 1). 

Target pCO2 gases were formulated by mixing compressed CO2 with compressed CO2-free air (low-OA treatment) or with compressed air (moderate- and high-OA treatments) using solenoid-valve-controlled mass flow controllers (Aalborg mass flow controllers, Model GFC17, precision = 0.1mL/min) at flow-rates proportional to the target pCO2 conditions. The CO2-free air used in the low-OA treatment was generated by scrubbing CO2 from compressed air with a Parker Balston FT-IR Purge Gas Generator. The mixed gases were then bubbled into the experimental seawater treatments with 45-cm flexible microporous tubes at rates sufficient for the pCO2 of aquaria seawater to equilibrate with the pCO2 of the mixed gases. Filtered seawater was introduced to each aquarium at a flow rate of 150 mL min-1. Temperature of all experimental tanks was maintained at 17°C with 1/4HP chillers (Aqua Euro USA, precision = 0.1°C). Seawater temperature of all aquaria was slowly increased to 19.5°C (Figure S1.2) between days 39-51 of the experiment in an effort to stimulate gonad development for a companion experiment. Oysters were fed 1% Shellfish Diet 1800® twice daily following best practices outlined in Helm and Bourne (2004). 

Temperature, pH, and salinity of all tanks were measured three times per week (M,W, F) for the duration of the experiment. Seawater pH was measured with an Accumet solid state pH electrode (precision = 1mV) and adjusted to the total scale, salinity was measured using a YSI 3200 conductivity probe (precision = 0.1 ppt), and temperature was measured using a NIST-standardized glass thermometer (precision = 0.1 °C). To empirically correct for the liquid junction bias of the pH electrode, we obtained the slope of the calibration using NBS buffers (pH 7.01 and pH 10.01) and adjusted the intercept using Dickson seawater certified reference material. Seawater samples were collected every two weeks from each tank for analysis of dissolved inorganic carbon (DIC) and total alkalinity (TA) on a VINDTA 3C coupled alkalinity gram titration and coulometric DIC analyzer system. In brief, samples were collected in 250 ml borosilicate glass bottles and immediately poisoned with 100 ul saturated HgCl2 solution, then refrigerated until analyzed. DIC, TA, salinity, and temperature were used to calculate calcite saturation state, pH, CO32-, HCO3-, aqueous CO2, and pCO2 of each sample using CO2SYS version 2.1 (Pierrot, Wallace, and Lewis 2011), using the seawater pH scale with K1 and K2 values from (Roy et al., 1993), a KHSO4 value from (Dickson 1990), and a [B]T value from (Lee et al., 2010).

Calcification Rate

Net calcification rate was calculated for oysters surviving to either 50 or 80 days (n = 35) by buoyantly weighing oysters prior to exposure (BW1) and on day 33 or 34 of the exposure (BW2) following the methods of Ries et al. (2009). Buoyant weight was measured in a 27.65 liter tank (48 cm long, 24 cm wide and 24 cm deep) filled with seawater from the flow-through system and was maintained at treatment temperature by an Aqua Euro USA Model MC-1/4HP aquarium chiller. Buoyant weight was measured by completely submerging the oyster on a flat platform suspended from a bottom-loading scale (Cole Parmer Symmetry S-PT 413E, precision = 0.001 g). Care was taken to ensure no bubbles were trapped inside oyster shells.  Oysters were weighed three times in the weighing basket, removing the oyster from the basket between each measurement. If replicate measurements varied by more than 0.01 g, oysters were re-weighed. A standard of known weight was weighed every 20 oysters to ensure that no drift was occurring in the scale.

            To establish a relationship between buoyant weight and dry weight for the purpose of estimating net calcification rate, shells of oysters sampled for tissue within four days of a buoyant weight measurement were soaked in 10% ethanol to remove salts, dried, and weighed. The dry weight was then regressed against the buoyant weight measurement to establish an empirical dry-buoyant weight relationship (Figure S2.1):

DryWgtBWi(mg)=1.87BWi−2.74

This empirical relationship was then used to calculate dry shell weight at each buoyant weight time point via linear regression. 

Calculated dry weights were then used to calculate daily calcification rate:

CalcificationRate(%) = (DryWgtBW2-DryWgtBW1)n* 100DryWgtBW1;

where dry weight (DryWgtBWi) was calculated for each individual using the buoyant weight pre-exposure (BW1) and 33-34 days into the exposure (BW2) and n was the number of days between the two measurements. Lastly, the average daily change in dry weight was divided by initial dry weight to standardize calcification rate for allometric effects, and multiplied by 100 to convert that fraction into a percent. 

Extrapallial fluid chemistry 

Oyster pHEPF was measured by removing each oyster from its tank, inserting a 5 mL syringe with a flexible 18-gauge polypropylene tip through the luer-lock port into the oyster’s extrapallial cavity, and extracting approximately 0.5-2 mL of extrapallial fluid. Care was taken to avoid puncturing the mantle tissue and inadvertently sampling either the hemolymph or stomach fluid. The pHEPF was measured immediately after extraction with an Orion 91’10DJWP Double Junction micro-pH probe calibrated with pH 7.01 and 10.01 NBS buffers (for slope) and Dickson seawater Certified Reference Material (for intercept).

Tissue collection 

Patterns observed in pHEPF over time informed the decision to focus molecular analyses on mantle edge tissue extracted on days 9 and 80 of the exposure in the control and high OA treatments (ntissue = 6 per treatment per time point). For consistency and to minimize potential effects on mantle tissue arising from drilling of the oysters’ right valves when installing the EPF extraction ports, mantle tissue was only collected from the left valve of the oysters (i.e., the side opposite of the EPF extraction port). Immediately following extraction of EPF at each time point, the oysters were shucked and the edge of the mantle tissue was sampled. Tissue was immediately flash frozen in liquid nitrogen before being transferred to a -80°C freezer for storage prior to DNA/RNA extraction.

DNA methylation processing and quantification

DNA was isolated from the oyster mantle edge tissue samples using the E.Z.N.A. Mollusc Kit (Omega) pursuant to the manufacturer’s instructions. Isolated DNA was quantified, sonicated, and fragmentation was verified using a 2200 TapeStation System (Agilent Technologies, Santa Clara, CA, USA). Samples were enriched for methylated DNA using MethylMiner kit (Invitrogen, CA, USA), then sent to Zymo Research (Zymo,CA,USA) for bisulfite-conversion, library preparation, and sequencing. Libraries were sequenced on an Illumina HiSeq1500 platform with paired-end 100bp sequencing. 

Raw sequences were trimmed to remove adapters and low quality sequences by removing 10 bp from both 5’ and 3’ ends of each fragment using the command trim_galore in the program TrimGalore! (Martin 2011) with the --clip_r1, --clip_r2, --three_prime_clip_R1, and --three_prime_clip_R2 flags set to 10, and the remaining flags left as the defaults. Quality of sequences was assessed with FastQC (Leggett et al. 2013) within the trim_galore command using the --fastqc_args flag. Next, a bisulfite converted version of the C. virginica genome (NCBI Accession GCA_002022765.4) was prepared with the bismark_genome_preparation command in Bismark (Krueger and Andrews 2011) and the --bowtie2 flag. The trimmed sequences were then aligned to the prepared reference using the bismark command with --non-directionality specified and the alignment score set to -score_min L,0,-0.8. Duplicates were removed from the data using the command deduplicate_bismark. All cytosines in the genome with coverage were extracted from the aligned deduplicated data for each individual using the command bismark_methylation_extractor with default settings. Next, a coverage report for all CpG associated cytosines within the genome was generated for each individual with the command coverage2cytosine using the output from the previous step. This report contains all cytosine loci (even those with no coverage) located within a CpG motif and includes separate columns for methylated vs. unmethylated coverage. Files were then processed using a standard methylKit pipeline, which included normalization of methylation counts among samples, destranding CpGs loci, and filtering loci that did not have at least 5x coverage for each sample (Akalin et al. 2012). Here, destranding involves combining the cytosine calls within a single CpG from either strand (i.e., top or bottom) and the proportion methylation was determined by dividing the number of methylated cytosine calls by the total coverage (regardless of methylation status). Before proceeding to downstream analyses samples were first visualized in a PCA to identify and remove outlier samples.

Gene expression processing and quantification 

RNA was extracted from the same mantle tissues used for DNA methylation using TRI Reagent Protocol (Applied Biosystems) following manufacturer’s instructions. A cDNA library was constructed for each sample using a TruSeq Stranded mRNA kit (Illumina, Cat# RS-122-2101) following manufacturer’s instructions. Library fragment size was determined using a BioAnalyzer (RNA 6000 Nano Kit; Agilent, CA, USA) and concentration was quantified on a Qubit (Qubit High Sensitivity RNA Kit; Invitrogen,CA,USA). Samples were sent to Genewiz (NJ ,USA) for library preparation and sequencing. Libraries were sequenced on two lanes (12 individuals/lane) of an Illumina HiSeq 3000 platform with paired-end 200bp sequencing. 

RNA reads were first trimmed to remove adapters and for quality control using Trimmomatic (Bolger, Lohse, and Usadel 2014), implemented in the dDocent pipeline (v.2.2.20; Puritz, Hollenbeck, and Gold 2014) and following recommendations by Puritz et al. (2014). This was done using the dDocent command with default settings, which included removing adapters and performing a quality control step that trimmed the leading and trailing bases with phred quality scores < 20, with additional trimming of 5bp windows with mean phred quality scores < 10 (see github repository for complete description). Reads were then aligned to the C. virginica genome (Accession: GCA_002022765.4) using a two-step mapping approach with the program STAR (v.2.7.0; Dobin et al. 2013). In the first step, a preliminary alignment to the reference genome was performed for each sample to identify novel splice junctions. In the second step, reads were realigned, but with all the splice junctions discovered from step one included along with the reference genome. Both steps were executed with the command STAR using default settings, with the exception that the number of alignments retained was adjusted by setting both the --outFilterMatchNminOverLread and --outFilterScoreMinOverLread flags to 0.17, and increasing the number of discoverable sequence junctions by setting the flag --limitSjdbInsertNsj to 1500000. In the second step, the --sjdbFileChrStartEnd flag was also included to specify the inclusion of splice junctions identified in step one. Expression levels for each gene were quantified for each sample with the program RSEM (Li and Dewey 2011) using the rsem-calculate-expression command with --star flag, default settings, the bam files generated from the mapping step with STAR, and the same C. virginica gene annotations used during mapping (Accession: GCA_002022765.4). 

References

Akalin, Altuna, Matthias Kormaksson, Sheng Li, Francine E. Garrett-Bakelman, Maria E. Figueroa, Ari Melnick, and Christopher E. Mason. 2012. “methylKit: A Comprehensive R Package for the Analysis of Genome-Wide DNA Methylation Profiles.” Genome Biology 13 (10): R87.

Bolger, Anthony M., Marc Lohse, and Bjoern Usadel. 2014. “Trimmomatic: A Flexible Trimmer for Illumina Sequence Data.” Bioinformatics  30 (15): 2114–20.

Dickson, Andrew G. 1990. “Standard Potential of the Reaction : AgCl(s) + 12H2(g) = Ag(s) + Hcl(aq), and and the Standard Acidity Constant of the Ion HSO4− in Synthetic Sea Water from 273.15 to 318.15 K.” The Journal of Chemical Thermodynamics 22 (2): 113–27.

Dobin, Alexander, Carrie A. Davis, Felix Schlesinger, Jorg Drenkow, Chris Zaleski, Sonali Jha, Philippe Batut, Mark Chaisson, and Thomas R. Gingeras. 2013. “STAR: Ultrafast Universal RNA-Seq Aligner.” Bioinformatics  29 (1): 15–21.

Helm, Michael M., and Neil Bourne. 2004. Hatchery Culture of Bivalves: A Practical Manual. Food & Agriculture Org.

Intergovernmental Panel on Climate Change. 2007. Climate Change 2007 - Mitigation of Climate Change: Working Group III Contribution to the Fourth Assessment Report of the IPCC. Cambridge University Press.

Krueger, Felix, and Simon R. Andrews. 2011. “Bismark: A Flexible Aligner and Methylation Caller for Bisulfite-Seq Applications.” Bioinformatics  27 (11): 1571–72.

Lee, Kitach, Kim Tae-Wook, Robert H. Byrne, Frank J. Millero, Richard A. Feely, and Yong-Ming Liu. 2010. “The Universal Ratio of Boron to Chlorinity for the North Pacific and North Atlantic Oceans.” Geochim. Cosmochim. Acta 74 (6): 1801–11.

Leggett, Richard M., Ricardo H. Ramirez-Gonzalez, Bernardo J. Clavijo, Darren Waite, and Robert P. Davey. 2013. “Sequencing Quality Assessment Tools to Enable Data-Driven Informatics for High Throughput Genomics.” Frontiers in Genetics 4 (December): 288.

Li, Bo, and Colin N. Dewey. 2011. “RSEM: Accurate Transcript Quantification from RNA-Seq Data with or without a Reference Genome.” BMC Bioinformatics 12 (August): 323.

Martin, Marcel. 2011. “Cutadapt Removes Adapter Sequences from High-Throughput Sequencing Reads.” EMBnet.journal 17 (1): 10–12.

Pierrot, D. E., D. W. R. Wallace, and E. Lewis. 2011. “MS Excel Program Developed for CO2 System Calculations.” Carbon Dioxide Information Analysis Center. https://doi.org/10.3334/cdiac/otg.co2sys_xls_cdiac105a.

Puritz, Jonathan B., Christopher M. Hollenbeck, and John R. Gold. 2014. “dDocent: A RADseq, Variant-Calling Pipeline Designed for Population Genomics of Non-Model Organisms.” PeerJ 2 (June): e431.

Ries, J. B., A. L. Cohen, and D. C. McCorkle. 2009. “Marine Calcifiers Exhibit Mixed Responses to CO2-Induced Ocean Acidification.” Geology 37 (12): 1131–34.

Roy, Rabindra N., Lakshimi N. Roy, Kathleen M. Vogel, C. Porter-Moore, Tara Pearson, Catherine E. Good, Frank J. Millero, and Douglas M. Campbell. 1993. “The Dissociation Constants of Carbonic Acid in Seawater at Salinities 5 to 45 and Temperatures 0 to 45°C.” Marine Chemistry 44 (2-4): 249–67.

 

General Notes

Note 1: For a complete step-by-step description of the RNA and DNA Methylation QC, mapping, and quantification pipelines see the manuscript github repository - `https://github.com/epigeneticstoocean/AE17_Cvirginica_MolecularResponse`.

Note 2: Using scripts will require modifying data pathways. Statistical analyses use relative paths that correspond to the manuscript github repository - `https://github.com/epigeneticstoocean/AE17_Cvirginica_MolecularResponse`.

Note 3: Raw sequence files corresponding to this study are located on NCBI with the BioProject ID - PRJNA594029. Reference files and annotations based on NCBI C. virginica genome files (002022765.2_C_virginica-3.0).

Note 4: Additional intermediate products can be found on the github repository - `https://github.com/epigeneticstoocean/AE17_Cvirginica_MolecularResponse`.

Files

Seawater Chemistry

1. `/AE17_WaterChemistry _carbChem.csv` : Biweekly water samples (Total Carbon and Alkalinity) and carbonate chemistry calculations using CO2Sys.
2. `/AE17_WaterChemistry _carbChem.xlsx`: Biweekly water samples (Total Carbon and Alkalinity) and carbonate chemistry calculations using CO2Sys.
3. `/AE17_WaterChemistry_weekly.csv` : Triweekly water samples of `Temperature`, `Salinity`, and `pH`.
4. `/AE17_WaterChemistry_weekly.xlsx` : Triweekly water samples of `Temperature`, `Salinity`, and `pH`. Phenotypes `/AE17_metadata_ExperimentalExposureSamples_Simple.xlsx` : Contains sample meta data (treatment, time sampled, etc).
5. `/AE17_metadata_ExperimentalExposureSamples_Simple.csv` : Contains sample meta data (treatment, time sampled, etc).
6. `/AE17_CalcificationInfo_ExperimentalExposureSamples_Simple.xlsx` : Contains calculated calcification data.
7. `/AE17_CalcificationInfo_ExperimentalExposureSamples_Simple.csv` : Contains calculated calcification data.
8. `/AE17_CalcificationComplete.csv` : Calcification with associated water chemistry data generated with `/AE17_mergePhenotypeAndSW.R`.
9. `/AE17_EPFpHData_ExperimentalExposureSamples_Simple.xlsx` : Contains pH measurements of extra-pallial fluid.
10. `/AE17_EPFpHData_ExperimentalExposureSamples_Simple.csv` : Contains pH measurements of extra-pallial fluid.
11. `/AE17_EPFpHComplete.csv` : EPF pH measurements with associated water chemistry and calculated delta pH generated with `AE17_mergePhenotypeAndSW.R`.

RNA-Seq

1. `/RSEM_gene_Summary.Rdata` : RSEM transcript quantification `RData` file (all counts).
2. `/RSEM_gene_EstCount.csv` : Estimated counts from RSEM.
3. `/RSEM_gene_FPKM.csv` : FPKM counts from RSEM.
4. `/RNAseq/RSEM_gene_TPM.csv` : TPM counts from RSEM.
5. `/Target_BiomineralizationGenes.csv` : Target list of biomineralization genes.
6. `/RNA_gene_postVoomAndNormalization_DGEListObj.RData`: Normalized gene count matrix (edgR and Voom pipeline).
7. `/RNA_Limma_Expression_Data_forWGCNA.RData`: WGCNA expression object (Figure 7).
8. `/RNA_Limma_Expression_networkConstruction_WGCNA.RData`: WGCNA network object (Figure 7).

MBD-BSseq

1. `/methylKitObj_cov5Filtered_medianNormalized_united.RData` : MBD-BSseq summary counts in RData file (all counts).
2. `/countMatrix_cov5Filtered_medianNormalized_methylCCounts.csv` : Methylated Cytosines counts.
3. `/countMatrix_cov5Filtered_medianNormalized_totalCounts.csv` : Total Cytosines counts (coverage).
4. `/countMatrix_cov5Filtered_medianNormalized_unMethylCCounts.csv` : Unmethylated Cytosines counts.
5. `/countMatrix_cov5Filtered_medianNormalized_beta.tar.xz` : DNA Methylation proportion, beta (compressed).
6. `/CpGCoordinates.csv` : List of CpG genomic coordinates.
7. `/20200202_geneBeta_cov5_byFeature.RData`: Beta values per gene.
8. `/DNAm_methylkit_dm50_q01_sigSummaryTable.csv`: Differentially methylated loci (DML).

References

1. `/CDS_wGoTerms_GeneLOC.RData`: Coding regions with full GO term descriptions.
2. `/STAR_gnomon_tximportGeneFile.RData`: Gene description file.

Scripts

(see https://github.com/epigeneticstoocean/AE17_Cvirginica_MolecularResponse for additional notes on using these scripts)

1. RNA-Seq Pipeline
   1. `02_STAR_genomeCreate.sh` : Script for creating index for STAR mapping.
   2. `03A_STAR_1Pass.sh` : Script 1st pass of STAR mapping.
   3. `03B_STAR_2Pass.sh` : Script of 2nd pass of STAR mapping.
   4. `04A_RSEM_createRefFromStar.sh` Script for creating genome reference for RSEM quantification.
   5. `04B_RSEM_calcExp.sh` : Script for RSEM quantification.
   6. `05_filtering_CreatingDGEListObj.R` : Script for creating standardized DGEList object of RNAseq data.
   7. `06_CreatingWGCNAObj.R` : Script for performing WGCNA clustering and object creation.
2. MBD-BSseq Pipeline
   1. `00_concateRename.R` : Used to merge data from sequence facility (not needed when downloading files from NCBI).
   2. `01_seqQualityTrim.sh` : Script for performing QC and trimming step on raw reads.
   3. `02_bismarkMapping.sh` : Script for performing mapping using bismark.
   4. `03_cytosineReport.sh` : Script for generating cytosine reports (estimates of methylation for all CpGs) for each individual.
   5. `04_methylationMatrix_diffMethylation_methylKit.R` : Script for generation methylation matrices and differential methylation analyses using MehtylKit.
   6. `05A_CpGIntersectionByFeature.sh` : Script for intersecting CpGs by genomic feature.
   7. `05B_CpGCountByFeature.sh`: Script for counting CpGs per genomic feature.
3. Statistics Scripts
   1. `AE17_diffExpression.R` : Script for performing differential expression using DGEList object.
   2. `AE17_fig1_EPFtimeseries.R` : Script for analyzing complete EPF pH timeseries data (Figure 1).
   3. `AE17_fig2_EPFfinal_calcification.R` : Script for analyzing long-term EPF pH and calcification data (Figure 2).
   4. `AE17_fig3_DNAm.R` : Script for performing PERMANOVA, dapc, and other genome-wide analyses with the MBD-BSseq data (Figure 3).
   5. `AE17_fig4_geneExpression.R` : Script for performing PERMANOVA, dapc, and other genome-wide analyses with the RNA-seq data (Figure 4).
   6. `AE17_fig5_DNAmvsGE.R` : Script for examining association between DNA Methylation and gene expression data (Figure 5).
   7. `AE17_fig6_diffDNAmvsGE.R` : Script for evaluation association between differential gene expression and difference in DNA methylation (Figure 6).
   8. `AE17_WGCNAmultiComp.R` : Script for examining correlation between gene body DNA methylation, gene cluster expression, and EPF pH (Figure 7).
4. Additional
   1. `AE17_mergePhenotypeAndSW.R` : Script used to merge meta data (`/AE17_metadata_ExperimentalExposureSamples_Simple.csv`), water chem data (`/AE17_WaterChemistry_weekly.csv`) and phenotype data (`/AE17_EPFpHData_ExperimentalExposureSamples_Simple.csv`;`/AE17_CalcificationInfo_ExperimentalExposureSamples_Simple.csv`) and create merged table (`/AE17_EPFpHComplete.csv`;`/AE17_CalcificationComplete.csv`).
   2. `basicR_functions.R` : List of simple functions used in other scripts.