Remote coral reefs can sustain high growth potential and may match future sea-level trends

Climate-induced disturbances are contributing to rapid, global-scale changes in coral reef ecology. As a consequence, reef carbonate budgets are declining, threatening reef growth potential and thus capacity to track rising sea-levels. Whether disturbed reefs can recover their growth potential and how rapidly, are thus critical research questions. Here we address these questions by measuring the carbonate budgets of 28 reefs across the Chagos Archipelago (Indian Ocean) which, while geographically remote and largely isolated from compounding human impacts, experienced severe (>90%) coral mortality during the 1998 warming event. Coral communities on most reefs recovered rapidly and we show that carbonate budgets in 2015 average +3.7 G (G = kg CaCO3 m−2 yr−1). Most significantly the production rates on Acropora-dominated reefs, the corals most severely impacted in 1998, averaged +8.4 G by 2015, comparable with estimates under pre-human (Holocene) disturbance conditions. These positive budgets are reflected in high reef growth rates (4.2 mm yr−1) on Acropora-dominated reefs, demonstrating that carbonate budgets on these remote reefs have recovered rapidly from major climate-driven disturbances. Critically, these reefs retain the capacity to grow at rates exceeding measured regional mid-late Holocene and 20th century sea-level rise, and close to IPCC sea-level rise projections through to 2100.


