On January 22, 2012 (the new moon), eight colonies of adult P. damicornis were collected from a depth of ~1–3 m on a fringing reef (17.4803 S, 149.7989 W) on Moorea, French Polynesia. Average seawater temperature on the fringing reef was 28.19 ± 0.50°C during the week of coral collection (Rivest, 2014). Colonies were individually held in aquaria that received indirect natural sunlight and coarsely filtered seawater. Daytime temperatures in the aquaria varied between 26.0 and 29.5°C. Larvae were collected from adult colonies following Rivest and Hofmann (2014). The “peak” larvae, determined by the predicted relationship of larval release according to the lunar cycle (January 30, 2012; e.g., Fan et al., 2006; Rivest and Hofmann, 2014, 2015), were used in the experiment. Due to the smaller larval production by some colonies, all larvae released were pooled. A portion of larvae from this pool were immediately photographed and preserved to quantify larval quality metrics and gene expression (see below) that represented the condition of freshly-released larvae. The remaining larvae in the pool were randomly assigned to experimental treatments immediately after collection.

Two pCO₂ levels were prescribed: ~465 μatm pCO₂ and ~1,030 μatm pCO₂. The control (ambient) treatment approximated the present-day environmental conditions for the adult coral (verified by environmental data; Rivest and Hofmann, 2015), while the high treatment represents a level of pCO₂ expected in the surface ocean by year 2100 (RCP8.5 scenario; IPCC, 2013). Two experimental temperatures (T) were crossed with the pCO₂ levels: 28 and 31°C. The control (ambient) temperature approximated the multi-year average temperature (20-min sampling frequency) for the fringing reefs close to the collection sites for adult P. damicornis (Leichter, 2015). The elevated temperature represents an average surface ocean temperature projected by year 2100 (RCP8.5 scenario, IPCC, 2013). The four treatments created by this experimental set-up are defined as ambient temperature–ambient pCO₂ (ATAC), ambient temperature–high pCO₂ (ATHC), high temperature–ambient pCO₂ (HTAC), and high temperature–high pCO₂ (HTHC).

Treatment conditions were generated in a Moorea Coral Reef Long-Term Ecological Research program (MCR LTER) facility as described previously (Rivest and Hofmann, 2015), with two aquaria for each combination of pCO₂ and temperature treatment levels (8 aquaria total). Photosynthetically-active radiation (PAR) was provided at 1,800 ± 140 μmol m⁻² s⁻¹ on a 13h:11h: light:dark cycle (Sol LED Module, Aquaillumination, 75 W, Ames, IA, USA), an ecologically relevant level for shallow reef environments studied by MCR LTER (Carpenter, 2016). To verify and monitor the physical parameters of the pCO₂ × temperature treatments, the carbonate chemistry (pH and total alkalinity), salinity, and temperature of the seawater in the aquaria was measured during the incubations. For a description of the methods, see Rivest and Hofmann (2015). pH and pCO₂ of the treatments were calculated from pH at 25°C, total alkalinity (A_T), temperature, and salinity using CO₂calc (Robbins et al., 2010). CO₂ constants K1 and K2 from Lueker et al. (2000) were used, and pH was expressed on the total scale (mol kg SW⁻¹).

Larvae were incubated under experimental conditions for 24 h in two 400 mL containers per aquarium at 0.15 larva mL⁻¹. Larval containers had 100 μm mesh sides and a PAR-transparent lid. Incubations were staggered by 1 h per aquarium, with the order randomized, due to the time needed to process the larvae at the end of the incubation. At 24 h, larvae within tanks were pooled, and 10 larvae were randomly selected for size measurements (n = 20). The remaining larvae were aliquoted and preserved at −80°C for downstream analyses of gene expression (25 larvae per aliquot), total protein content (5 larvae per aliquot), and symbiont density (5 larvae per aliquot). Therefore, for comparisons between treatments, n = 4 for gene expression, total protein, and symbiont density.

