Thermal stress-related spatiotemporal variations in high-latitude coral reef benthic communities

High-latitude coral reef communities have been postulated as the first areas to undergo reorganisation under climate change. Tropicalisation has been identified in some high-latitude communities and is predicted in others, but it is unclear how the resident benthic taxa are affected. We conducted a long-term (2007–2016) assessment of changes to benthic community cover in relation to thermal stress duration on the Southeast Florida Reef Tract (SEFRT). Thermal stress events, both hot and cold, had acute (thermal stress duration affected benthic cover that year) and chronic (thermal stress duration affected benthic cover the following year) impacts on benthic cover. Chronic heat stress was associated with declines in cover of the reef-building coral families Acroporidae, Montastraeidae, Meandrinidae, Mussidae and Siderastreidae, which coupled with the absence of cold stress and rising annual temperatures boosted macroalgae cover. Cover of smaller, weedy coral families, Poritidae, Agariciidae and Astrocoeniidae, was either unaffected or positively related to heat stress duration and rising mean temperature. Thermal stress was related to spatiotemporal variations in benthic cover on the SEFRT, likely enhanced by local stressors, such as elevated nutrients and sedimentation. Coral and octocoral cover declined within four of six sub-regions, sponge cover increased in half of the sub-regions, and macroalgae cover increased in four sub-regions during the study. Under current conditions, increased macroalgae and weedy coral cover are anticipated to inhibit reef recovery.


