Depth variation in benthic community response to repeated marine heatwaves on remote Central Indian Ocean reefs

Coral reefs are increasingly impacted by climate-induced warming events. However, there is limited empirical evidence on the variation in the response of shallow coral reef communities to thermal stress across depth. Here we assess depth-dependent changes in coral reef benthic communities following successive marine heatwaves from 2015 to 2017 across a 5–25 m depth gradient in the remote Chagos Archipelago, Central Indian Ocean. Our analyses show an overall decline in hard and soft coral cover and an increase in crustose coralline algae, sponge, and reef pavement following successive marine heatwaves on the remote reef system. Our findings indicate that the changes in benthic communities in response to elevated seawater temperatures varied across depth. We found greater changes in benthic group cover at shallow depths (5–15 m) compared to deeper zones (15–25 m). The loss of hard coral cover was better predicted by initial thermal stress, whilst the loss of soft coral was associated with repeated thermal stress following successive warming events. Our study shows that benthic communities extending to 25 m were impacted by successive marine heatwaves, supporting concerns about the resilience of shallow coral reef communities to increasingly severe climate-driven warming events.

Methods
Study Sites
The Chagos Archipelago is a remote archipelago located in the Central Indian Ocean. Situated at the southern end of the Laccadives-Maldives-Chagos ridge, it covers a total area of 640,000 km2 and consists of five islanded atolls. The archipelago has been uninhabited since the early 1970s, except for the US military base on the southern atoll, Diego Garcia (81). Benthic community data were collected with repeat sampling, in March–April 2013/2014 (before) and 2018/2019 (after) following successive severe marine heatwaves in 2015–2017. The peaks of marine heatwave events occurred between May–June 2015 (4.2 ± 0.4 – 6.1 ± 0.3 °C-weeks) and April–June 2016 (12.1 ± 0.1 – 15.8 ± 0.1 °C-weeks) (Figure S1). Surveys in each period were carried out at 16 sites on forereef slopes across four atolls: Peros Banhos (PB), Salomon (SA), Great Chagos Bank (GCB), and Egmont (EG) (Figure 1), and conducted at four depth zones: 5–10 m, 10–15 m, 15–20 m, and 20–25 m. For ease of interpretation, shallower depth zones in the text refer to benthic groups occurring between 5–15 m and deeper zones to benthic groups occurring between 15–25 m.

Data collection and benthic community assessment
Benthic community assessment was carried out using 10-minute continuous video swims conducted at a constant speed (0.1 ms-1) at each site and depth zone. Equipped with two spotlights and two red laser pointers set at 10 cm apart, the camera set-up provided a scale of measurement of the benthos and adjustment for lower light levels at greater depth. During the survey, the camera was held approximately 0.5 m above the substrate and at a 45° angle to capture benthic organisms and substrate types under overhangs and canopies (82). Each video was converted into a sequence of still images (25 frames per second in Pinnacle Studio, v22.2.0). Thirty of these images were selected for analysis using Matlab (R2018a.Ink) (total across all sites, n = 5610). Selected images were separated by 80–100 frames in each sequence to avoid re-sampling the same area of the reef.

Benthic composition within each image was quantified using Coral Point Count with Excel extensions (CPCe) (83). To reduce observer bias, image analysis was equally distributed among SSP, RR, and LR (10 images per person per site). The percentage cover of benthic organisms and other substrate types was quantified at 15 randomly allocated points on each image, using a stratified random design (84). The broad classification of benthic groups was adapted from Denis et al. (2017) and the NOAA Coral Reef Information System (86). This study examined five benthic groups that have been found to undergo shifts in dominance and contribute to habitat building on coral reefs following thermal stress events: hard coral, soft coral, sponge, crustose coralline algae (CCA), and reef pavement (48,87–90). Despite being an important competitor to corals (47), macroalgae were not included as we only found small changes in macroalgal cover (turf and fleshy macroalgae) following marine heatwaves (pre-heatwaves: 6.5 ± 0.2 % and post-heatwaves: 6.0 ± 0.1%).

Estimating change in benthic groups following 2015–2017 bleaching events
The change in benthic group cover which occurred following the 2015–2017 warming events was calculated using:

δ = (log(Vf) - log(V0))/t

where δ is the geometric logarithmic change in cover, Vf  is the percentage cover of the benthic group collected at the end of the time series, i.e. 2018 or 2019 and V0  is the percentage cover of benthic groups at the beginning of the time series, i.e., 2013 or 2014 and t is the duration of the time in years between benthic surveys, i.e., before and after the successive marine heatwaves (91). This approach optimised the dataset by including surveyed sites that were monitored in different years before and after the heat stress event, i.e., to include pre-disturbance community data collected in 2013 or 2014 (hereinafter pre-heatwaves), and post-disturbance community data sampled in 2018 or 2019 (hereinafter post-heatwaves). The geometric logarithmic change in cover metric allows the quantification and comparison of non-linear time series by taking into consideration the exponential decline and increase in cover of benthic communities over time (sensu 91) (for example 63,92). Compared to other metrics that measure ecological change over time (92), the geometric logarithmic rate of change preserves proportionality of change in cover (e.g. a 50% to 5% loss is proportionately similar to a 10% to 1% loss) and is symmetrical, i.e., produces similar rates of change for declines and matching increases (e.g. a 50% to 5% loss and a 5% to 50% increase yields the same rates (91).To summarise the geometric logarithmic change in cover metric, a negative value represents a loss in benthic group cover from pre- to post-heatwaves and a positive value indicates an increase in benthic group cover following thermal stress. Variation in the change in benthic group cover across depth zones was visualised using a forest plot (ggplot: ggplot2 package, (93)). For ease of interpretation, the log change in cover was back-transformed to the percentage change in cover (Figure S2), using:

i = (e^δ  - 1) x 100

where δ is the geometric logarithmic change in cover and ἱ is the percentage change in cover.

Exposure to repeated thermal stress
To determine how benthic groups were affected by thermal stress across depth zones during the 2015–2017 marine heatwaves, maximum degree heat week (hereinafter maxDHW) was used to quantify thermal anomalies from 2015–2017. DHW is a proxy for accumulated thermal stress, represented by a 1°C increase above the local mean climatic temperature over a 12-week period at a given pixel and expressed as degree Celsius weeks (°C-weeks); and maxDHW is the annual maximum accumulated thermal stress in a year. Using the National Oceanic and Atmospheric Administration (NOAA) Coral Reef Watch 5 km resolution product (94), the annual maxDHW between January 2015 and December 2017 was extracted for each study site (using R packages: ncdf4 (95), raster (96), rgdal (97) and sp (98)). A DHW threshold of 4°C-weeks may be indicative of significant coral bleaching, and a DHW value of 8°C-weeks is a signal for severe and widespread bleaching with likely mortality (11,99–102). Local-scale hydrodynamics for e.g. upwelling and changes in mixed layer depth can influence temperature regime across depth gradients (103,104). However, in the absence of depth-specific temperature data from in-situ loggers across all study sites,  DHW that were derived from SST were used as a proxy of thermal stress experienced along the depth gradient at the study sites during the warming events (Figure S1). Here, initial thermal stress refers to maxDHW recorded in 2015 in the first year of the successive severe heatwaves (4.1–6.4 °C-weeks) and repeated thermal stress refers to the cumulative maxDHW recorded in all years from 2015–2017 (17.1–21.1 °C-weeks, Table S1). Cumulative maxDHW has been shown to be a good predictor of bleaching in the Western Indian Ocean (105) and therefore the change in benthic group cover following repeated marine heatwaves.