Larval size, protein content, and Symbiodinium density were measured as described in Rivest and Hofmann (2015). Briefly, larval size was quantified through image analysis using ImageJ (Rasband, 1997). Protein content was measured using the Bradford assay (Bradford, 1976; Rivest and Hofmann, 2014). Symbiont density was quantified from larval homogenates using a haemocytometer.

All data were analyzed using R v3.0.1 (R Core Team, 2013). Effects of pCO₂ and temperature on tank conditions and larval quality response variables were estimated using linear mixed-effect models (nlme package in R; Pinheiro and Bates, 2000) with pCO₂ and temperature treatments as fixed factors and “tank” as a random factor (following Rivest and Hofmann, 2015). Statistical assumptions of normality and homogeneity of variance were met, sometimes following a power-based transformation of the response variable.

RNA was extracted from aliquots of 25 P. damicornis larvae homogenized in 0.5 mL TRIzol® using glass beads and a bead-beater. Homogenates were brought to a total volume of 2 mL TRIzol® and split into two microcentrifuge tubes. Proceeding with one half of the homogenate, RNA was isolated using the guanidinium thiocyanate method (Chomczynski and Sacchi, 1987). To improve the quality of the RNA extracts, a Qiagen RNeasy® MinElute® Cleanup Kit was used, following the manufacturer's instructions. After the clean-up, RNA quality was assessed using an Agilent RNA 6000 Pico Kit and a Bioanalyzer Agilent 2100 instrument, following the manufacturer's instructions. The electropherograms showed clean, resolved peaks with no obvious signs of RNA degradation and consistent quality across all extracts.

cDNA libraries were prepared from at least 1 μg total RNA from each extract using a TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA USA), following the manufacturer's instructions. From the two tubes of larvae preserved from each aquarium, one was used to prepare a library for paired-end RNA sequencing and the other was used to prepare a library for single-read RNA sequencing. A unique adapter index was ligated to the cDNA fragments in each library. The libraries contained cDNA sequences of insert size ~260 bp.

