Molecular basis of parental contributions to the behavioural tolerance of elevated pCO2 in a coral reef fish

Knowledge of adaptive potential is crucial to predicting the impacts of ocean acidification (OA) on marine organisms. In the spiny damselfish, Acanthochromis polyacanthus, individual variation in behavioural tolerance to elevated pCO2 has been observed and is associated with offspring gene expression patterns in the brain. However, the maternal and paternal contributions of this variation are unknown. To investigate parental influence of behavioural pCO2 tolerance, we crossed pCO2-tolerant fathers with pCO2-sensitive mothers and vice versa, reared their offspring at control and elevated pCO2 levels, and compared the juveniles' brain transcriptional programme. We identified a large influence of parental phenotype on expression patterns of offspring, irrespective of environmental conditions. Circadian rhythm genes, associated with a tolerant parental phenotype, were uniquely expressed in tolerant mother offspring, while tolerant fathers had a greater role in expression of genes associated with histone binding. Expression changes in genes associated with neural plasticity were identified in both offspring types: the maternal line had a greater effect on genes related to neuron growth while paternal influence impacted the expression of synaptic development genes. Our results confirm cellular mechanisms involved in responses to varying lengths of OA exposure, while highlighting the parental phenotype's influence on offspring molecular phenotype.


Introduction
As atmospheric CO 2 levels increase, so does the amount of CO 2 taken up by the ocean, causing a decrease in seawater pH (ocean acidification (OA)), with potentially broad-ranging effects on the physiology and ecology of marine organisms [1]. However, the effects of OA on marine ecosystems will depend on the relative sensitivity and tolerance of different species, and their ability to acclimatize and adapt to the environmental changes brought about by rising CO 2 levels [2]. Recent studies show that increased pCO 2 levels can affect the growth, survivorship and physiology of some marine fishes [3][4][5][6][7]. Elevated CO 2 has also been shown to alter behaviours in a wide variety of fish [8][9][10] and invertebrates [11][12][13]. However, most experimental studies focus on short-term exposure to elevated pCO 2 and do not account for phenotypic variation that may enable populations to adapt over the time scale at which OA will occur [2,14].
Fish use chemical alarm cues (CACs) from conspecifics to detect predation threats and respond to CAC presence by moving away from the cue and reducing activity [15]. Under elevated pCO 2 conditions, juvenile fish can exhibit an impaired response to CAC, showing a decreased avoidance to CAC, and in some cases failing to associate CAC with the threat of predation [16][17][18][19]. The underlying cause of this behavioural change is thought to be impaired function of GABA A receptors (GABA is the major inhibitory neurotransmitter in the brain) caused by changes in the concentration of acid-base relevant ions to maintain pH homeostasis in a high CO 2 environment [17,[20][21][22]. Elevated CO 2 has also been demonstrated to affect learning ability, decision making, turning preference, auditory preferences, visual acuity, shoaling, boldness and escape responses in a range of different fish species [16,18,[23][24][25]. Such alterations in behaviour could affect individual performance and survivorship, with potential implications for community structure and population replenishment [8,26]. However, there is also individual variation in behavioural sensitivity to elevated CO 2 , with some individuals exhibiting impaired behaviours, whereas others do not, especially at CO 2 levels projected to occur in the ocean this century (i.e. approximately 700 µatm CO 2 ) [23,27,28]. Such individual variation could be the raw material for population-level adaptation to rising CO 2 levels [14,29].
Genetic variation among parents, combined with nongenetic effects from the parental environment, drive an offspring's phenotype through natural selection, epigenetic inheritance and changes in molecular mechanisms [30]. These parental effects, genetic and non-genetic, are key to understanding how populations will adapt and survive in the face of climate change. Both mothers and fathers make important contributions to the success and development of their offspring, yet recent studies have focused singularly on either maternal or paternal effects, rarely examining how they may interact to shape their offspring's performance [31]. In the marine stickleback, maternal exposure to high temperatures had a clear impact on the growth of their offspring in the same environmental conditions [32]. Maternal inheritance of increased metabolic capacity led to larger juveniles, suggesting transgenerational plasticity can mediate short-term impacts of ocean warming [32]. Another recent study, on wild salmon, found a correlation between telomere length in juveniles and paternal time spent in seawater that appeared unaffected by increased temperatures [33]. Short telomeres can indicate poor biological health, suggesting that fathers with an increased time at sea produce healthier offspring even in the face of environmental variability [33]. There is a significant lack of knowledge about sex-specific parental contributions in fish, especially under elevated pCO 2 conditions, and the interactions of genetic and nongenetic effects add to the complexity [34]. Identifying the parental contributions to specific traits and the molecular mechanisms behind them will allow us to better predict how fish will respond and adapt to rapid climate change.
Recent studies on the spiny damselfish, Acanthochromis polyacanthus, indicate that individual variation in sensitivity of the behavioural reaction to elevated pCO 2 is heritable and could thus facilitate adaptation [28,35,36]. A portion of fish from the natural population respond normally to CACs under predicted end-of-century pCO 2 levels of 700-800 µatm and thus appear to be tolerant to the behavioural effects of OA. Moreover, there is a strong correlation between the behavioural tolerance of juvenile spiny damselfish to elevated pCO 2 and the tolerance of their fathers in both wild and captive populations [28], indicating that behavioural tolerance to elevated pCO 2 has a genetic basis and is heritable. Furthermore, molecular differences were identified between offspring of tolerant versus sensitive parental pairs, with those under transgenerational exposure to elevated pCO 2 showing differential gene and protein expression in the brain [35,37]. Specific genes related to this signature of tolerance include those that regulate the circadian rhythm processes and could be creating a phase shift in the circadian clock to better control their acid-base regulation [35]. Recent evidence suggests that this behavioural variation may be passed on paternally through generations [28]; however, it is still unknown if this heritability is identifiable in transcriptional changes.
In this study, we investigate maternal and paternal influence on the molecular signature of behavioural tolerance to elevated pCO 2 in A. polyacanthus. Using a unique transgenerational experimental design, we cross-bred behaviourally tolerant mothers with behaviorally sensitive fathers and behaviourally sensitive mothers with tolerant fathers. We exposed breeding pairs and their offspring to ambient and elevated pCO 2 conditions. As described in previous studies [28,35,36], offspring were exposed to ambient pCO 2 (control), elevated pCO 2 from hatching (transgenerational high pCO 2 treatment), or for only 4 days before sampling (acute high pCO 2 treatments). Due to previous research implicating the role of neurotransmitters in behavioural changes under elevated pCO 2 , we continued to study the brain. We measured genome-wide gene expression in 68 juvenile fish from the different combinations of parental pairs (CO 2 -tolerant fathers paired with CO 2 -sensitive mothers versus CO 2 -sensitive fathers paired with CO 2 -tolerant mothers) and CO 2 treatment conditions, exploring the role of the maternal and paternal phenotype in influencing individual offspring expression profiles under various pCO 2 levels and exposure durations.

Material and methods (a) Experimental design
In this study, we sampled the brains of juvenile A. polyacanthus from the laboratory experiment described in a previous paper [28]. Full details of the experimental design are provided in that paper [28]. Briefly, adult A. polyacanthus were collected from the northern Great Barrier Reef in Australia and transported to the aquaria facilities at James Cook University (Townsville, Australia; JCU ethics permit A1828). Adults were placed in elevated pCO 2 (754 µatm) (electronic supplementary material, table S1) for 7 days then tested for their reaction to CAC in a twochannel flume. Their sensitivity to high pCO 2 was determined by the amount of time spent in the CAC: spending less than 30% of the time in the CAC was considered tolerant to elevated pCO 2 , and individuals spending greater than 50% were considered non-tolerant. Adults were then sexed and paired for breeding based on their behavioural sensitivity. For the purpose of this study, we focused only on the tolerant male × sensitive female (T♂S♀) and sensitive male × tolerant female (S♂T♀) pairs; however, all parental combinations were made for the larger experimental design [28,36]. Adults were then placed in 40 L aquaria at control (414 µatm) and elevated pCO 2 royalsocietypublishing.org/journal/rspb Proc. R. Soc. B 288: 20211931 (754 µatm) conditions for three months prior to breeding season (electronic supplementary material). For this specific study, we used clutches from 3T♂S♀ pairs and 3 S♂T♀ pairs in control and 2T♂S♀ pairs and 2 S♂T♀ pairs at elevated pCO 2 . Breeding pairs were checked daily for the presence of egg clutches. Consistent with previous observations for this species, there was no apparent effect of elevated pCO 2 on reproductive output, although this was not directly quantified in the current study. On the day of hatching (approx. 10 days post-fertilization), the offspring were transferred to either control or elevated pCO 2 conditions (figure 1). This led to two long-term conditions referred to as control (control parents, offspring reared in control pCO 2 ) and transgenerational (elevated pCO 2 parents, offspring reared in elevated pCO 2 ). Furthermore, we exposed half of the offspring from the control condition to elevated pCO 2 (754 µatm) 4 days prior to sampling at the end of the experiment. These offspring came from parents exposed to either control conditions or elevated pCO 2 creating two separate offspring acute treatments named for their parental conditions: control-acute and high pCO 2 -acute (figure 1). No significant juvenile mortality was recorded throughout the experiment across all conditions. After five months, up to four fish per parental pair were euthanized resulting in approximately nine individuals from each parental type in each of the four conditions for a total of 68 fish (n = 16 or 18 from each of the four conditions, electronic supplementary material, table S2). Brains were dissected out and immediately flash frozen in liquid nitrogen and then kept at −80°C for further processing. Different fish were euthanized for this study than those used by Welch & Munday [28] to measure behavioural response to CAC. In the previous study, offspring were tested after six weeks in a two-channel flume to determine their behavioural phenotype [28]. For this study, we sampled fish after five months to ensure sufficient brain tissue was available for the molecular analysis. We did not perform any behavioural testing prior to sampling at five months to avoid potential influence of handling disturbance and the testing regime on the brain transcriptome.

(b) RNA extraction and preparation
Whole-frozen fish brains were homogenized in RLT-plus buffer for 30 s using a bead beater (BioSpec) with single-use silicone beads and total RNA was extracted using the AllPrep DNA/ RNA mini kit following the provided protocol (Qiagen). RNA quantity was determined using a Nanodrop (Thermo Scientific) then reanalysed for both the quality and quantity with the Agilent 2100 Bioanalyzer. Extracted RNA was processed at the King Abdullah University of Science and Technology (KAUST) Bioscience Core Lab for library preparation and sequencing. Conversion to cDNA and preparation for Illumina sequencing was done using the TruSeq RNA Illumina Library Prep Kit. Samples were sequenced on the Illumina HiSeq 2500 paired end at 100 bp in KAUST, Saudi Arabia.

(c) Gene expression analysis
Raw reads were quality checked with FastQC [38] and sequences were trimmed using Trimmomatic v. 0.36 [39] with set parameters including a phred score of 33, a minimum length of 40 and removal of the first 13 bp of each sequence as well as the Illumina adapter sequences. Reads were rechecked in FastQC for quality to verify correct adapter trimming. Sequences were then aligned to the Acanthochromis polyacanthus genome (ENSEMBL ASM210954v1) using hisat2 v. 2.1.0 [40] with an average successful alignment rate of approximately 80%. SAM output files from mapping were then sorted and converted into bam files using SAMtools v. 1.5. Feature counts (Subreads v. 1.5.3) [41] was used to create a matrix of raw exon read counts and the gene-ids provided in the genome annotation. Minimum mapping quality was set to 20 and both fragment pairs were required to align to the reference. The data matrix was then imported into R v. 3.4.0 and differential expression between conditions and parental pairs were statistically determined with DESeq2 [42]. A likelihood ratio test (LRT) was used to determine the interactions between the multiple variables as well as the best design formula for further analysis ( p-adj < 0.05). Significant differentially expressed genes (DEGs) were identified using a q-value < 0.05 and a minimum log twofold change of 0.3. All graphs were created using the R packages Vegan and ggplot2 [43,44].
Gene ontology of each gene was annotated using Blast2Go Pro v. 4.19 (Gene Ontology Consortium [45]. Enrichment was performed in RStudio v. 1.1.463 with the GO-MWU script that uses a Mann-Whitney U test and a Benjamini-Hochberg correction to determine statistical significance [46], (scripts available at: https://github.com/z0on/GO_MWU). The enrichment test was carried out separately for each of the different GO domains using stringent filtering including collapsing redundant categories, only selecting those genes represented by more than 5 GO-Terms, and employing an FDR cutoff of 0.1.

Results (a) Condition-specific expression patterns
Duration of exposure to elevated pCO 2 had the greatest effect on expression profiles of A. polyacanthus offspring. The number of DEGs was higher when comparing the two acuteelevated pCO 2 treatments to control than in the comparison of the transgenerational-elevated pCO 2 treatment to control (figures 2 and 3). The influence of length of exposure to elevated pCO 2 on patterns of gene expression in the brain was also seen by the high number of overlapping DEGs between offspring of both acute treatments (1642), whereas offspring There were 264 DEGs when comparing different environmental pCO 2 levels (control versus transgenerational-elevated pCO 2 ) for offspring of tolerant mothers, but only 68 DEGs when comparing these treatments for offspring of tolerant fathers (figure 3). Only 11 of these genes were shared between offspring of the two parental groups, including upregulation of a gene known to protect cells from reactive oxygen species (GSTA1) and a cytoskeleton-related gene (CKAP2). Genes significantly differentially expressed in the offspring of tolerant mothers included several related to growth (TNNT3, MYSS and MLRV). Further analysis on total body weight revealed that offspring of tolerant mothers were on average smaller than those of tolerant fathers, with a significant difference between the two in the transgenerational-elevated pCO 2 treatment (electronic supplementary material, figure S1 and table S3).
GO enrichment revealed upregulation of the circadian rhythm pathway as well as circadian regulation of gene expression in offspring from tolerant mothers, including the circadian repressors PER3, CIART and HT7R ( figure 4). Unique differentially enriched pathways in offspring of tolerant fathers when comparing control to transgenerational included downregulation of reactive oxygen species biosynthesis and the H4 histone complex as well as upregulation of axon development ( figure 5).

(c) Influence of parental identity in acute-elevated CO 2 treatments
Fewer DEGs due to parental identity were identified in both acute-elevated pCO 2 treatments compared with the control and transgenerational-elevated pCO 2 treatment. Two hundred and seventy-two DEGs were found between offspring of T♂S♀ (tolerant fathers) and S♂T♀ (tolerant mothers) parents in the control-acute treatment, and 270 DEGs between offspring of the two parental pairs in the high pCO 2 -acute treatment ( figure 3). There was a strong common response to the acute-elevated pCO 2 treatments versus control, regardless of the parental condition or parental identity, which consisted of 1642 DEGs. The enriched functions for these shared genes include collagen binding, neurotransmitter receptor activity and GABA receptor activity (figure 5). GABA-related genes included KCC1, NXPH1 and NXPH2, several solute carriers specializing in GABA transport, and GAD2 involved in GABA catalysis. Carbonic anhydrase (CA) genes (CAH4, 8, 10 and 15) were also differentially expressed between control and both acute-elevated pCO 2 conditions. Circadian rhythm activators (RORa and b) involved in stability of the clock were highly upregulated in all acute treatments when compared to control. Comparing control to control-acute and high pCO 2 -acute led us to identify 604 shared DEGs from offspring of tolerant fathers and 567 shared DEGs from offspring of tolerant mothers. In offspring of tolerant mothers, SC6A1 and S6A11 were both upregulated in the acute-elevated pCO 2 treatment. Offspring of tolerant mothers exhibited highly upregulated glutamate-related genes (NMD3B and GRIK2), as well as genes related to neuron plasticity, growth and development (AMIGO1, GFRA2 and CDK5R1). Several genes downregulated in offspring of tolerant mothers exposed to acuteelevated pCO 2 conditions were related to development and growth (TNNI2, TNNC2, MYSS, ACTA1 and ACTN3). Pathways unique to the offspring of tolerant mothers were nervous system development, serine family biosynthesis and RNA methylation (figure 5). In offspring of tolerant fathers, solute carrier genes involved in GABA transport were both upregulated (S6A11 and S6A13) and downregulated (SC6A1) in acute-elevated pCO 2 conditions. These offspring also showed upregulation of genes related to synaptic plasticity (MPDZ) and inhibition of synaptic development (IGSF9B) as well as a cytoskeleton gene (ANK1) and growth factor genes (FGF6 and FGF14) involved in cell proliferation and nervous system development. Downregulated genes including transcriptional repressors (HES5 and 6) as well as NEUROG1 involved in transcriptional regulation and cell differentiation were identified in offspring of tolerant fathers under acute-elevated pCO 2 conditions. Also identified as downregulated genes were a histone-binding gene (NASP), a gene involved in immunity and cell death (PERF), and one involved in DNA repair (PAF15). These genes corresponded to unique pathways such as histone binding, ATPase activity, hydrogen transport and microtubule polymerization enriched in offspring of tolerant fathers ( figure 5).

Discussion
In this study, we show that parental phenotype shapes the transcriptomic response in the brain of spiny damselfish offspring, regardless of the duration or level of pCO 2 exposure. Within the control and transgenerational-elevated pCO 2 treatments, we found an almost 10-fold increase in the number of significantly differentially expressed genes in offspring from CO 2 -tolerant fathers versus CO 2 -tolerant mothers compared to the number of DEGs between the control and transgenerational-elevated pCO 2 treatments. A notable response uniquely upregulated in offspring of tolerant mothers was circadian rhythm regulation. This included circadian repressors that were previously found upregulated in offspring of tolerant parents when both parents exhibited tolerance in the same condition [35]. Circadian rhythm genes are closely linked to homeostatic ion-regulatory adjustments and have been identified as an important component in the response of fish to elevated pCO 2 [35,[47][48][49]. Previous research linked this circadian rhythm shift to adaptive potential in the face of OA that could be inherited across a generation [35]. The ability of offspring from tolerant mothers to make changes in the circadian clock could lead to greater flexibility of ion-control and thus avoid maladaptive behavioural reactions to elevated pCO 2 . Research on parental effects on circadian regulation especially in fish are lacking; however, several studies in mice have confirmed both the role of parental identity and maternal effects on circadian rhythm regulation of their offspring [50,51]. For example, Borengasser et al. [50] identified the effects of maternal obesity on circadian disruptions via suppression of key genes leading to detrimental metabolic effects in the offspring. This is consistent with our findings that the mother's phenotypic identity influences circadian rhythm regulation in juvenile fish exposed to elevated CO 2 . While there is an effect of parental identity on gene expression profiles with the transgenerational and control exposures, this identity weakens under acute exposure to elevated pCO 2 . The decrease in significant differential gene expression between parental phenotypes is likely to be due to an intense cellular response to the acute exposure to elevated pCO 2 [52]. This population-wide response could minimize the observed molecular variation based on parental phenotype. Importantly, we found a common response to acute exposure to elevated pCO 2 , regardless of parental condition or phenotype, which consisted of the upregulation of GABA receptor activity and neurotransmitter receptor pathways. GABA receptors play a pivotal role in maintaining the inhibitory-excitatory balance in the brain; under elevated CO 2 these receptors are thought to become overwhelmed by altered neuronal gradients of Cl − and HCO 3− attempting to maintain pH homeostasis [20]. A recent study [21] highlights the impact of near-future pCO 2 levels on the creation of a vicious cycle caused by functional alterations to GABA receptors and ultimately leads to the behavioural impairments observed in fish [18,24,36]. This vicious cycle can quickly cause increased excitatory ion fluctuations in the presence of relatively small changes in pCO 2 levels, such as those experienced in the acute treatments. Elevated CO 2 could potentially have other neurological effects. For example, glycine receptors are also Cl − /HCO 3− channels and thus could potentially be affected by elevated CO 2 in a similar way to GABA receptors [53], although this has not been tested.
Indeed, we found several glycine receptors significantly upregulated in both acute-elevated pCO 2 conditions and offspring of both parental phenotypes. CA genes, another key player in acid-base regulation, were also identified as upregulated in both acute-elevated pCO 2 conditions when compared to control, regardless of parental phenotype. CAs are highly abundant in neurons and play an important role in the hydration of CO 2 in the cell [7,54,55]. Upregulation of these genes suggests an elevation of intracellular HCO 3 − that would further excite the GABA receptors and contribute to the vicious cycle, thereby leading to an ion imbalance in the brain. Previous studies only identified CA genes to be downregulated in fishes; however, these studies focused on the gill tissues which could have different regulatory mechanisms in response to elevated pCO 2 than the brain [56,57]. Identification of these upregulated genes in all offspring regardless of parental phenotype suggests tolerance or sensitivity to elevated pCO 2 cannot protect these fish from this vicious cycle under short-term exposures. However, the absence of differential expression of these genes and pathways in offspring exposed transgenerationally to elevated pCO 2 suggests that acclimation via stabilization of acid-base regulation occurs with long-term exposure [58].
Welch & Munday [28] found that heritability of behavioural tolerance to increased pCO 2 in the fish from this experiment was only detectable in the acute exposure conditions; therefore, we endeavoured to find a molecular response correlated with the observed behavioural response between offspring of different parental pairs in the sampled acute treatments. Upregulation of genes and pathways associated with neuron plasticity and development, as well as glutamate-related genes were identified in offspring of tolerant mothers (S♂T♀). Glutamate is a key part of GABA synthesis, as well as other neurotransmitter pathways, and upregulation of these genes could be part of a compensatory mechanism along with neuronal plasticity in the face of short-term increases in pCO 2 [18,24,36]. The identification of downregulated growth-associated genes suggests a negative impact of elevated pCO 2 on development in juvenile fish and is consistent with previous research [7,59,60]. This suggests that the maternal line could be primarily responsible for neural plasticity used to compensate for changes in seawater chemistry.
In offspring of tolerant fathers (T♂S♀), we identified upregulated genes and pathways directly related to neuronal plasticity, which suggests these offspring also have compensatory neural mechanisms in response to acute increases in pCO 2 . Histones and transcriptional repression mechanisms were found to be downregulated in this study, and in the previous study by Schunter et al. [36] that investigated offspring of fully tolerant parents, under acute-elevated pCO 2 conditions. Histones control chromatin dynamics and also play a part in gene expression by regulating transcription [61]. Decreases in repression of transcription and histone binding lead to increases in gene expression, suggesting both a genetic and epigenetic control of parental phenotype on the offspring's expression profile in response to acutely elevated pCO 2 . Epigenetic changes can be heritable through the germ line and can be strongly influenced by environmental factors that induce specific phenotypes in parents and their offspring [62,63] A recent study on medaka identified paternally influenced changes in methylation patterns across multiple generations of male fish after initial exposure to hypoxia [62], whereas studies on zebrafish suggest that maternal methylation is reprogrammed in the embryo [66]. Therefore, it is possible that epigenetic influence on the offspring is more likely to come from fathers than from mothers [62,66]. Previous research on A. polyacanthus has already shown that epigenetics plays a key role in transgenerational acclimation to increased temperatures via histone regulation and selective DNA methylation [67,68]. This leads us to speculate an important role for epigenetic inheritance from fathers in acclimation to elevated pCO 2 .
By combining a multi-generational experimental design with genome-wide gene expression measurements, we were able to identify effects of parental phenotypes on expression patterns of their offspring to end-of-century pCO 2 levels. Our results suggest that there are both maternal and paternal contributions involved in OA tolerance and that their relative importance may differ depending on length of exposure to elevated CO 2 . Although maternal and paternal influence on specific traits is a growing research field, few studies have attempted to distinguish the separate parental contributions in transgenerational experiments [69]. In the marine stickleback, maternal influences were shown in growth rates and genetic covariance under elevated temperatures [32], but in wild salmon a correlation between telomere length and paternal growth was identified [33]. Therefore, it appears that parental contributions may vary between species as well as between specific traits. This is supported by our results; for example, maternal tolerance impacted circadian rhythm genes in their offspring while epigenetic modifications were influenced in the tolerant paternal line. The role of epigenetics in transgenerational acclimation to future climate change scenarios is also seen in recent research on DNA methylation of A. polyacanthus [67,68].
Further research is necessary to determine the exact mechanisms of maternal and paternal effects. Transgenerational studies spanning further generations could give insight into whether inheritance of pCO 2 tolerance is maintained across multiple generations, as well as clarifying the developmental time point when these expression changes occur. Also, more complex experimental designs could help clarify the relative importance of genetic versus non-genetic effects and tease apart the maternal versus paternal components of these different mechanisms. Nevertheless, it is clear from this study that variations in parental phenotype can be highly important in shaping the response of fish to future OA conditions. Data accessibility. RNAseq data for all individuals can be found under the BioProject ID PRJNA690095. All other data are provided in the electronic supplementary material [70].