Life-Long Coral Skeletal Acclimatization at CO2 Vents in Papua New Guinea Reveals Species- and Environment-Specic Effects

The responses of corals and other marine calcifying organisms to ocean acidication (OA) are variable and span from no effect to severe responses. Here we investigated the effect of long-term exposure to OA on skeletal parameters of four tropical zooxanthellate corals living at two CO 2 vents in Papua New Guinea, namely in Dobu and Upa Upasina. The skeletal porosity of Galaxea fascicularis, Acropora millepora, and Pocillopora damicornis was higher (from 17% to 38%, depending on the species) at the seep site compared to the control only at Upa Upasina. Massive Porites showed no differences at any of the locations. Pocillopora damicornis also showed a ~ 7% decrease of micro-density and an increase of the volume fraction of the larger pores, a decrease of the intraskeletal organic matrix content with an increase of the intraskeletal water content, and no variation in the organic matrix related strain and crystallite size. The fact that the skeletal parameters varied only at one of the two seep sites suggests that other local environmental conditions interact with OA to modify the coral skeletal parameters. This might also contribute to explain the great deal of responses to OA reported for corals and other marine calcifying organisms.


Introduction
Tropical coral reefs support the livelihoods of hundreds of millions of people around the world, harbor 25% of all marine species, and protect thousands of kilometers of shoreline from waves and storms 1 . However, coral reefs face a wide and intensifying array of threats deriving from pollution and overexploitation on centennial scales which is leading to a decline in their health 2 . In addition, global climate change compounds these threats in multiple ways. Declines in seawater pH and associated decreases in carbonate ion concentration driven by ocean acidi cation (OA) is projected to have profound implications for marine calci ers, as carbonate ions are essential for biotic calci cation 3 . Many biological features may affect coral responses to OA, such as colony morphology, skeletal mineralogy and structure, body size, tissue thickness, symbiont types, and/or the mechanisms of nutrient acquisition 4 . Moreover, the discrepancy among responses could derive from different experimental designs and analytical methods (e.g., addition of acid vs CO 2 bubbling to mimic OA), co-limiting environmental conditions (e.g., temperature, light intensity, ow, feeding, etc.), and exposure times (days to months or even life times) 5 .
Most studies to date, both under controlled conditions in aquaria and under natural conditions in the eld (e.g., CO 2 vents), support predictions of decreased rates of calci cation and increased rates of dissolution and bioerosion as seawater pH decreases 6 . However, coral calci cation rates recorded for several decades within skeletal cores indicate there hasn't been a constant decline as ocean pH decreased and temperatures warmed throughout the 20th century. On the contrary, at some locations calci cations rates have remained stable and in others they have increased over this time period 7,8 . Even where declines in calci cation have occurred, many other factors such as ocean warming, sea level rise, changes in surface ocean productivity, as well as many localized anthropogenic disturbances that co-occur with OA and also in uence coral growth obscure our ability to attribute such changes solely to OA 9 .
Most of the available knowledge about OA effects on marine organisms derives from short-term laboratory or mesocosm experiments on isolated organisms 10 , which can substantially underestimate full organism acclimatization 11 . In fact, taxa that result apparently unaffected by high CO 2 under controlled conditions may be: 1) vulnerable in the long-term 12 , 2) affected during life stages that were not considered during the experiment 13 , or 3) be indirectly affected by OA-driven ecological changes (e.g., food webs, competition, diseases and/or community structures, habitat properties such as microbial surface bio lms) 14 . Likewise, other taxa that respond negatively to OA under controlled conditions may be capable of acclimatizing in the longer term. Thus, eld experiments, where organisms are naturally exposed to OA for their entire life, as found around submarine CO 2 vents, could provide important new insights. However, vent systems are not perfect predictors of future ocean ecology owing to temporal variability in pH, spatial proximity of populations unaffected by acidi cation, and the unknown effects of other changing parameters (e.g., temperature, currents) 15 . Nonetheless, vents acidify sea water on su ciently large spatial and temporal scales to integrate ecosystem processes such as reproduction, competition and predation 16 . Field-based studies conducted at volcanic CO 2 seeps in Italy [16][17][18] , Japan 19 , Mexico 20 , and Papua New Guinea (PNG) 14 provide a unique opportunity to investigate long-term effects of OA on marine ecosystems that have been naturally exposed to chronic low pH and concomitant altered carbonate chemistry parameters for years/decades. These studies have already demonstrated substantial changes in community structure and functional biodiversity 21 of benthic species, as well as an array of responses to OA spanning from sharp decrease to no effect on calci cation rate 22 . Studies conducted on corals at volcanic CO 2 vents in PNG have supported the mixed effects observed in laboratory experiments 14,23 . Hard coral cover is similar at acidi ed and control sites (33% versus 31%). However, the cover of massive Porites corals doubled under OA, whereas the cover of more structurally complex corals is reduced by one third 23 . Some species are signi cantly less common or even absent under OA. For instance, while the coverage of Pocillopora damicornis decreases by 43% in acidi ed sites, in situ growth measurements have found small differences in linear extension rate 14 , but large differences in recruitment success 24 . Population reductions in situ, combined with observations of negative physiological impacts, including declines in calci cation under OA, strongly suggest that low pH imposes selection pressure on less resilient taxa within the PNG system 22 .
The aim of this study was to assess the effects of long-term exposure to OA on the skeletal parameters (micro-density, porosity, bulk density) of four tropical zooxanthellate coral species Galaxea fascicularis (Linnaeus, 1767), Acropora millepora (Ehrenberg, 1834), massive Porites Link, 1807, and P. damicornis (Linnaeus, 1758), living at PNG CO 2 vents 14 . The study was conducted at two locations in Milne Bay Province, PNG, namely Upa Upasina and Dobu ( Fig. 1). At each location, corals were collected from a shallow water (1-5 m) cool volcanic CO 2 vent (seep hereafter) and a control site. The algal endosymbionts of these corals did not differ between seep and control sites, nor between the two seep locations 25 . Seep and control sites have been characterized in detail 14,23 .