Pooled, barcoded libraries were sequenced at the UC Davis Genome Center DNA Technologies Core (http://dnatech.genomecenter.ucdavis.edu/) on the HiSeq 2500 Illumina sequencing platform. Nine libraries were pooled in each of 2 lanes. One lane was used to generate 100 bp paired end reads, and the second lane was used to generate 50 bp single reads, resulting in 187,925,419 paired-end and 157,325,043 single-end reads, respectively. Adapter contamination and bases with Phred quality scores < 20 were removed with Trim Galore v3.23 (Trim Galore!, Bioinformatics Group, Babraham Institute) with Cutadapt v1.2.1 (Martin, 2011) and FastQC v0.10.1 (Andrews, 2014). The reads generated were uploaded to the European Nucleotide Archive under the project accession PRJEB26650. The individual RNA-seq libraries are accessioned with numbers ERS2472791-ERS2472808.

A de novo transcriptome was assembled using Trinity v2013-8-25 (Haas et al., 2013) implemented on the Knot computing cluster at the UCSB Center for Scientific Computing (http://csc.cnsi.ucsb.edu/clusters/knot). Subsequent bioinformatics analyses were performed using high performance computing resources provided by Louisiana State University (http://www.hpc.lsu.edu). Because P. damicornis harbors Symbiodinium zooxanthellae, host and symbiont mRNA was present in our samples. A holobiont transcriptome was assembled and then filtered to separate the host coral and symbiont zooxanthellae contigs. Contigs in the holobiont transcriptomes were collapsed using CD-HIT v4.6 (Fu et al., 2012) with a 0.98 cutoff, and contigs <500 bp long were removed. The host and symbiont contigs were filtered from the holobiont transcriptome following the procedure outlined in Davies et al. (2016). Briefly, a series of blast searches were performed against four databases created from cnidarian and Symbiodinium genome and transcriptome sequences publicly available as of July 2016. The databases contained i) sequences from cnidarians (n = 28) that may contain Symbiodinium contamination, ii) sequences from cnidarians (n = 10) known to be free of Symbiodinium, iii) sequences from Symbiodinium that may contain cnidarian contamination (n = 19), and iv) sequences from Symbiodinium known to be free of cnidarian sequences (n = 13). The symbiont-free cnidarian database contained sequences obtained from aposymbiotic larvae and adult individuals (for example, deep sea corals). The cnidarian-free Symbiodinium databases contained sequences from free-living zooxanthellae or from lines cultured in the laboratory. A list of all genomes and transcriptomes used to construct the blast databases and their citation information is available in Supplementary Table 1. The databases were indexed, and translated blast searches (tblastx v2.4.0, Camacho et al., 2009) were performed on the coral holobiont transcriptome against all four databases. If a contig from the holobiont transcriptome had an amino acid overlap >100 bp with 80% identity cutoff to any cnidarian sequence, it was considered a host contig. However, if the same contig had a 100 bp amino acid overlap with 80% identity to a clean Symbiodinium sequence, it was removed from the coral host transcriptome. Symbiodinium contigs were identified in the same manner as described above, removing any contigs that were assigned to aposymbiotic cnidarian sequences with 100 bp amino acid overlap and 80% identity from the Symbiodinium transcriptome. After filtering, 1,200 contigs (1.5 and 2.6% of the coral and symbiont transcriptomes, respectively) that had equally good matches to both contaminant-free cnidarian and contaminant-free Symbiodinium blast databases were removed because it was not possible to determine whether they were of coral or symbiont origin. This conservative approach increases confidence that the findings for P. damicornis are based unambiguously on the response of the coral host and those for Symbiodinium are based unambiguously on the response of the symbiont.

To measure changes in gene expression associated with pCO₂, temperature, and the interaction of these environmental stressors, reads were mapped against the coral and symbiont transcriptomes and counted using RSEM v1.2.31 (Li and Dewey, 2011). Read counts for each technical replicate (one set of 50 and 100 bp reads each per tank) were added together before testing differential gene expression among treatments. Effects of temperature, pCO₂, and their interaction were tested via a likelihood ratio test in DESeq2 v1.12.4 (Love et al., 2014) in R v3.3.1 (R Core Team, 2013), and genes differentially expressed among treatments were identified with a p-value cutoff corresponding to an FDR of 0.05 (Benjamini and Hochberg, 1995). Overall patterns of gene expression were assessed by plotting a two-dimensional principal components analysis (PCA) of log-transformed gene counts. Two PCAs for both the coral and symbiont transcriptome data were performed: one comparing larvae among the experimental treatments and one comparing larvae from the experimental treatments to the initial sample of larvae collected immediately after release from the parents. Including the freshly-released larvae in our analysis allowed us to test the hypothesis that High-T and High-pCO₂ treatments resulted in developmental delay, causing larvae from these treatments to exhibit a gene expression pattern more similar to freshly-released larvae.

A rank-based gene ontology (GO) analysis was performed to identify GO terms that were significantly enriched by up or down-regulated genes following the method of Wright et al. (2015). First, coral and symbiont transcriptomes were translated into protein sequences in TransDecoder v3.0.1 (http://transdecoder.github.io/) using default parameters, and GO terms for the protein sequences were retrieved using InterProScan v5.22-61.0 (Jones et al., 2014). Next, the uncorrected p-values for differentially expressed genes in the DESeq2 analysis were converted by taking the log of the p-value and multiplying by −1 if the gene was down-regulated. After this conversion, highly down-regulated genes had highly negative values, and highly up-regulated genes had highly positive values. Using the GO terms list from InterProScan and the positive and negative log p-values for genes expressed under each condition, the Mann-Whitney U (MWU) test in R was used to determine functional enrichment of up and down-regulated genes. Scripts for this procedure and example files are available at https://github.com/z0on/GO_MWU.

Because P. damicornis can host multiple clades of symbionts (Putnam et al., 2012), and the Symbiodinium community could shuffle in response to changes in seawater conditions (Jones et al., 2008), the clade(s) of Symbiodinium present were identified for freshly-released larvae and larvae exposed to experimental treatments. Symbiodinium form eight distinct clades, which can be definitively characterized based on chloroplast 23S-rDNA and nuclear 28S-rDNA sequences (Pochon et al., 2006). To identify which clade(s) of Symbiodinium were present in the coral larvae sampled in this experiment, a local blast search (Zhang et al., 2004) of the symbiont transcriptome was performed against a database created from 23S and 28S genes from each of the eight Symbiodinium clades (see Table 1 in Pochon et al., 2006).

To quantify changes in gene expression between high-temperature, high-pCO₂ conditions and ambient conditions, a de novo reference transcriptome for P. damicornis larvae was assembled. The holobiont assembly yielded a transcriptome of 422.3 Mb with 169,472 contigs. The mean contig length was 1,279 bases, the N50 (the shortest sequence length such that 50% of the total sequence output is contained in sequences that are shorter) was 2,317 bases, and the GC content was 43.39%. The filtered reference transcriptomes contained 49,941 and 23,479 contigs for the coral and symbiont, respectively. The coral transcriptome had a N50 of 3,531 and a GC content of 40.65%, while the symbiont transcriptome had a N50 of 1,931 and a GC content of 55.41%.

After filtering, the total number of reads mapped per biological replicate ranged from 9.4 to 22.7 million for the coral transcriptome and from 1.7 to 5.4 million for the symbiont transcriptome. For the coral transcripts, only a small number of genes were differentially expressed (17 genes, FDR = 0.05, Supplementary Table 4), and a LRT indicated that the best-fit model for coral gene expression was one that included the interactive effects of temperature and pCO₂, but no main effects of either treatment (0 differentially expressed genes for either temperature or pCO₂ alone, FDR = 0.05). All 17 of these genes showed elevated expression in the ACAT and HCHT treatments, relative to the ACHT and HCAT treatments. Gene ontology categories for these genes (as determined by an InterProScan search) included transcription factor activity, ion and protein binding, nuclease and receptor activity, and regulation of apoptotic process (Supplementary Table 4).

Despite having fewer mapped reads (and thus lower power of detection), many more differentially expressed genes were identified in the symbiont transcriptome (89 genes, FDR = 0.05, Supplementary Table 5), consistent with a greater separation among pCO₂ treatments in the principal components analysis (Figure 2). A LRT indicated that the best-fit model was for symbiont gene expression was one that included only the effects of CO₂ (89 genes, FDR = 0.05), with no T or interaction effects (0 differentially expressed genes for either T or TxpCO₂, FDR = 0.05).

When the freshly-released larvae were included in the PCA (Supplementary Figure 1), the first principal component separated the freshly-released larvae from the eight experimental samples. This axis accounted for 64 and 52% of the variance among samples for the coral and symbiont, respectively, with 5,772 transcripts being differentially expressed between time points in the coral and 796 transcripts being differentially expressed between time points in the symbionts. However, there was no evidence that any of the experimental treatments resulted in a gene expression pattern that was more similar to the freshly-released larvae than to any other treatment.

For the coral transcriptome, only the Molecular Function category included enriched gene ontologies (Figure 3, Supplementary Table 6). Here, gene expression affected by the interaction of pCO₂ and temperature was enriched for enzyme and binding activities, which were upregulated in the ACAT and HCHT treatments, relative to the others, while structural molecule and oxoreductase activities were downregulated in the ACAT and HCHT treatments.

Several GO categories were functionally enriched for symbiont genes that were up- and down-regulated in response to High-pCO₂ (Figure 4, Supplementary Table 7). Under the Biological Function class, the up-regulated genes were in metabolic process and ion transport categories, and the down-regulated genes had biosynthetic, metabolic process, and RNA modification categories. Under the Molecular Function class, up-regulated genes were largely in transporter and binding categories. Down-regulated genes were enriched for several enzyme-related categories (i.e. isomerases, peptidases, hydrolases). Under the Cellular Component class, up-regulated genes included cellular complexes, and down-regulated genes included ribosomal subunit complex and photosystem.

Blast searches of the symbiont transcriptome against databases of Symbiodinium 28S nuclear and 23S chloroplast rDNA genes identified the symbiont clade present in coral larvae from this study as D1 (Pochon et al., 2006). Symbiont genes comp97910 and comp84390 had 99% and 96% sequence identity to clade D1 28S Genbank sequence AJ311948.1 and clade D1 23s Genbank sequence JN558001.1, respectively.

When examining compartments of the holobiont separately, conditions of OA and temperature had little effect on the gene expression profile of the coral animal, given the small number of differentially expressed genes and the greater power to detect differential expression given the higher number of mapped reads. These results provide insight for previous work showing that larvae released on the peak day of the spawning period, as were used in this study, are often larger in size with higher protein content and have different response phenotypes to high-pCO₂, high-temperature conditions (i.e., changes in oxygen consumption rates, changes in lipid utilization) than larvae released on other days (Cumbo et al., 2012, 2013a; Rivest and Hofmann, 2014, 2015). The ability of these larvae to tolerate high temperature, high pCO₂ conditions may be rooted at different levels of control other than gene expression: using existing enzymatic machinery or post-translational modification. For example, different oxygen consumption rates of P. damicornis larvae have been observed under different pCO₂ levels despite the same maximum activity of aerobic metabolism enzymes (Rivest and Hofmann, 2014). Indeed, in Acropora millepora embryos, key components of DNA/RNA metabolism are down-regulated while heat shock proteins are upregulated in response to elevated temperatures, but responses to long-term exposures (>24 h) differ (Rodriguez-Lanetty et al., 2009; Meyer et al., 2011). Larvae released close to the peak day, which are generally larger and contain more symbionts and protein that were vertically transmitted from the parent colony (Cumbo et al., 2012), may also receive more mRNA, proteins, and hormones that enhance physiological plasticity immediately after release (McCormick, 1999; Place and Smith, 2012). Future studies should pair measurements of larval gene expression with enzyme and whole-organism measurements of phenotype to provide a comparative investigation of mechanisms of plasticity.

The subdued response of the coral animal to multistressor conditions is unique and may inform its relationship with the symbiont. In contrast to previous studies on the animal compartment of other species (Voolstra et al., 2009; Moya et al., 2015), here there was no evidence that exposure to either increased temperature or increased pCO₂ led to upregulation of a generalized cellular stress response in P. damicornis larvae. We also observed no transcriptomic evidence for either depressed metabolism (as in Voolstra et al., 2009; Moya et al., 2012) or elevated metabolism (as in Davies et al., 2016) in response to high pCO₂. These results support the above hypothesis that the coral larvae can tolerate high temperature and high pCO₂ using mechanisms other than changes in gene expression. Additionally, the coral's subdued response in gene expression may indicate a lack of stress response on the part of the symbionts. Many previous studies have shown that photosynthetic dysfunction of symbionts during thermal stress produces reactive oxygen species (ROS), which elicit a stress response in the coral animal in the form of upregulation of antioxidant genes (Griffin et al., 2006; Merle et al., 2007; DeSalvo et al., 2008; Richier et al., 2008; Sunagawa et al., 2008). Our observed lack of change in expression in antioxidant genes may indicate that the symbionts in the coral larvae are not stressed by the elevated temperatures imposed during the experiment and are not producing a greater number of ROS. This hypothesis is supported by Rodriguez-Lanetty et al. (2009), which suggested that ROS from symbionts are required to stimulate upregulation of host oxidative stress genes.

Ours is the first study to compare transcriptome-wide expression patterns between host and symbiont compartments in coral larvae in response to temperature and pCO₂ conditions. Although other studies have examined the expression patterns of candidate genes in larvae (e.g., Putnam et al., 2013), the full transcriptome provides a global view of all genes involved in a particular stress response. Previous studies examining gene expression in animal and symbiont compartments of adult coral provide an interesting comparison to our results with larvae. For adult Acropora millepora near CO₂ seeps in Papua New Guinea, expression changes in Symbiodinium were also greater in response to pCO₂ than in the coral host (Kenkel et al., 2017). In adult A. hyacinthus, heat stress did not change symbiont transcriptomes (Barshis et al., 2014) but did affect expression of hundreds of host coral genes (Barshis et al., 2013). Fold change differences in Symbiodinium gene expression were smaller than for coral host genes in A. aspera following heat exposure (Leggat et al., 2011). The relative role of symbiont physiology and performance, as influenced by gene expression, may differ among life stages and species of coral. While host maternal effects and host plasticity at levels other than transcription are likely important, our results indicate that the symbiont may be playing a relatively large role in shaping the overall performance of symbiotic larvae post-release.

The gene expression profile in the symbionts of P. damicornis larvae changed dramatically in response to conditions of OA and warming. Parsing out the symbiont-specific gene expression response is important in the effort to accurately predict effects of global ocean change on coral-dinoflagellate symbioses, as Symbiodinium and other dinoflagellates can display muted or dramatic responses in gene expression to elevated temperature and other stressors (Rosic et al., 2010, 2011; Guo and Ki, 2012). Our results—greater expression changes in Symbiodinium vs. host—agree with a recent study in Papua New Guinea using adult A. millepora (Kenkel et al., 2017). As ocean conditions continue to change, coral animal populations will need to acclimatize or locally adapt to changes in environmental conditions (direct effects on animal host) as well as changes in symbiont performance (Kenkel et al., 2017) and in composition of hosted Symbiodinium genotypes/clades/types (DeSalvo et al., 2010).

Although symbiont density was significantly affected by the interaction of pCO₂ and temperature, symbiont gene expression was only affected by pCO₂. One of the most prominent effects of high pCO₂ observed in this study was the upregulation of ion transport genes in symbionts. Other studies of gene expression responses of coral holobionts to elevated pCO₂ have observed increased expression of ion transporters in high pCO₂ treatments (Kaniewska et al., 2012; Vidal-Dupiol et al., 2013; Davies et al., 2016). However, most of these previous studies have observed changes in ion transporter expression on the part of the coral host and have hypothesized that these changes in gene regulation allowed the host to maintain conditions necessary for calcification. Our results represent the first observation (of which we are aware) of altered ion transporter expression on the part of the symbiont. Observations of isolated Symbiodinium have shown that compatible intracellular ion composition with the host improves carbon assimilation by zooxanthellae (Seibt and Schlichter, 2001). Our observations suggest that symbiont gene regulation may be playing an important role in the maintenance of holobiont homeostasis. However, it is also important to note that ion regulation is often energetically costly (Seibt and Schlichter, 2001), so if symbionts are required to increase active ion regulation under high pCO₂ as part of physiological maintenance, it may come at the cost of long-term carbon assimilation. Alternatively, increased ion regulation under high pCO₂, may reflect the higher capacity of symbionts to enhance carbon metabolism under increased carbon availability.

The larvae in this experiment contained Clade D1 Symbiodinium, known for the high thermotolerance of its symbiosis with corals (Stat and Gates, 2011). Although symbiont density varied by temperature as a function of pCO₂ level, the lack of significant differential expression of symbiont genes among temperature levels supports the idea that the symbionts can tolerate the range of temperatures imposed with minimal adjustment from the level of gene expression. In Moorea, Pocilloporid corals are symbiotically flexible, known to harbor 5 symbiont types from clades A, C, and D (Putnam et al., 2012). While Clade D1 has a global distribution in tropical corals, it is relatively rare and is found in higher abundance in reefs exposed to high sea surface temperatures (Stat and Gates, 2011). The exclusive presence of Clade D1 in the larvae from 8 P. damicornis colonies suggests that these colonies may have a history of exposure to elevated temperatures in their shallow fringing reef habitat. Acclimatization to an environmental history of elevated temperature exposure may also explain the muted gene expression response in the coral animal. Additionally, the clade or genotype of symbiont harbored by the larvae in this study may have influenced the gene expression response of the coral, as previous work with adult Montipora faveolata showed that host gene expression profiles were significantly distinct based on the clade of symbiont hosted, above and beyond differential effects of thermal stress on gene expression (DeSalvo et al., 2010). Consequently, P. damicornis larvae harboring other symbiont clades may have produced different responses to conditions of high temperature and pCO₂. Indeed, the magnitude of host gene expression response to thermal stress differs between corals hosting different symbiont clades (DeSalvo et al., 2010).

Our observation of the plethora of differentially expressed genes in Symbiodinium as compared to the coral host offers important insight for future efforts to examine the response of the coral holobiont to future ocean conditions. Measurements of symbiont density alone are not enough to adequately describe or predict the functional change of symbiosis, as symbiont density and gene expression were affected by different environmental factors in our study. Furthermore, while a change in symbiont density would affect the coral host through changes in photosynthate delivery, ROS production, etc. (Cunning and Baker, 2014), clearly, the physiological function of symbionts, as controlled by their gene expression, is responsive to external seawater conditions. Corals whose symbiont densities do not change in response to environmental stress may still experience a fundamental functional shift in their symbiotic relationships.

Our study on genome-wide changes in expression of P. damicornis larvae in response to multiple stressors provides new insight for functioning of the coral-algal symbiosis under high-pCO₂, high-temperature conditions as well as the role of the gene expression patterns of each partner that correlate with responses observed at the level of the whole organism. Firstly, our results indicate that many of the whole-organism responses previously observed for P. damicornis larvae in scenarios of OA and warming may reflect physiological plasticity at levels of control other than gene expression. For example, multiple studies have observed significant changes in oxygen consumption rates in P. damicornis larvae in response to elevated temperature and pCO₂ (Cumbo et al., 2013a,b; Putnam et al., 2013; Rivest and Hofmann, 2014). However, in this study, coral genes involved in metabolic processes were not among the few differentially expressed genes. Additionally, Rivest and Hofmann (2014) found that citrate synthase activity did not differ among larvae cultured under combinations of temperature and OA conditions. Taken together, these results suggest that plasticity in metabolic rate may be accomplished through changes in fluxes of pathways using existing biochemical machinery, rather than the synthesis of new enzymes, at least for the population and seawater conditions evaluated.

Much of the work to date on responses of coral to environmental conditions has measured biological and physiological traits at the level of the whole-organism holobiont. For example, rates of growth/calcification and respiration are measured for a nubbin, colony, or intact larva. These measurements are made on the performance of the holobiont, and it is assumed that both coral and Symbiodinium contributed in some way to the response measured. The results of this study reveal strong potential merit in examining the gene expression of both partners, rather than just focusing on the gene expression of the coral animal host as a way to inform control of broader biological performance. While symbiont density and clade/genotype are known to cause differences in whole-organism performance like calcification and oxygen fluxes (e.g., Jokiel and Coles, 1977; Goreau and Macfarlane, 1990; Goulet et al., 2005), our results indicate that the performance of an existing symbiont population can change significantly with potential ramifications for host physiology and holobiont performance. Future work is needed to link differential expression in host and symbiont compartments with changes in processes of larval or colony biology.

Finally, it is clear that the responses of symbionts to future ocean conditions could play a large role in shaping success of coral larval stages. The relative role of symbiont function could depend on the clade/genotype hosted. Robust differential gene expression was observed in Symbiodinium clade D1 in P. damicornis larvae in response to experimental treatments. Pocilliopora damicornis is able to host other clades of Symbiodinium (Putnam et al., 2012), so corals may respond differently to environmental conditions based on unique transcriptomic responses of their symbiont clades. If the clade or genotype of symbiont is an important factor in holobiont function (through transcriptomic responses), then the brooding reproductive mode (vs. broadcast spawning) may be advantageous, as it allows parents to curate the symbionts the offspring receive. Symbiont clades that performed well and proliferated under the environmental conditions experienced by the parent can be vertically transmitted to brooded larvae and may allow them to perform better under those same environmental conditions immediately upon release.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest. The reviewer SF and handling Editor declared their shared affiliation.