Study area
The Chagos archipelago is located in the central Indian Ocean (Lat: 5-7 o S, Long: 71-73 o E; and comprises five main atolls (see main Fig. 1C): Great Chagos Bank (~18,000 km 2 ) which is mostly submerged, but which has eight reef islands developed on its western and northern rims, Peros Banhos (~463 km 2 ) which has 24 reef islands around its rim, Salomon (~38 km 2 ) with 8 rim islands; Egmont (~40 km 2 ) with 2-3 islands on its rim (several islands having merged over the last decade or so), and Diego Garcia (~200 km 2 ) which has only 4 individual islands, but has the largest land area (2,733 ha) of any of the atolls. There are, in addition, a further ten submerged atolls and banks. Most islanded atolls have now been uninhabited for nearly 40 years, and thus the Chagos reefs are almost entirely free of direct anthropogenic impacts. Exceptions include the impacts of low (though poorly quantified) levels of illegal fishing on the outer atolls, impacts from recreational fishing around Diego Garcia (which has lower fish biomass than the other atolls 1 ), and the effects of the terrestrial military development on Diego Garcia. The Chagos reefs are thus amongst the remotest in the Indo-Pacific and, in April 2010, the entire archipelago and its associated exclusive economic zone (an area of ~640,000 km 2 ) was officially established as a no-take marine protected area. Impacts linked to climatic change and to broad scale oceanic and meteorological disturbances thus currently represent the only serious threats to reef health and ecosystem function.
Our ecological and carbonate budget surveys were conducted during March to April 2015 on 28 reefs across each of the five islanded atolls (numbers of survey sites in brackets): Diego Garcia (5); Peros Banhos (7), Salomon (6), Great Chagos Bank (8), and Egmont (1), as well as at one site on the submerged Blenheim Reef (see main Fig. 1C). The atolls and platforms of Chagos are characterised by marked differences in wave energy regimes around their margins (see main Fig. 1C), driven by the seasonally-shifting wind regime, with the predominant wind direction being from the south-east (see main Fig. 1B). Our site selection strategy around atolls was driven by a desire to survey sites on both the more sheltered (the south-west, western and northern margins of the atolls, and those on the more exposed margins (the north-east, east and south-eastern margins), as well as integrating sites that had been the focus of earlier ecological surveys. At each site, surveys were conducted at a depth of ~10 m i.e., a little above the upper shelf break, and with replicate transects established running parallel to the reef crest, with a spacing of ~5 m between transects. With only two exceptions we collected data along 4 replicate transects at each site (the exceptions being Middle Island (n = 5) and Cannon Point (n = 3) on Diego Garcia).

Wave exposure modelling
To enable us to classify these sites on the basis of their wave exposure regime, spatially explicit estimations of wave exposure were modelled as a function of wind speed and direction, and fetch length (i.e. the distance over open ocean that wind can travel in a specific direction unobstructed by land or reefs). To do this we followed the same protocols previously applied in reef environments 2,3 . Fetch lengths were calculated using the USGS model 4 which uses the procedure recommended by the Shoreline Protection Manual 5 . A binary raster representing the distribution of land masses and reef crests was generated using the outputs of the Millennium Coral Reef Mapping Project at a spatial resolution of 30 m 2 ref. 6 . Fetch lengths were generated for 16 compass directions (every 22.5⁰) by calculating the arithmetic mean of 5 radials spread at 3⁰ increments around the desired wind directions. 650 km was used as a maximum limit for fetch length as this is the distance required for maximal wave conditions 7 . Hourly wind measurements from 1973 to 2001 were obtained from Diego Garcia airport (n = 219,943) and used to calculate both the probability of wind blowing from each of the 16 compass directions and also the mean velocities for each direction. Fetch lengths were then converted into wave energies based on linear wave theory and using established equations ref2 . Based on model outputs (see main Fig. 1C) we have classified our sites into "exposed" (>1000 J m -3 ) or "sheltered" (<1000 J m -3 ), this division being based on a natural break in the rank order of the data across all sites. We refer to these groupings in the text for descriptive purposes.

Quantifying carbonate production and erosion rates
To quantify gross carbonate production and erosion and thus to determine net carbonate production (G, where G = kg CaCO3 m 2 yr -1 ) we used a version of the ReefBudget methodology of Perry et al. ref8 adapted and parameterised for the Indian Ocean. However, there are some importance differences: carbonate production by corals and coralline algae is calculated using geometric relationships derived from individual colony morphology, rather than calculated using rugosity at the level of the transect; and clinoid sponge bioerosion is calculated using published rates and the proportion of hard substrate under the transect line available for erosion. Sediment production by macro-bioeroders (urchins and parrotfish) is also estimated, but other aspects of sediment production and post-depositional lithification are not. Within each site, two different sets of transects were conducted, one set to census benthic carbonate producers and bioeroders, and one set to quantify parrotfish abundance and size.
Along each benthic assessment transect we measured the distance within each linear 1 m covered by each category of benthic cover beneath a 10 m guide line using a separate flexible tape. All overhangs, vertical surfaces and horizontal surfaces below the line were surveyed (i.e., if the guide line crossed over a table coral, the upper and lower surfaces of the coral, plus the benthos under the canopy, were recorded). The following groups were recorded: scleractinian corals to the genera and morphological level (e.g., Acropora branching, Porites massive etc.); crustose coralline algae (CCA) including CCA below macroalgal or soft coral cover; turf algae; fleshy macroalgae; non encrusting coralline algae (e.g., Halimeda sp., articulated coralline algae); sediment; bare substrate (e.g., granitic rock, limestone pavement); sediment; rubble; and other benthic organisms. Substrate rugosity was calculated as total reef surface divided by linear distance (a completely flat surface would therefore have a rugosity of 1).

Calculating carbonate production
In contrast to the original Reefbudget approach 8 , we used the morphology and size of individual coral colonies in combination with genera specific skeletal density (g cm -3 ) and linear growth rates (cm year -1 ) across each transect to estimates carbonate production rates in kg CaCO3 m -2 year -1 (where m 2 refers to the planar surface of the reef). Where possible, we used published growth rates and skeletal densities from the Western/Central Indian Ocean, but when these were not available means of published growth rates and densities for each coral genera were used instead. These data were then combined with geometric transformations based on colony morphology to give a growth rate for each colony for the area under the transect line (taking a transect line width of 1cm): massive colonies were assumed to be hemispherical in cross-section; encrusting, foliose and plating colonies, as well as colonies of crustose coralline algae (CCA) were assumed to be growing primarily at the edge of the colony (and at 10% of this growth rate across the remainder of the colony); for branching colonies, the proportion of the colony area of growing branch tips was assumed to be growing at published rates, and the remainder of the colony at 10% of these rates. The equations used in our calculations are thus as follow: Massive: Encrusting etc.: Where CPi = carbonate production for colony i, g = growth rate, x = surface length of colony, d = skeletal density and ca = proportion of colony that are growing axial branches.
Measuring the linear surface of growing tips on branching corals is time-consuming. Therefore, in order to calculate the amount of each colony growing as axial branch tips, we used data previously obtained on the size of branching colonies and the length of growing tips from 405 coral colonies in northern Mozambique (293 Acropora, 62 Pocillopora, 26 Porites, 24 other) (ref. 9). We conducted linear regressions between colony size and length of growing tips for each genera/morphology combination for which we had greater than 20 replicates in order to calculate ca. To calculate the production for a single transect over a year, the following equation was used: Where CPj is the total carbonate production of both corals and crustose coralline algae for transect j in kg CaCO3 year -1 . To estimate the production rate of the reef, we then used the following equation: Where Gprodj is the carbonate production rate of both corals and crustose coralline algae for transect j in kg CaCO3 m -2 year -1 , and l is the transect length in centimetres.

Micro-(endolithic) and Macro-(clionaid sponge, polychaetes, bivalves etc.) bioerosion
In general, macro-borer communities are less well-characterised in the Indo-Pacific than in the Caribbean. This is particularly true of clionaid sponges, which are generally cryptic and difficult to identify in the field, particularly for a non-expert, and are typically a less dominant part of the macroborer community than in the Caribbean. To this end, instead of conducting an intensive search of the substrate for clionaid sponges as described in Perry et al., (2012), we instead utilized published rates of total macrobioerosion measured at Indo-Pacific sites, alongside a census of substrate available for bioerosion from the benthic line-intercept transects. This comprises of all dead carbonate substrate available to bioeroding sponges, including that covered by macroalgae or algal turf and live coral cover and soft corals. While both live and soft corals can prevent settlement of most bioeroding sponges, live corals are often colonised by other bioeroders (particularly polychaete worms), while soft corals are often ephemeral, and removal of soft coral cover, either by predation or storms could allow settlement and establishment of bioeroders prior to regrowth. All substrate not available to bioeroders was excluded, and the following equation used:

Macrobioerosion = S*R*Ea
Where S is the percentage of surface area of the transect available for erosion, R is the rugosity of the transect, Em is the erosion rate of microborers, and Ea is the erosion rate of macroborers in kg m -2 yr -1 .

Urchin bioerosion
The main agents of echinoid bioerosion on reefs belong the family Diadematidae, (Diadema spp., and Echinothrix spp.) and the genera Echinometra, Echinostrephus and Eucidaris. We censused urchin abundance, size and species composition within 10 x 2 m belt transects along the main benthic transect lines at each site. These census data were then combined with carbonate ingestion rates from the literature and adjusted to account for the rugosity of the substrate. A variety of techniques have been used to estimate bioerosion rates in these urchin species: including CaCO3 content of the gut (e.g., ref. 10) or of their faecal pellets (e.g., ref. 11), both with or without estimations of reworked sediment, spine abrasion and gut turnover (e.g., ref. 12,13). This makes it difficult to compare bioerosion rates derived from different studies. However, evaluating the published data on erosion rates against test size across all urchin species suggest a relatively tightly correlated plot. SI Figure 1A shows the pooled bioerosion rates relative to test size for six species of urchins across 10 studies in the Indian and Pacific Oceans, with SI Figure 1B showing the rates for Echinometra mathaei and the Diadematidae. This allows us to calculate the erosion rate (kg urchin -1 year -1 ) for each individual urchin using one of the following equations:

Parrotfish bioerosion
To determine the species-size-life phase abundances of bioeroding parrotfish at each site we used underwater census data for Chagos collected in 2010 and 2012 ref. 22 , augmented by additional data collected in 2015. At each site four 50 m x 5 m belt transects were used (covering 1000 m 2 of reef). Time to complete transects was not held constant but adjusted according to the numbers of fish in the sampled areas. All fish surveys were completed by the same experienced observers (N.A.J.G., S.K.W or G.N.M). Biomass of individual fish was then calculated using estimated length data and length-weight relationships 23,24 ) and multiplying by abundance of the species or family of the fish. To calculate parrotfish bioerosion rates by each individual fish we then used a model based on total length and life phase to predict the bite rates (bites hr -1 ) for each species. For species for which data was not previously available, we have made use of new observations) 9 or have used data from similar sized species with the same feeding functional group. Daily bite numbers were calculated using diurnal feeding activity reported in Bellwood 25 . Bite volume data has been obtained from published studies 25,26 . Where no bite volume data exists, we have either measured the size of individual bites in situ using Vernier callipers to obtain width and length of bite, or obtained from published studies estimates of bite area (e.g., Ref 27). We have then used a bite depth of 0.1 mm to obtain a conservative estimate of bite volume for Scarus sp. Not all bites on the substratum remove material and we have assumed, following Bellwood & Choat 28 , that only bites that leave visible scars are eroding the substrate. In order to estimate the proportion of bites resulting in scars on the substratum for each size class, we extrapolated from published literature 26,29 and used the same method as above where data was missing. We used the following equation to calculate species specific erosion rates for the median value within each size class: Bioerosion rate (kg.ind -1 yr -1 ) = v.sprop.br.d*365 Where v is bite volume (cm 3 ), sprop is the proportion of bites leaving scars, br is bite rate (bites day -1 ) and d is substratum density (kg cm -3 ), here taken to be 1.49*10 -3