Materials And Methods
Page 5/18

Study sites and coral sampling
The study was conducted at two shallow-water (1-5 m) volcanic CO 2 seeps at ambient temperature and adjacent control sites at Milne Bay Province, PNG, namely Dobu and Upa Upasina (Fig. 1). Almost pure CO 2 (~ 99%) has been streaming from the seabed for an unknown period of time (con rmed for approximately 70 years, but possibly much longer) 14 , resulting in localized acidi ed conditions. Environmental parameters were measured across a 4-year period (2010-2013) at 1-5 m depth in both control and seep sites at Dobu and Upa Upasina, as previously reported 14,23 . Two-cm coral fragments were collected at 1-5 m depth from adult colonies of P. damicornis, G. fascicularis, A. millepora, and Massive Porites at control and seep sites in Dobu and Upa Upasina in August 2010 (N = 6-15 fragments per site, each fragment from a different colony; Supplementary Table S1) in Dobu and Upa Upasina 25 . Tissue from the coral fragments was totally removed following standardized protocols 26 .

Skeletal parameters determination
The skeletal paramaters of 192 fragments from the control and seep sites at Dobu and at Upa Upasina (6-15 for each site) were obtained by buoyant weight measurements with a hydrostatic balance (Ohaus Explorer Pro balance ± 0.0001 g) equipped with a density determination kit. After determining the dry mass the fragments were placed inside a dryer chamber connected to a vacuum pump to evacuate air and water from the pores. After 3 h, water was gently introduced to fully saturate the samples which were then weighed in air. The buoyant weight was then obtained by applying the density determination kit, and the skeletal parameters were calculated by means of standard calculations (details in Supplementary Methods; 26 ).