Introduction
High-latitude coral reef communities are of substantial importance as climate change strengthens. Communities towards the edge of a species range are generally limited by environmental conditions near their physiological threshold, but are also subjected to more extreme climatic events (Busch et al. 2011;Rehm et al. 2015). With range shifts predicted under ocean warming, high-latitude communities may play a crucial role in species persistence, enabling species to maintain or expand their distribution (Gibson et al. 2009;Vergés et al. 2014;Rehm et al. 2015). Understanding recent spatiotemporal changes to high-latitude benthic communities and their relationship with temperature is thus critical.
Globally, the benthic composition of tropical reef ecosystems has changed significantly in just a few decades (Gardner et al. 2003;Bellwood et al. 2004;Hoegh-Guldberg and Bruno 2010;Hughes et al. 2018), with a shift from scleractinian corals (hereafter referred to as corals) to turf algae, macroalgae, octocorals or sponges in some locations (Hughes 1994;Norström et al. 2009;Jackson et al. 2014;Pawlik et al. 2016). The decline in coral cover has been linked to changes in environmental and physical parameters, such as ocean temperature, acidification, storm prevalence and local anthropogenic stressors (Walther et al. Topic Editor: Morgan S.Pratchett Electronic supplementary material The online version of this article (https://doi.org/10.1007/s00338-020-01994-8) contains supplementary material, which is available to authorized users. 2002;Hughes et al. 2003;Banks et al. 2008;Hoegh-Guldberg et al. 2017;Hughes et al. 2018). Concurrently, the reorganisation of coral reef communities has occurred following intense thermal stress events (Walther et al. 2002;Hughes et al. 2003;Hoegh-Guldberg et al. 2017), defined as periods where temperatures rise above or fall below threshold levels found to be detrimental to marine organisms (Bruno et al. 2007). As climate change strengthens, high-latitude reef systems are becoming increasingly valuable study areas to predict the future of coral reefs worldwide. The immigration of tropical taxa with ocean warming has been predicted to lead to the tropicalisation of sub-tropical coral reefs, resulting in macroalgae cover reduction and range expansions for coral species (Yamano et al. 2011;Baird et al. 2012;Vergés et al. 2014). High-latitude reefs have also been proposed as refugia (Riegl 2003;Cacciapaglia and Woesik 2016), acting as resilient areas from temperature stress (Pandolfi et al. 2011) and as tropical/temperate transition zones that may function as areas of heightened evolutionary activity (Kawecki 2008;Gibson et al. 2009;Beger et al. 2014). With these possibilities in mind, understanding changes to the resident taxa on high-latitude coral reefs and their drivers is fundamentally important.
One such high-latitude reef system, which can be used to investigate spatiotemporal changes to the benthic community and the impact of temperature on benthic assemblages, is the Southeast Florida Reef Tract (SEFRT). The SEFRT, the northern portion of the Florida Reef Tract and within the Southeast Florida Coral Reef Ecosystem Conservation Area, spans 170 km on a latitudinal gradient (Walker 2012) and is the northern limit of coral distribution in the Western Atlantic. Distinct sub-regions have been defined, based on known biogeographic boundaries (Walker 2012) and major port channels (Fig. 1). Latitude, biogeographic boundaries and major divides in the reef have all been identified as delineators for changes in the distribution and diversity of marine organisms (Stevens 1989;Veron 1995;Walker 2012). Here, we used long-term annual monitoring data (10 years, from 2007 to 2016) to assess spatial and temporal variations in benthic community cover on the SEFRT in relation to temperature and thermal stress. Benthic cover was calculated from annual photographic data spanning 31 permanent sites over six sub-regions running north to south along the SEFRT. Temperature data, collected via bihourly in situ measurement, were analysed concurrently with benthic cover data to assess: (1) how benthic community cover varied spatially and temporally on the SEFRT from 2007 to 2016 and (2) how temperature fluctuations and thermal stress duration affect variations in benthic community cover and within the coral community.

Study location
The Southeast Florida Reef Tract (SEFRT) extends in a northerly direction from Biscayne Bay to St Lucie Inlet along the east coast of the south Florida peninsula (Fig. 1). It comprises three parallel linear reefs (west to east from the coastline; Fig. 1 inset centre), separated by sand channels (Banks et al. 2007), with predominately low coral cover benthic communities, which lie within 3 km of a heavily urbanised coastline. Geologically, the three linear reefs are relict reef structures on the southeast Florida shelf. They initiated as fringing reefs at different points in time over the last 10,000 years, before transitioning to barrier reefs with sea level rise (Banks et al. 2007(Banks et al. , 2008. Additionally, a nearshore hard bottom ridge complex, consisting of karst deposits, is found closest to shore along the SEFRT. North to south, the SEFRT was split into six sub-regions ( Fig. 1) based on previously identified biogeographic boundaries (Walker 2012), and by the major port channels, Port Everglades and Port of Miami. The Bahamas Fault Zone, which marks the northern terminus of the outer reef, separates Martin (sub-region 1) and South Palm Beach (sub-region 2). The middle reef ends just north of Boca Inlet, which separates South Palm Beach and Deerfield (sub-region 3). Deerfield stretches North from Boca Inlet to Hillsboro Inlet, where the inner reef ends. Broward (sub-region 4) stretches from Hillsboro Inlet to Port Everglades, and Miami (sub-region 5) south to Government Cut at the Port of Miami. The final sub-region is Biscayne (sub-region 6), which extends south from Government Cut to central Biscayne Bay. For this study, sites on the nearshore hard bottom ridge complex and inner reef, both at depth ranges between 4 and 10 metres, were grouped. Inner reef sites are in all sub-regions, middle reef sites in Deerfield, Broward and Miami and outer reef sites in South Palm Beach, Deerfield, Broward, Miami and Biscayne. Benthic survey sites were positioned to represent each reef habitat (inner, middle and outer reef) in each subregion. Therefore, sampling effort varied between sub-regions to reflect differences in potential reef habitat. Sampling effort also varied by year, with no surveys in Deerfield in 2009 or Biscayne before 2010 (Table S1).

Monitoring projects
Data from two long-term annual monitoring projects, collected from 2007 to 2016 by the Coral Reef Restoration, Assessment and Monitoring (CRRAM) lab at Nova Southeastern University (NSU), were used in this study. Photographic data collected as part of the Southeast Florida Coral Reef Evaluation and Monitoring Project (SECREMP) were used to calculate benthic cover at four permanent 22 9 0.4 m belt transects at 22 sites, ranging from 5 m to 17 m depth. Photographic data from nine sites, each consisting of a single permanent 20 9 1.5 m transect surveyed from 2007 to 2016 as part of the Broward County Biological Monitoring Project (BC Bio), were used to supplement SECREMP and increase spatial and temporal coverage. BC Bio site depths ranged from 4 to 20 m.

Image analysis
SECREMP images were analysed using PointCount'99 (Dustan et al. 1999), with the benthic taxa directly underneath 15 randomly placed points per frame identified, for a total of 900-1000 points per transect (60-67 images per transect), to determine percentage cover. Photographic quadrats for BC Bio sites were analysed using Coral Point Count with excel extensions (CPCe 4.1) research software (Kohler and Gill 2006). The benthic taxa underneath 25 random points per image, totalling 1000 per site, were identified to determine the per cent cover (40, 1 9 0.75 m images per transect). Benthic taxa identified included corals (to species level), octocoral, sponges, macroalgae, crustose coralline algae, turf algae/substrate and cyanobacteria. A power analysis was conducted to estimate the minimum necessary number of points to identify per image to ensure comparability between BC Bio and SECREMP transects and justify using sites from both longterm monitoring projects. Three sites were monitored over   Ruzicka et al. 2013). A second power analysis was conducted to ensure that the minimum detectable level of change was suitable for the cover changes considered in this study. The power analysis for minimum detectable level of change was performed on both data sets. Both power analyses were calculated using the power.t.test function in R (R Core Team 2016), based upon a two-tailed t test from Cohen (1988) (Eq. 1).
Equation 1 Power Analysis calculating minimum detectable difference in cover (Cohen 1988). n = number of samples; d = standardised mean difference, calculated from l1 to l2 with l being the mean; d = minimum detection level; t = 2-tailed t test statistic; a = level of significance; v = degrees of freedom; b = power of the test.

Temperature analysis
Water temperature data were compiled using two methods to cover the entire SEFRT and enable analysis of Sea Surface Temperature (SST) and temperatures at reef depth, the latter of which are reflective of in situ conditions experienced by the benthic community. In situ temperatures were recorded at all SECREMP sites between February 2007 and June 2016 using HOBO temperature loggers (HOBO Pro V2). Hindcast, SST data from the Hybrid Coordinate Ocean Model (HYCOM; https://www. hycom.org) were used to analyse trends in maximum, mean and minimum SST over the entire SEFRT. Thermal stress was assessed in relation to duration above or below threshold levels only. Thermal thresholds (thermal high or thermal low) for each sub-region were calculated using HYCOM data. Only HYCOM point data near sites and within the bounds of the SEFRT were used. HYCOM SST was compared to in situ temperature data for errors, (e.g. slightly elevated annual SST indicative of the Gulf Stream). Thermal stress duration was calculated from in situ temperature loggers, located at site depth. Thermal highs were considered as 1°C above the maximum of the mean summertime SST, and thermal lows 1°C below the minimum of the mean wintertime SST over the study period (2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016), defined using NOAA's coral reef watch limits (Liu et al. 2013), which reflect thermal stress found to significantly affect coral health over a prolonged period. Thermal stress thresholds were calculated independently for each sub-region to account for environmental acclimation of the benthic community with latitude. The maximum mean summertime temperature was defined as the maximum of the mean SST during summer months (July-September) during the study. Heat stress duration was defined as the number of days where average daily in situ temperatures at a specific site are 1°C above the maximum of the mean summer temperatures. The minimum mean wintertime temperature was defined as the minimum of the mean temperature during winter months (December-March) during the study. Cold stress duration was defined as the number of days where in situ temperatures at a specific site are 1°C below the minimum of the mean winter temperatures.

Statistical analysis
Statistical data analysis was conducted in (R Core Team 2016). Initial evaluation of the spatial and temporal similarities in benthic cover between sub-regions was conducted using principal component analysis (PCA; Wold et al. (1987)). Benthic community data were consolidated into the major benthic groups: coral, octocoral, sponge, macroalgae, zoanthid and ''other'' fauna (anemones, corallimorphs, hydrozoa, tunicates, polychaetes, crustose coralline algae and cyanobacteria). Coral community data were assessed at the family level, due to low cover of certain species and because species identification can be difficult using imagery (Table S2). Kendall-Theil nonparametric regression was used to assess temporal change in SST maximum, mean, minimum and range during the study period to compare to changes in benthic cover. Kruskal-Wallis test (KW) was used to assess spatiotemporal change in thermal stress duration (n = 220).
Generalised linear mixed-effects models (GLMMs) were used to analyse potential spatiotemporal drivers of change in cover of the major benthic taxa and within the coral community. A single model for each benthic taxa was fitted to assess whether cover significantly varied over time and space and if temperature, particularly thermal stress duration, accounted for some of this variation using the glmmTMB function in the ''glmmTMB'' package (Brooks et al. 2017) and included all observations from all years (n = 769; Eq. 2). Multiple transects (observations) were surveyed at a site, and there were multiple sites within a sub-region. We accounted for the hierarchical structure of the data by applying a mixed-effects model with the random effect site nested within the random effect sub-region. Binomial GLMMs were primarily used because the response variable (in this case the benthic cover) was a percentage, to factor the non-Gaussian distribution of data, differences in sample size between years and to account for the hierarchical data structure. Spatial, temporal and temperature-related predictors were incorporated into each model (Table 1). Sample year, sub-region and reef habitat were included in the model as fixed categorical predictors. Time was included as a categorical variable so that results were not skewed by the first or last year of the study and to allow for between year comparisons. However, we verified that equivalent results would have been obtained if time had been considered as a continuous predictor. We considered including the interaction between sub-region and year, but were unable to due to lack of convergence, likely caused by sample size changes over time, specifically no surveys in Deerfield in 2009 and no surveys in Biscayne before 2010. Supplementary analysis of temporal trends within each sub-region was therefore conducted using individual sub-region GLMMs for the major benthic taxa.
Eight temperature-related variables (Table 1) were calculated for their effect on the major benthic taxa. These were derived from two main metrics-in situ temperature and thermal stress-and were selected to reflect current and potential future impacts of climate change and thermal stress duration. A lag effect of thermal stress on the benthic community the following year (for example, the effect of 50 days of heat stress in 2010, on the benthic community in 2011) was included in each model to assess for the latent impact of temperature stress on benthic cover. The full model (Eq. 2) was fitted to the data. Prior to model selection, correlation of variables was assessed using pairs plots and tested using Spearman's rank correlation coefficient in R using the ''cor'' function. Minimum temperature correlated strongly with mean and temperature range, heat stress duration with maximum temperature and cold stress duration with mean, range and latent cold stress duration. When testing for multicollinearity using Variance Inflation Factor (VIF) based upon the check_collinearity function in the ''performance'' package (Lüdecke et al. 2019), minimum temperature was confirmed to be highly collinear in all models, as was sample year. Variables were removed from the model one at a time, starting with the highest VIF and then retested. Sample year and minimum temperature were omitted from all models, and temperature range omitted from several models. Then a backward stepwise approach was taken to select the minimum adequate model using Akaike Information Criterion (AIC). The model with the lowest AIC was selected as the minimum adequate model. In the event of two models having an AIC score within two, the simplest model was chosen as the best fit for the data.
Equation 2 GLMM examining variation in benthic taxon, where Benthic cover ijk refers to the jth observation at site i in sub-region k. Random intercept site is nested within sub-region and is assumed to be normally distributed with mean 0 and variance r 2 .
Model validation was assessed by plotting deviance residuals against fitted values and deviance residuals against each significant variable in the minimum adequate model. Temporal autocorrelation was assessed by plotting residuals against time for each survey site (Figs. S1-13). This approach was repeated for each of the four main benthic taxa (coral, macroalgae, octocoral and sponge) and for each coral family (Table S2), resulting in a specific minimum adequate model for each benthic taxon. Post hoc analysis was performed on fixed categorical predictors in

Temperature change
Water temperature on the SEFRT fluctuated seasonally with an annual mean temperature of 26.38°C (± SE 0.02°C) and annual temperature range of 8.46°C (± SE 0.11°C; Fig. 2a). Annual mean SST and minimum SST significantly increased over the study period (Kendall-Theil nonparametric regression, NPR, p \ 0.001). Heat stress duration and cold stress duration (Fig. 2b)

Temporal trends
The taxonomic composition of benthic cover on the SEFRT changed over time, with a 63% decline in coral cover and a concomitant 92% increase in macroalgae cover from 2007 to 2016. Within this ten-year period, the greatest changes occurred from 2015 to 2016: coral cover declined 43% and macroalgae cover increased by 94%. Additional changes in the other major benthic taxa were seen from  (Figs. 3 and 4) suggest there would be a significant interaction between sub-region and year, but we were unable to appropriately test that interaction in our model (see methods). We conducted additional analysis, modelling temporal trends within each sub-region.

Sub-region benthic cover change
Significant declines in coral, octocoral and macroalgae were seen in sub-regions with the highest cover of each taxa over the study (Fig. 3)

Drivers of benthic community variation
Initial assessment of the spatiotemporal variations between sub-regions using PCA revealed that macroalgae had the largest influence on spatiotemporal similarity between subregions and years in principal component 1, with octocoral,  Table S29). Modelling of these taxa was therefore conducted to assess the sources of cover variation in relation to time, space and temperature.
For each taxon, sample year was removed from the minimum adequate model. Even though year had a significant effect on the cover of all benthic taxa (Tables S5-S28), VIF testing revealed high correlation with the significant temperature variables caused multicollinearity. The relationship between sample year and each benthic taxon was subsequently tested to assess temporal trends. Coral cover was negatively related to increasing maximum temperature, mean temperature, latent heat stress duration and latent cold stress duration in the minimum adequate model (Table S4). Coral cover differed among reef habitat, with coral cover on the inner reef significantly higher than the middle and outer reefs (p = 0.01). Macroalgae cover was positively related to maximum temperature, mean temperature and latent heat stress duration in the minimum adequate model. It was negatively related to heat stress duration, cold stress duration and latent cold stress duration (Table S6). Macroalgae cover differed by reef habitat and sub-region, with cover higher on the outer reef than the inner reef (p = 0.03) and in Martin than in all other subregions (p \ 0.01). Octocoral cover was positively related to mean temperature, cold stress duration and latent cold stress duration and negatively related to latent heat stress duration in the minimum adequate model (Table S8).
Octocoral cover varied by sub-region and reef habitat, with cover on the outer reef higher than on the inner reef (p = 0.001) and cover in Martin lower than in all other regions (p \ 0.001). Sponge cover was positively related to maximum temperature, latent heat stress duration and latent cold stress duration in the minimum adequate model (Table S10). Sponge cover differed by reef habitat, with sponge cover significantly lower on the inner reef than on the middle reef and outer reef (p \ 0.001). It should be noted that multiple temperature variables within the sample year were omitted from models due to collinearity with other variables. These variables may be as good a predictor for benthic cover as those variables included.

Drivers of coral community variation
Coral families on the SEFRT can be broadly divided into those considered reef-building, Montastraeidae, Merulinidae, Acroporidae, Mussidae, Meandrinidae and Siderastreidae and those thought to be weedy, Poritidae, Agariciidae and Astrocoeniidae. Multiple coral families considered reef-building were negatively related to latent heat stress duration (Tables S12-S22). In the minimum adequate model, Montastraeidae cover was positively related to heat stress duration and negatively related to maximum temperature and latent heat stress duration. Acroporidae cover was negatively related to mean temperature, latent heat stress duration and latent cold stress between sub-regions. PC1 and PC2 account for 85% and 9% of variation, respectively. Colours and shapes are consistent for sub-regions (Fig. 1); years are labelled. Arrows reflect influence of taxa driving dissimilarity in community composition between samples duration. Merulinidae cover was positively related to maximum and mean temperature. Mussidae cover was positively related to temperature range and negatively related to heat stress duration, latent heat stress duration and latent cold stress duration. Mussidae cover differed by reef habitat, with cover significantly higher on the inner reef than middle and outer reefs (p \ 0.01). Meandrinidae cover was positively related to heat stress duration and negatively related to latent heat stress duration. Siderastreidae cover had a positive relationship with mean temperature and a negative relationship with latent heat stress duration temperature.
Coral families considered weedy were either unrelated or positively related to latent heat stress duration and rising mean temperature (Table S24-S28). Poritidae cover was positively related to mean temperature. Agariciidae cover was negatively related to latent cold stress duration and varied by reef habitat. Astrocoeniidae cover was positively related to latent heat stress duration. Reef habitat was found to help explain Astrocoeniidae cover variation, with cover on the outer reef higher than the inner reef (p \ 0.05).

Discussion
Ocean warming and thermal stress events have contributed to multiple spatiotemporal changes in benthic community composition on the SEFRT. Latent heat stress duration, rising maximum annual temperature and a significant increase in mean temperature were related to a decline in coral cover. From 2007 to 2016, significant changes in coral community composition were found, with declining trends in the cover of reef-building families and increasing trends in smaller, weedy families associated with latent heat stress duration. Negative feedback loops associated with increased macroalgae cover and coral community change to weedy species, coupled with continued increases in water temperature and thermal stress events, are likely to hinder the recovery potential of reef-building corals on the SEFRT.
Coral cover significantly declined on the SEFRT over the past decade, driven by declining coral cover in the subregion with the highest coral cover, Broward. Our data suggest that, as has been proposed elsewhere (Loya et al. 2001;Côté and Darling 2010;van Woesik et al. 2011), coral community change results in ''winners'' and ''losers'', based upon the resilience of coral species to, among other factors, thermal stress. The cover of reefbuilding coral families, Montastraeidae, Acroporidae, Meandrinidae, Siderastreidae and Mussidae, was negatively related to latent heat stress duration, which peaked following the consecutive El Niño years in 2014 and 2015. Conversely, the cover of weedy coral families, Agariciidae, Astrocoeniidae and Poritidae was either positively related or unrelated to latent heat stress duration. Coral cover decline on the SEFRT largely resulted from disease, which coincided with substantial thermal stress. The primary reefbuilding coral families on the SEFRT, Acroporidae (Acropora cervicornis) and Montastraeidae (Montastraea cavernosa), suffered significant declines, particularly in areas where cover was previously highest. From late 2014, an unprecedented disease, termed stony coral tissue loss disease (SCTLD), spreads along the Florida Reef Tract (FRT), resulting in mortality of multiple species on the SEFRT, including key reef builders M. cavernosa, Meandrina meandrites and Siderastraea siderea (Walton et al. 2018). Elevated disease prevalence during warm conditions has also been documented separately for A. cervicornis (Goergen et al. 2019). The increased prevalence and susceptibility of corals to disease have previously been found to result from an increase in the duration of temperatures above the bleaching threshold (Bruno et al. 2007). Cooling in winter months may previously have negated a portion of the impact of heat stress on corals on the SEFRT and suppressed disease, but annual minimum and mean temperatures have significantly increased. Our findings suggest that latent thermal heat, coupled with rising annual temperatures, may have exacerbated the impacts of additional stressors on the SEFRT, resulting in declining cover of many reef-building species. The comparatively resilient poritids, astrocoeniids and agariciids, which increased in cover at specific sites, appear as ''winners'' on the SEFRT, but do not provide the equivalent level of cover or structural complexity that was lost. Their apparent heightened resistance to disease and thermal stress (Côté and Darling 2010) will likely lead to further increases in cover and larval output, in turn increasing benthic competition and contributing to the negative feedback loop predicted to impact the recovery of reef-building coral populations.
Heat and cold stress duration had both acute (thermal stress duration affected benthic cover that year) and chronic (thermal stress duration affected benthic cover the following year) effects on benthic cover on the SEFRT, with increased macroalgae in 2016 following coral cover decline from heat stress being most conspicuous. While macroalgae blooms can be ephemeral, macroalgae was positively related to increasing maximum temperature and latent heat stress duration. An increasing trend in macroalgae cover was found in half of the sub-regions from 2013 to 2016, culminating in an increase of over 60% in four sub-regions from 2015 to 2016. This appears to have resulted from a combination of the additional space available after the decline in coral cover (Mumby and Steneck 2008), prolonged warm water temperatures and an absence of any major disturbance (hurricane or extreme cold stress) during 2015 and early to mid-2016, which remove macroalgae (Mumby et al. 2005). The macroalgae cover increase found on the SEFRT does not support evidence of tropicalisation on high-latitude reefs found elsewhere, where macroalgae have been consumed by herbivores following poleward species migration (Vergés et al. 2014). Our data suggest that increased cold stress duration, such as those felt in 2010 and 2011, may have caused macroalgae cover to decline, facilitating subsequent increases in sponge cover. Sponges likely benefited from macroalgae suppression and the newly available substrate, as sponges which lack chemical defence have previously been shown to grow rapidly after disturbance (Loh and Pawlik 2014;Wooster et al. 2017). The positive association of octocoral cover with cold stress duration is also likely related to macroalgae suppression, with a reduction in smothering by macroalgae resulting in increased recruitment and juvenile survival, as has been found for scleractinian corals (Box and Mumby 2007). While previous studies in the Florida Keys found intense cold stress caused sponge and octocoral decline when water temperature fell below 16°C , water temperature never dropped near that level in any of the sub-regions where octocoral or sponge cover was high. Octocoral cover was instead negatively associated with increased latent heat stress duration supporting results from studies in the Pacific (Sammarco and Strychar 2013), which have shown zooxanthellae symbionts in octocorals to be as susceptible to heat stress as those in corals. Equally, studies in the Mediterranean and Caribbean have shown that heat stress causes thermal bleaching and mass mortality in gorgonians (Garrabou et al. 2009;Prada et al. 2010). The positive association of sponge and octocoral cover to cold stress duration, which declined over the study period, coupled with heat stress events is predicted to continue favouring macroalgae growth over coral, octocoral and sponges. Macroalgae can reduce fecundity, inhibit recruitment, impair growth and cause tissue mortality in corals (Lirman 2001;Box and Mumby 2007). We therefore anticipate a negative feedback between macroalgae and coral which hinders recovery under present conditions. While thermal stress and temperature variability explained a portion of the observed variance in benthic cover found on the SEFRT, multiple stressors, including coastal development, nitrogen enrichment, dredging and hurricane damage, have also been suggested as causes for coral decline on the FRT (Ruzicka et al. 2013;Walton et al. 2018;Cunning et al. 2019;Lapointe et al. 2019). Analysis of the specific impacts of local stressors or biological/ecological factors was beyond the scope of this study, but initial differences in benthic cover within sub-regions and reef habitats were present and a significant relationship between each benthic taxon and sub-region and/or reef habitat found, suggesting these factors account for some spatiotemporal variation. Coastal development and marine construction activities, which reduce irradiance through turbidity and sedimentation, are widespread on the SEFRT. Multiple ports, inlets and sewage outfalls are found along the coastline, contributing to elevated nutrients and contaminants (Futch et al. 2011). Chronic local stress from poor water quality and increased sedimentation is often cited as contributing to coral disease and affecting recovery potential (Ortiz et al. 2018;Lapointe et al. 2019;Otaño-Cruz et al. 2019), but refugia potential has been predicted in other turbid, nearshore areas (Cacciapaglia and Woesik 2016;van Woesik and McCaffrey 2017). Additionally, we predict that some correlations between temperature variables and benthic taxa, specifically the positive association of Montastraeidae, Mussidae and Meandrinidae with heat stress duration and the negative relationship between heat stress duration and macroalgae cover, are likely influenced by other biological or environmental factors. Heat stress duration was predominately higher at shallow sites, but shallower reefs often have higher irradiance, which can increase coral calcification rates (Tomascik and Sander 1987), and an increased number of herbivorous fishes and wave action, which could suppress macroalgae and enhance coral survival and growth (Burkepile and Hay 2009). Understanding the interactive effects of biological control and local and global stress will be vital to assess the prospects of resilience for coral reef communities worldwide.
Local stressors have also likely affected the wider benthic community. Macroalgae and octocoral cover significantly varied by sub-region, driven by high macroalgae cover and low octocoral cover in the northernmost subregion, Martin, which lies north of the Bahamas fault zone. Conditions there more closely resemble those associated with a temperate ecosystem, particularly with lower winter water temperature; future studies should consider separating this sub-region from those to the south. Octocoral and sponge cover was higher than coral cover in most of the sub-regions. Sponge cover was positively related to increasing maximum temperature, as well as heat stress and cold stress duration, while octocoral cover was positively related to mean temperature and cold stress duration. While this suggests they may have increased resilience to thermal stress, sponge cover also varied by reef and octocoral cover by sub-region and reef. Resilience to a range of stressors has been posited as a cause of octocoral increase in the Florida Keys (Ruzicka et al. 2013), while current flow has been suggested for an increased abundance of cold water octocoral in temperate regions (Yesson et al. 2012). Sponge cover was highest between two inlet channels and two sewage outfalls. The location may contribute to an increased influx of nutrients and dissolved organic carbon (DOC) which have been shown to support sponge growth (Ward-Paige et al. 2005;Pawlik et al. 2016). The impact these local stressors have on the benthic community should continue to be investigated.
At the start of this study, differences in the major benthic taxa cover were identified between sub-regions on the SEFRT. While there is still variability between sub-regions and they can still be broadly characterised by these taxa, biotic homogenisation is occurring throughout the SEFRT, with declines in coral and octocoral cover and replacement by macroalgae or sponges. Coral cover significantly declined in two-thirds of the sub-regions in 2016, following a similar trend to that reported in parts of the Caribbean and Florida Keys (Gardner et al. 2003;Ruzicka et al. 2013;Jackson et al. 2014;Estrada-Saldívar et al. 2019). Octocoral cover declined in four sub-regions, in contrast to elsewhere in the Caribbean where increasing (Ruzicka et al. 2013) or stable trends have been reported (Lenz et al. 2015). Declining coral and octocoral cover in sub-regions where it was highest at the start of the study (Broward and South Palm Beach, respectively), loss of reef-building species and increased macroalgae cover, which contributed most to benthic cover in four sub-regions, suggests that resilience is severely compromised. The potential loss of uniqueness of ecosystem regions is a major concern. The combination of environmental conditions driving benthic cover change and the negative feedback loops associated with community changes suggest recovery may be difficult.
Coral reef community changes found on the SEFRT highlight the multiplicative, often interacting, impacts that ocean warming and thermal stress events have on resistance in multiple taxa and suggest that other high-latitude reef communities will also be negatively affected by climate change. Our long-term monitoring data identified specific drivers of change within benthic communities, which can be used to predict changes globally. In the past, proposed trajectories and drivers of community changes in coral reefs have often been theoretical exercises based on expertise or short-term data (Bellwood et al. 2004). Here we identify the specific thermal stress drivers of coral reef community change based on long-term data collected over a large spatial scale. The schematic proposed (Fig. 5) reveals how different environmental stress scenarios, which may arise from climate change, can result in multiple different community change trajectories depending on the type, frequency and severity of thermal stress they experience. While normal disturbance regimes with heat and cold stress did not always result in reef-building coral decline, latent heat stress from repeated thermal events and increasing temperatures were related to changes in communities away from reef-building coral. Macroalgae and the weedy coral families, Agariciidae, Poritidae and Astrocoeniidae, were positively related to or unaffected by heat stress duration and increased or remained stable in cover following thermal stress on the SEFRT. However, other environmental factors may create additional trajectories for benthic community cover, which may impact reef resilience and the potential for coral reefs to recover in the Anthropocene. Fig. 5 Schematic depicting the drivers of benthic community cover on coral reefs. Curved arrows represent normal disturbance regimes following thermal stress, with reefbuilding and weedy coral cover remaining stable. Straight arrows represent likely directional change in cover following stated thermal stress scenario; red = heat stress, blue = cold stress, grey = increasing mean in situ temperature, increasing maximum in situ temperature or absence of winter cooling the Florida Department of Environmental Protection Coral Reef Conservation Program and in Broward County. This study was funded in part from contracts and agreements from the Florida Fish and Wildlife Conservation Commission (Grant Nos G0099, RM085 and RM143) and the Florida Department of Environmental Protection, Office of Coastal and Aquatic Managed Areas, as amended, through National Oceanic and Atmospheric Administration Awards NA05NOS4261187, NA08NOS4260327 and NA13NOS4820015. The views, statements, findings, conclusions and recommendations expressed herein are those of the author(s) and do not necessarily reflect the views of the State of Florida or the Department of Commerce, NOAA or any of its sub agencies. Partial funding was also through an award from the Broward County Board of County Commissioners to Nova Southeastern University. Funding for the development of HYCOM has been provided by the National Ocean Partnership Program and the Office of Naval Research. Data assimilative products using HYCOM are funded by the US Navy. Computer time was made available by the DoD High Performance Computing Modernization Program. The output is publicly available at http://hycom.org.

Compliance with ethical standards
Conflict of interest On behalf of all authors, the corresponding author states that there is no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.