Statistical treatment of data
Variation in bioerosion, and in net and gross production of carbonate, were analysed with respect to reef, atoll and exposure using analysis of variance, where reef was nested within atoll. Data were plotted and tested for normality and homogeneity of variance prior to analyses and transformed to meet these assumptions when appropriate. Significant results were further investigated with a post-hoc Tukey test. The interaction between atoll and exposure was investigated by pooling data at the reef level. Differences in the influence of different variables on net carbonate budgets between reefs were examined using correlation-based principle components analysis (PCA). Data were log(x + 1) transformed and normalized to account for skewness and to ensure all metrics were on a common scale. We overlaid eigenvectors to identify the direction and contribution of different variables to the patterns. The PCA was conducted using the prcomp function in R 3.1.1.

Reef accretion rates and comparisons to recent and future sea level trends
To assess the accretion potential of reefs, and to explore their capacity to respond to future projected regional sea-level rise rates, we converted our net production rate estimates to potential accretion rates (mm yr -1 ). To do this we used an approach previously applied to Caribbean reefs 30 , but modified to also factor for variations in accumulating framework porosity as a function of between-site variations in reef community composition. Specifically, we estimated the maximum accretion potential of each reef as a function of the net carbonate production rate of the site (calculated as gross production less gross erosion rate) and assumed that a proportion of the bioeroded framework (that is converted to sediment) is also reincorporated back into the accumulating reef structure. This proportion was calculated as the sum of 50% of the parrotfish-derived sediment (as a highly mobile bioeroder which defecates randomly over the reef), as well as all sediment produced by urchins and by macroborer erosion. To keep our estimates conservative we worked on the assumption that only 50% of this bioerosional sediment yield is actually incorporated back into the reef (based on data from Hubbard 31 ), and excluded any sediment generation by other benthic sediment producers (e.g Halimeda -although this calcareous green algae was extremely rare within our survey sites). Finally, we made an allowance for variations in the porosity of the accumulating reef framework as follows: 30% for head and massive coral dominated assemblages, 70% for branched and tabular dominated assemblages, and 50% for mixed coral assemblages (based on data in Kinsey & Hopley 32 ). A loss factor to account for natural framework removal through physical processes was also included based on the data of Morgan and Kench 33 , such that we assumed that 20% of the annual framework produced was removed from more sheltered reef settings, and 50% from exposed settings.
To explore the accretion potential of the Chagos reefs in relation to sea level we compared our potential accretion rate estimates to both measured recent past and future projected sea-level rise rates. Recent sealevel trends around Chagos, based on tide-gauge records, are difficult to precisely constrain due to an apparent tide-gauge offset in the Diego Garcia record around 2002 (ref. 34). However, research that has integrated these records with satellite altimetry data suggests rise rates of ~2.2 mm yr -1 around Diego Garcia and 1.7 mm yr -1 around Peros Banhos and Salomon between 1993 to 2011 (ref. 34). These rates are above those reported (~1 mm yr -1 ) from the central Indian Ocean region over the period 1950-2000 (ref. 35), but similar to those reported for Diego Garcia (~1.75 mm yr -1 ) over the period 1950-2009 (ref. 36). However, a defining feature of Indian Ocean sea-level trends over time, and certainly over recent decades, is one of significant intra-and inter-annual variability, in part driven by the influence of the Indian Ocean Dipole (IOD). This causes marked variations in sea-surface temperatures and wind regimes across the region, and drives considerable cross-basin changes in sea-level. Caution is thus needed in interpreting long-term trends from the relatively short term datasets presently available. Whilst Dunne 34 argue that because of such inter-annual variability the long-term regional trend is not significantly different from zero (2.22 ± 1.80 SE mm yr -1 (ref. 34), it is perhaps more appropriate to conservatively assume that the local rate of rise over the last half century or so has probably not exceeded the global average of ~1.7 mm yr -1 (ref. 37) and ~2.2 mm yr -1 as measured over the last ~20 years in Chagos 34 . To compare contemporary reef accretion potential to future sea level rise trajectories we used the IPCC AR5 report projections for the period 2081-2100 (ref. 37), but also accounting for the impacts of future wind-stress 38 .