Time-Domain Nuclear Magnetic Resonance for pore size distribution determination
This technique was used to investigate the 'pore-size' distribution of P. damicornis coral skeletal fragments from each control and seep site (details in Supplementary Methods). Each coral fragment was saturated with distilled water, then removed from the water and placed on a wet paper to dry the excess of water on its surface. Then, every fragment was put inside a glass tube, sealed and measured. A home-built relaxometer based on an electromagnet JEOL C-60 operating at 0.5 T with a coil ≈ 8 mm in diameter, and equipped with a Spinmaster portable console (Stelar, Mede, Pavia, Italy) was used. The Carr-Purcell-Meiboom-Gill (CPMG) sequence with 200 µs echo time was used for T 2 measurements. The measured multi-exponential relaxation curves, affected by unavoidable measurement noise, were transformed into T 2 distributions by the algorithm UPEN (Uniform-Penalty inversion algorithm) 27 , implemented in UpenWin. The ratio between the NMR signal under a particular portion of the distribution and the total NMR signal will correspond to the ratio of the volume of the pores with a particular pore size to the total pore volume, giving, as an example, the macro-scale pore fraction. A total of 60 fragments from the control and seep sites at Dobu (15 for each site) and at Upa Upasina (15 for each site) were analyzed.

Thermogravimetric Analysis for organic matrix content determination
Thermal gravimetric measurements were performed on using a TA Instruments thermobalance model SDT-Q600 with 0.1 µg of balance sensitivity. Powdered samples (5 to 10 mg), held in alumina pans, were heated under a linear gradient from ambient (ca. 20 °C) up to 600 °C with an an isotherm at 120 °C for 5 min to remove the adsorbed water; heating rate: 10 °C/min under an N 2 atmosphere, with ux xed to 100 ml/min. Two main weight loss regimes were identi ed: a rst one in a range around 125-250°C (related to the loss of structured water molecules) followed by another thermal region between 250 °C and 470° C (generally associated with organic matrix pyrolysis) 28 . A total of 32 fragments from the control and seep sites at Dobu (10 for each site) and at Upa Upasina (6 for each site) were analyzed. Before the analysis the baseline and temperature were calibrated.

Synchrotron high-resolution X-ray powder diffraction
The coral fragments were measured at the ID22 beamline of the European Synchrotron Radiation Facility (Grenoble, France) using a wavelength of 0.4 Å (details in Supplementary Methods). The fragments were air-dried and ground with an agate mortar and pestle to a ne powder which was then loaded into borosilicate glass capillaries of 0.7-1 mm in diameter and measured at room temperature and after ex-situ heating at 300 °C for 2 h. A Rietveld re nement was used to calculate the unit-cell parameters of the diffraction pattern pro les. The line pro le analysis was applied to a speci c diffraction peak to obtain the coherence length (nm) along various crystallographic directions, which was achieved by tting the pro le to a Voigt function and deconvoluting the Lorentzian and Gaussian widths. Analyses were conducted on fragments of P. damicornis from Upa Upasina in the control (N = 3) and seep (N = 3) site.

Statistical analyses
Skeletal parameters (including macro-scale pore fraction, water content and organic matrix content analysed in P. damicornis) were compared between control and seep sites using the non-parametric Kruskal-Wallis test (χ 2 ), due to deviations from parametric ANOVA assumption being veri ed (Normality: Shapiro-Wilk's test; equal variance: Bartlett's test). Statistical analyses were performed using SPSS 20.0.
Data visualization, and graphics were obtained with the ggplot2 R package in R 29 . Statistical differences were accepted when p < 0.05.
Permutation multivariate analysis of variance (PERMANOVA) were perfomed using PRIMER v6 30 and based on Euclidean distances (999 permutation) to test for (i) variations of environmental parameters amongst locations and sites, and across seasons; (ii) variations of skeletal parameters amongst locations, sites, and species. When the main tests revealed statistical differences (p < 0.05), PERMANOVA pairwise comparisons were carried out. Distance-based redundancy linear modeling (DISTLM) with a test of marginality in PRIMER was also performed to account for the contributions environmental parameters in explaining the total observed variance in the skeletal parameter datasets of each species. DISTLM used the BEST selection procedure and adjusted R 2 selection criteria. BEST BIOENV routine in PRIMER 6 was also carried out to explain the observed patterns of nano-, mico-and macro-scale paraemetrs of P. damicornis (999 permutations). BEST BIOENV method performs permutation tests on environmental variables for determining which subset of variables produced the highest correlation with biological data. The environmental and biological data matrices were used to calculates a Spearman's rank correlation coe cient between these two matrices.

Environmental parameters
The complete dataset of environmental parameters was analyzed to test for the differences between sampling locations, and, within each location, for control vs seep sites. Effects of seasonality were also considered. PERMANOVA analyses (Supplementary  Fig. S1).
Controls at the two locations displayed lower pH and higher DIC, total alkalinity, salinity, and temperature at Upa Upasina compared to Dobu, while pCO 2 and Ω AR did not change. Seep sites at the two locations showed higher pH, Ω AR , total alkalinity, and temperature, and lower pCO 2 , DIC at Upa Upasina compared to Dobu, while salinity did not change (Table 1).
DISTLM analysis and the related tests of marginality indicate that carbonate chemistry parameters (i.e., pCO 2 , pH TS , and Ω AR ) explained most of the variance observed in the environmental dataset ( Fig. 2A).
Indeed, PCO analyses show that these parameters are those that mostly explained the differences between controls and seeps within the two locations (Fig. 2B), while temperature and DIC seem to explain the differences between the two locations (Fig. 2B).

Skeletal parameters in corals from control and seep sites at the two locations
Variations of the skeletal parameters bulk density, micro-density, and porosity are reported in Fig. 3 and in Supplementary Table S1. In all four species, none of these skeletal parameters were affected by CO 2 at Dobu (p > 0.05). At Upa Upasina, A. millepora and G. fascicularis showed signi cantly higher porosity (+ 38% and + 25%, respectively) and consequently lower bulk density (− 13% and − 15%, respectively) in the seep compared to the control site (Fig. 3A,B; Supplementary Table S1). At the same location, P. damicornis showed signi cant variation in all investigated skeletal parameters, with lower micro-density (− 7%), higher porosity (+ 17%), and lower bulk density (− 12%) at the seep compared to the control site ( Fig. 3C; Supplementary Table S1). Skeletal parameters in massive Porites did not change between seep and control sites in either locations ( Fig. 3D; Supplementary Table S1). PERMANOVA analyses con rmed Site as a signi cant factor determining the observed variations in porosity and bulk density, while also accounting for the signi cant species-speci c variations in all skeletal parameters (Supplementary Table S3).
Furthermore, porosity and bulk density also showed a signi cant interaction between Location and Site, while micro-density showed a signi cant interaction between Site and Species (Supplementary Table S3).
The decreased micro-density in P. damicornis between seep and control sites at Upa Upasina was further explored to assess additional macroscale and microscale skeletal changes. Speci cally, Time-Domain Nuclear Magnetic Resonance (TD-NMR), Thermogravimetric Analysis (TGA), and synchrotron highresolution powder X-ray diffraction (HRPXRD) analyses were performed. NMR measurements were performed on P. damicornis also from Dobu. In the T 2 distributions obtained by NMR (Supplementary Fig.   S2), it was possible to identify a cut-off at 3 ms. This allowed to divide the pores containing water into two classes, distinguishing the smaller pores (smaller volumes, estimated pore sizes < 1 µm) from the remaining larger ones (larger volumes, estimated pore sizes > 1 µm). For the sake of simplicity, the two classes was named micro-scale and macro-scale pores 12,31 . The macro-scale pore volume fraction showed a 7% increase (p < 0.05) at the seep site compared to the control site only in Upa Upasina ( Fig. 4; Supplementary Table S4). Figure 4 and Supplementary Table S5 summarize the intraskeletal organic matrix (OM) and water content (% mass loss) evaluated by TGA. In Dobu no signi cant difference between the control and seep sites was found neither for OM nor for water content (p > 0.05). In Upa Upasina, OM content showed a signi cant decrease in the seep compared to the control (p < 0.01), while a signi cant increase in water content was observed in the seep compared to the control (p < 0.01).
Three coral skeleton fragments of P. damicornis each from the control and seep sites in Upa Upasina were analysed by HRPXRD. All HRPXRD patterns were well indexed as aragonite and no additional diffraction peaks were detected. Then, they were re ned using the Rietveld method 32 and lattice parameters and strain (Supplementary Table S6), and microstructural data 33 , crystallite size, and microstrain (Supplementary Table S7), were calculated. No signi cant differences were found between the control and seep site. To test the in uence of the OM on the mineral strain, ex-situ heat treatment prior the HRPXRD measurements were performed. The heat treatment removes the OM effects on the strain 34 . The data showed that the OM induced a positive strain on aand c-axis and a negative one on the b-axis, but no signi cant differences were found between the control and seep site (Supplementary Table S6). We also measured crystallite size after the thermal annealing together with the transition to calcite (Supplementary   Tables S6 and S7). These latter parameters did not show any signi cant difference between the control and seep sites.

Relations between variations of skeletal parameters and environmental variables
The DISTLM analysis exploring possible sources of variation in skeletal parameters of the single species related to the environmental parameters con rmed that P. damicornis was the species in which the environmental variables showed the strongest in uence on skeletal parameters, as indicated by the higher proportion of variance explained, compared to the other two species (Supplementary Fig. S3). Massive Porites was excluded from the analysis since no signi cant changes in any of the measured skeletal parameters were observed (Supplementary Fig. S3). The carbonate chemistry parameters explained almost all signi cant variations of skeletal parameters between locations and sites (p < 0.05; Table 2; Supplementary Fig. S3). In P. damicornis, also salinity contributed to explain some of the variance (Supplementary Fig. S3). For this species, the in-detail BEST BIO/ENV analyses on the multi-scale skeletal parameters (Supplementary Table S8) showed that Ω AR and, in general, carbonate chemistry parameters, signi cantly correlated with most skeletal parameters. Moreover, salinity signi cantly correlated with changes in macro-scale pore volume fraction and intraskeletal water content, while temperature correlated with changes in intraskeletal OM content (p < 0.05; Supplementary Table S8).

Discussion
While measurable acidi cation of the tropical oceans has been underway for several decades now, detection and attribution of the effects of OA on reef-building corals has been challenging because multiple environmental changes, including ocean warming, are co-occurring with OA, impacting coral growth 35 .
This study investigated the effects of long-term exposure to elevated CO 2 on skeletal properties in tropical zooxanthellate corals living at CO 2 vents.
Of the four species investigated, massive Porites was the only species showing unmodi ed skeletal parameters between seep and control sites at both locations. This is in agreement with a previous investigation performed at the same locations showing that the cover of massive Porites was twice as high at the seeps compared to the control sites, while net calci cation and skeletal bulk density were similar 14 .
Increased Porites abundance and unchanged skeletal extension, skeletal bulk density, and net calci cation rates with decreasing pH was also observed along a natural pH gradient in Palau 36 . Our results here corroborate the hypothesis that massive Porites spp have the capacity to acclimatize and thrive under persistent exposure to ocean acidi cation whereas any other and especially the structurally more complex branching corals appear to be more sensitive 23 .
Similar to Mediterranean 12 and other tropical coral species 37 , G. fascicularis, A. millepora, and P. damicornis showed a signi cant increase in porosity and decrease in bulk density at reduced pH TS in Upa Upasina. Although the Dobu seep site showed signi cantly higher CO 2 levels, with consequently lower pH and Ω AR compared to the Upa Upasina seep, a signi cant increase of porosity and decrease of bulk density in the seep compared to the control site was observed only at Upa Upasina. This nding suggests that low pH is not the main driver and that other local environmental conditions not considered in this study interacted with OA to modify the skeletal parameters of these three coral species.
The species-speci c DISTLM analyses and the BEST BIO/ENV analyses performed on P. damicornis pointed out that carbonate chemistry related parameters, in particular Ω AR , contributed to most of the observed changes for all species. This nding is in agreement with several studies showing a strong dependence of calci cation on seawater aragonite saturation state 38 . Among the investigated species of the current study, only P. damicornis appears sensitive to the minor variation in salinity, which may indicate a trade-off between osmoregulatory and acid-base regulatory systems (controlling H + and HCO 3 − uptake from seawater), as observed in other calcifying marine organisms 39 . For example, responses of acid-base balance of crabs to both hypercapnia and changes in seawater salinity are related to iono-regulation because both homeostatic processes share the same mediators, such as H + and HCO 3 -transporters 40 . Pocillopora damicornis exposed to hyposalinity conditions for 24 h showed decreased photosynthetic and respiratory activities 41 . Thus, this species may be more vulnerable to low pH under uctuating salinities, a condition that is likely to occur in coastal environments in the face of global climate change 42 . Furthermore, P. damicornis was the only species displaying a variation in micro-density, with lower values at the seep site compared to the control. Micro-density, which represents the mass per unit volume of the biogenic calcium carbonate composing the skeleton 43 , depends on the mineral composition of the skeleton and content of intraskeletal OM and water 26 . The additional analyses of macro-and micro-scale parameters performed in this species revealed an increase in macro-scale pore volume fraction and intraskeletal water content and a decrease in OM, and eventually strong linked water 28 . In particular, the observed increase in intraskeletal seawater content at the Upa Upasina seep can partially justify the observed decrease in skeletal micro-density. Changes in OM and water content with pH reduction are reported showing either increased content in the tropical Stylophora pistillata kept in aquaria at pH 7.2 for approximately one year 28,37 , or no variation in the temperate B. europaea naturally living at a CO 2 vent 44 .
However, it must be noted that different methodologies were used in the two studies. Transcriptomic data in aquaria experiments conducted on A. millepora, P. damicornis, and S. pistillata show that several genes encoding OM protein are up-regulated under reduced pH 45 . In particular, P. damicornis exposed to pH 7.8, 7.4 and 7.2 in aquaria for 3 weeks showed a 4 to 70-fold increase in up-regulation of genes encoding skeleton organic matrix proteins at all pH treatments 46 .
The decrease in intra-skeletal OM in the samples from the seep site was not associated with a signi cant change in the strain, micro-strain, and crystallite size. These observations may indicate that the amount of intra-crystallite OM does not change, in agreement with the fact that the crystallite sizes after the thermal annealing are the same for samples from the control and the seep sites. Thus, the observed decrease in OM is likely associated with a decrease in the inter-crystallite OM. In addition, the stability of aragonite through the transition to calcite did not show a signi cant difference between control and seep samples, as well as the lattice parameters of the calcite formed after thermal annealing. The crystallographic features of aragonite from coral skeletons have been previously investigated 47 47 . In addition to the decreased OM, P. damicornis at the seep site of Upa Upasina also showed an increase of intra-skeletal seawater content, which could be partially related to the observed decrease of skeletal microdensity, but whose role in coral phenotypic plasticity to OA still has to be investigated. Again, these different responses seem species and environmental-speci c.

Conclusions
This multi-species study showed that the combination of different environmental conditions can have a stronger effect on macro-scale skeletal parameters than average low pH values alone. Our ndings showed a common phenotypic response among three zooxanthellate corals which all displayed a more porous skeletal phenotype under OA but also highlighted that OA is not always the main driver and that other local environmental conditions likely interacted with OA to determine the observed responses. Additionally, these skeletal macromorphological adjustments in P. damicornis did not affect the measured crystallite features, suggesting that the fundamental structural components produced by the biomineralization process might be substantially unaffected by increased acidi cation 12,49 . Nonetheless, the extra porous phenotype here described in the branching coral species may render structurally complex corals more vulnerable to damage and bio-erosion under climate change compared to massive growth forms 37,50 , which in the future may lead to a weakening of the reef framework and subsequent degradation of the complex coral reef ecosystem. More generally, our ndings highlight the importance of using a multi-parameter and multispecies analysis when investigating the vulnerability of coral species to OA, to understand what induces corals in a certain environment to acclimatize, and whether other species under the same conditions have the same capacity to adjust to future changes. These species and environmental speci c differences highlighted in the present study contribute to explain the large range of responses to OA reported for corals and other marine calcifying organisms.

Declarations Data Availability
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.    Figure 4 Macro-scale pore volume fraction (in the gure simply macro pore fraction) and intraskeletal OM and water content for P. damicornis from control and seep sites at Dobu and Upa Upasina. (A) TD-NMR measurement of macro-scale pore volume fraction. (B) TGA measurements of microscale parameters, namely intraskeletal organic matrix (OM), intraskeletal water content, and total (the sum of OM and water). The box indicates the 25th and 75th percentiles and the line within the box marks the median. Whisker length is