Data from 1,024 bleaching surveys across the Main Hawaiian Islands were obtained from the Hawaiʻi Coral Bleaching Collaborative (Hawaiʻi Coral Bleaching Collaborative, 2021), a monitoring partnership created in response to the 2019 Hawaiʻi coral bleaching event and comprised of NOAA Fisheries Ecosystem Sciences Division; State of Hawaiʻi Department of Land and Natural Resources, Division of Aquatic Resources; Scripps Institution of Oceanography; Arizona State University; National Park Service, Pacific Island Inventory & Monitoring Network; Hawaiʻi Institute of Marine Biology; and The Nature Conservancy of Hawaiʻi. Survey and site selection methods differed slightly across organizations, and thus we included organization as a random factor in our analysis (see “Data Analysis” below). Coral bleaching (% bleaching of live coral cover) and % live coral cover were measured using either a rapid visual assessment or line transects (10, 25, or 30 m depending on the organization), and transects were quantified using non-overlapping photoquadrats analyzed through CoralNet, or a line-intercept method. Some islands were monitored using ordinal bins (bleaching levels 1-5), but we elected to use only continuous (% bleaching) data to achieve a more detailed metric, which limited our analysis to five islands: Hawaiʻi, Maui, Oʻahu, Lānaʻi, and Molokaʻi ( Figure 2A ). Each survey represents a single visit to a distinct location (“site”); some locations were spaced closer together than others (see discussion of spatial autocorrelation below). While surveys included some taxon-level data, these data were not available for all locations, and we therefore did not analyze bleaching by taxon. We also decided not to analyze taxa because when several different survey techniques are used, as in this dataset, species records can vary widely between survey techniques, while % coral cover is more comparable between survey methods (Jokiel et al., 2015). Depth was recorded for all sites but was missing for 10 observations; in these cases, we imputed depth from aircraft-based bathymetry from the Global Airborne Observatory (Asner et al., 2020a).

We extracted turbidity values from each site using a recently developed Shallow Water Turbidity (SWaT) method, which is designed to detect turbidity in nearshore shallow waters using the Planet Dove satellite constellation (3.7 m spatial resolution) (Li et al., 2022). SWaT accounts for bottom reflectance in total water-leaving radiance, and therefore improves upon historical remote-sensing turbidity algorithms, which do not address bottom reflectance and therefore may overestimate turbidity values in clear, shallow, nearshore waters where bottom reflectance is high (Lee et al., 2012). SWaT was previously validated against field turbidity monitoring data from the Pacific Islands Ocean Observing System and United States Geological Service, exhibiting an accuracy of R² = 0.83 and RMSE = 0.90 FNU (Li et al., 2022). We extracted weekly mosaics of Planet imagery for all weeks during the bleaching event (28 July 2019 to 1 December 2019) and used bathymetry processed from Sentinel-2 data (Li et al., 2021) to subtract depth-dependent bottom reflectance from water-leaving radiance. This produced water column reflectance, i.e., the in-water turbidity signal. We then used the red band of water column reflectance to derive turbidity in Formazin Nephelometric Units (FNU) based on Dogliotti et al. (2015).

We found the maximum turbidity value within a buffer of 100 m of each survey site, which was the lowest buffer that both maximized data availability while remaining appropriate to the scale of particle mixing. We used turbidity maxima because sediment runoff often occurs in sporadic events and, therefore, turbidity regimes are well-characterized by peaks (De Carlo et al., 2007). In addition, because the SWaT algorithm subtracts bottom effects it can over-correct for bottom reflectance (underestimate turbidity) in extremely shallow areas directly abutting the coast, and therefore, nearshore surveys were more accurately represented by maxima. A visual inspection showed limited evidence of artificial outliers in imagery (discrete, extreme turbidity values) near survey points, and we found high (> 0.70) Pearson's correlation between turbidity calculated from maximum and mean values, so determined that using maxima was appropriate. All analysis was conducted using Python and the gdal library on Linux command line.

To attribute turbidity to urban sources, we created a binary variable for urban coastline using a dataset of ‘urbanized areas and urban clusters’ from the U.S. Census Bureau (2015) TIGER database. Survey points were snapped to the nearest coastline and considered urban if listed as ‘urban’ by the Census; each survey point was then manually checked against true-color satellite imagery to ensure that classification based on 2015 Census data was still accurate.

We controlled for additional bleaching drivers determined from prior literature: thermal stress, sea-surface irradiance (hereafter “irradiance”), wave power, total effluent, nitrogen runoff, depth, and distance from shore. We also examined % live coral cover as a possible confounding variable, given that high coral cover may indicate low prior bleaching and, perhaps, low vulnerability.

We derived thermal stress as both 1) SSTa and 2) Degree-Heating Weeks (DHW) from NOAA Coral Reef Watch Operational Daily Near-Real-Time Global 5-km Satellite Coral Bleaching Monitoring data (2019), obtained via ERDDAP (Simons, 2020). Since coral survey dates were inconsistent, ranging from 20 August 2019 to 7 December 2019, we extracted daily thermal data between a start date of 1 August 2019 (1 week before the Main Hawaiian Islands mean SST surpassed the regional bleaching threshold of 28.5° C) and each unique survey date. SSTa was measured as the mean daily SSTa for that period, and DHW as cumulative DHW at the survey date. SSTa and DHW were linearly related, so we selected only SSTa for analysis. We measured irradiance as the maximum Photosynthetically Available Radiation (PAR, mol m^-2 day^-1) at each site during the bleaching event, derived from the MODIS Aqua satellite (NASA, 2019), which has a maximum spatial resolution of 4 km. As with temperature, we accounted for variable survey dates by fetching data between 1 August 2019 and each individual survey date. We then calculated PAR as the maximum pixel value at each point from eight-day PAR composites.

‘Wave power’ (kW m^-2) was the long-term mean of wave power (2000-2013) at each site, obtained from Kappel et al. (2017a). In this dataset, long-term means were based on hourly wave data from a SWAN (Simulating WAves Nearshore) model, which were converted to maximum daily wave power, then averaged by month from 2000-2013. While these data precede 2019, they represent the best available data for our study region and serve as a relative index of wave exposure between sites. ‘Total effluent’ (gal km^-2 d^-1) was likewise derived from Kappel et al. (2017b), based on effluent flux from Onsite Sewage Disposal Systems (OSDS). In this dataset, OSDS flux was attributed to Tax Map Key (TMK) parcels (Whittier and El-Kadi, 2014), which were aggregated into 500-m pixels and added as total flux within a 1.5-km radius of each cell. ‘Nitrogen runoff’ (gal km^-2 d^-1) was obtained from the same dataset.

Finally, we found ‘distance from shore’ as the orthogonal distance between each point and the nearest landmass. All geospatial processing was conducted using the sf package in R (Pebesma, 2018) and visualized in QGIS. Final independent variables for analysis were turbidity, urban area (binary), SSTa, distance from shore, depth, total effluent, nitrogen runoff, live coral cover, wave power, and irradiance ( Table 1 ). We also tested interactions between all variables and SSTa, and the three-way interaction of turbidity, urban area, and SSTa given our interest in the source (urban v. non-urban) of turbidity.

We used a Bayesian hierarchical spatial model to estimate the effect of turbidity and other environmental variables on coral bleaching. Given that our dataset contained a high number of zeros (18% of data showed no bleaching), we used a Zero-Altered approach (Zuur and Ieno, 2016) where bleaching was modeled first as a binary outcome of ‘bleaching occurrence’ (whether bleaching was present/absent in the survey), and data was then subset to only nonzero bleaching surveys to assess ‘bleaching severity’ (% bleaching of live coral cover). A Bernoulli distribution was used to fit bleaching occurrence and a beta distribution for bleaching severity, as these distributions are best suited to, respectively, presence-absence and proportional data (Zuur and Ieno, 2016). Bleaching was modeled using logit link functions as:

where Y i o c c and Y i s e v are bleaching occurrence and severity, E denotes fitted bleaching values, X_i is a covariate matrix of independent variables, β is a vector of regression parameters with α intercepts, μ is the mean of the beta component, π is the probability of bleaching absence vs. presence, ϕ summarizes shape parameters for the beta distribution, and W_i represents the models’ spatial structure, described below.

We standardized all continuous independent variables by subtracting the mean from each observation and dividing by the standard deviation. Beta distributions only take proportional data < 1, so we converted our response variable where Y_i = 1 to negligibly lower values using Y_i = (Y_i * (N-1) + 0.5)/N, following Zuur and Ieno (2016). We assessed all variables for collinearity using Pearson’s correlation index; a correlation higher than 0.70 existed between ‘nitrogen’ and ‘total effluent’ variables (p = 0.98, P < 0.001) so we retained only ‘total effluent’ for analysis ( Figure S1 ). We also used Point Biserial Correlation to test for correlation between categorical and continuous independent variables and confirmed that all p < 0.70. Preliminary plots of independent variables ( Figure S2 ) indicated high positive skewness in distance from shore, depth, and wave power. Distance from shore was log-transformed, which is advised when distributions display extreme positive skewness, and wave power and depth (milder skew) were square-root-transformed. We also removed two outliers in total effluent, which represented exceptionally high values ( Figure S2 ). Outlier removal and variable transformations had no impact on model fit > 2 DIC but improved homoskedasticity in residual-covariate plots ( Figure S3 ).

A variogram of Pearson’s residuals in initial models indicated that spatial autocorrelation existed in our dataset. We therefore modeled a spatial random effect using Stochastic Partial Differential Equations (SPDE) based on a Matern variance-covariance matrix (Lindgren and Rue, 2015; Zuur et al., 2017). To drive SPDE, we created a mesh ( Figure 2B ) with a nonconvex hull (convex = -0.05, concave = 0.05), which is best suited to well-separated islands, and specified maximum node distance between 1 and 10 km. We used a histogram of between-plot distances to determine range and sigma priors of, respectively, 5 km and 1 km, at 0.05 probability. We also conducted a sensitivity analysis of node distances, ranges, and sigmas, and found that variable influence was not substantially changed by priors.

To account for slight differences in survey method by organization, we included ‘organization name’ (identity of the organization collecting data) as a random factor in our analysis. Random effects are appropriate when correlation is expected between observations in the same group (organization), but when the relationship between independent and dependent variables is still expected to carry across all group levels (Harrison et al., 2018). Since we did not expect ‘organization’ to change the environment-bleaching relationship, but produce variable intercepts, a random factor was appropriate. We also included island as a fixed effect, given that islands represent different oceanographic, geomorphological, and human development regimes on reefs.

Final variables in the model were determined through backward stepwise selection, where we removed variables that did not increase model fit by 2 DIC, i.e., when only influential variables remained ( Table S1-S2 ). We tested for homogeneity of variance by comparing residuals to all covariates, and to fitted values ( Figure S3A–H ). We also tested for normality of residuals ( Figure S3I ), and compared fitted to observed values, visually determining that the relationship was linear but not overfitted. Following Donovan et al., 2020; 2021), we considered variables influential if their 80% credible interval did not cross zero, but also assessed the 95% credible interval to identify highly influential variables (note: frequentist p-values do not apply in Bayesian statistics, but we will use the term “significant” to refer to these influential variables, for the sake of interpretability). We used penalized complexity priors for random effects following Fuglstad et al. (2019); default vague Gaussian priors for fixed effects (mean = 0, precision = 0.001), and log-gamma priors for hyperparameters (mean = 1, precision = 0.1). All analysis was conducted using the R-INLA (Rue et al., 2009).

Turbidity was a negative driver of bleaching severity (80% credible interval -0.195 < β < -0.049), and thus, turbidity decreased percent bleaching in live coral cover ( Figure 4B ). In contrast, turbidity did not affect bleaching occurrence (binary variable of bleaching presence/absence) ( Figure 4A ). The interaction between turbidity and heat stress (SSTa) was positive for both bleaching occurrence (0.098 < β < 1.121) and severity (0.139 < β < 0.373). Thus, thermal stress suppressed the protective effect of turbidity on bleaching severity, i.e., turbidity had a weaker effect at higher temperatures ( Figure 5 ). In urban areas, turbidity had a largely flat effect on bleaching severity, though in non-urban areas, the relationship between turbidity and bleaching severity reversed to weakly positive (interactive effect of urban turbidity -0.480 < β < -0.143). Urban turbidity also had a negative effect on bleaching occurrence (-1.137 < β < 0.462), but this effect was not significant ( Figure 4A ). Urban areas, when treated in isolation without interaction, had a highly positive effect on bleaching occurrence (0.499 < β < 1.894) and no effect on bleaching severity. Therefore, the urban environment alone did not create a protective effect on corals during the bleaching event, but the presence of urban turbidity may have indicated higher coral resistance or resilience to high temperatures (note: covariates were standardized; see Table S3 for β relative to standard deviations).

SSTa had a positive effect on bleaching occurrence and severity, though this effect was highly significant (0.886 < β < 1.848) for occurrence only. Therefore, SSTa was a major driver of whether corals bleached, but not specific levels of bleaching. In contrast, bleaching severity showed a significant dependence on the interactions between temperature and other environmental covariates such as turbidity, urban turbidity, and wave power ( Figure 6 ; Figure S4 ). The interaction between wave power and SSTa had a negative effect on bleaching severity (-0.282 < β < -0.016) and therefore wave action reduced bleaching levels under increasing thermal stress.

Beyond temperature, variables associated with bleaching occurrence were % live coral cover (0.189 < β < 0.613) and distance from shore (-1.087 < β < -0.430). Areas with higher live coral cover showed a higher propensity to bleach, though this effect was weak. Corals located further from the coast were more likely to show low to no bleaching. In contrast, % live coral cover had a weakly negative relationship with bleaching severity (-0.113 < β - 0.009), as increasing live coral cover was correlated with a slight decrease in bleaching severity. Bleaching severity was higher in deeper reefs (0.145 < β < 0.286). Sea-surface irradiance was not a significant variable for either bleaching occurrence or severity.

Including island as a fixed effect dramatically reduced model DIC values, and therefore island had a strong influence on bleaching occurrence and severity. With the Big Island as a reference value, all islands showed higher bleaching occurrence and severity, although Oʻahu showed the strongest, positive effect, and Lānaʻi the weakest effect on both bleaching metrics ( Figure 6 ). In addition, including organization name (identifier for the organization conducting surveys, which was a proxy for survey method) as a random effect substantially reduced model DIC values (-40 DIC) and thus increased model fit ( Table S1–2 ).

Turbid coral bleaching refugia have been documented worldwide, but in contrast, so have heavily damaged reefs under severe sedimentation, toxin exposure, effluent, dredging, and other turbid-area stressors (Brodie et al., 2012; Cacciapaglia and van Woesik, 2016; Carlson et al., 2019; Sully and van Woesik, 2020). By coupling high-resolution remote sensing with in-situ field observation, we demonstrate that turbidity may provide protection from bleaching severity but increasing temperature stress may weaken or reverse this effect, particularly in rural areas. Across the Main Hawaiian Islands, turbidity ameliorated bleaching severity overall, but thermal stress counteracted this effect, as indicated by a positive turbidity-SSTa interaction. In rural areas in particular, turbidity failed to improve coral bleaching under increased SSTa, while in urban environments, turbidity still mildly suppressed bleaching despite climbing temperatures ( Figure 5 ). This pattern persisted when bleaching was measured as bleaching occurrence (presence/absence): turbidity decreased bleaching occurrence, though this relationship weakened as temperatures rose, especially in rural areas. Unlike for bleaching severity, these effects were not significant (80% credible interval crossed 0), likely because bleaching occurrence was overdetermined by SSTa, which was by far its most significant driver.

These findings suggest that temperature can override the protective benefits of turbidity, but that marginal reefs (chronically silty, urban areas) may have unique features that confer resilience. One possible reason for this pattern is that urban, turbid areas in Hawaiʻi may contain stress-tolerant corals. Urban hubs in Hawaiʻi have experienced prolonged disturbance, where surviving corals may have been selected over decades of intensive human contact. Stress-tolerant communities are common in urban, turbid areas, and have been observed in Singapore, Brazil, and Australia (Mora, 2008; Guest et al., 2016a; Heery et al., 2018; Loiola et al., 2019). In these regions, weedy or stress-tolerant communities have survived and/or adapted to urban trends like early grazer removal, heavy metal deposition, macroalgal proliferation, artificial habitat from urban structures, and genetic isolation (or, conversely, novel corridors) along the built coastal form (Goodkin et al., 2011; Darling et al., 2012; Heery et al., 2018). Rural areas, in contrast, represent the new development margin in Hawaiʻi, where earthmoving is recent and increasing. In these developing but not yet “urban” areas, coral community shifts may not have taken full course and stress-tolerant species may be uncommon or absent. Therefore, the combination of turbidity and heatwaves in rural areas may be particularly damaging. This aligns with findings by Asner et al. (2022), which linked coral mortality during Hawaiʻi’s 2019 bleaching event to new development (conversion of “natural landscape to human-built environments” over five recent years).

Coral taxa associated with urban landscapes in Hawaiʻi also may have experienced more sustained pressure that begets trophic plasticity over time. Heterotrophy has been noted in turbid reefs worldwide (e.g., De'ath and Fabricius, 2010; Morgan et al., 2016), conferring resistance to bleaching through higher photosynthate translocation during heatwaves (Tremblay et al., 2016) and carbon subsidies when symbionts are lost (Grottoli et al., 2006). However, heterotrophic taxa have been recorded in remote and urban areas alike (Zweifler et al., 2021), and further research is required to delineate the trait space of corals surviving in urban versus rural, sedimented habitats, as well as linkages between these traits and repeat bleaching susceptibility. Importantly, our results show that living in an urban environment, absent turbidity, does not confer bleaching resilience, suggesting that bleaching tolerance is limited to taxa thriving in murky segments of the urban landscape.

Previous research shows there are considerable coarse- and fine-scale morphological differences between species and phenotypes that impact coral response to sediment (Jones et al., 2016). For example, branching morphologies have been linked to higher sediment tolerance due to high surface inclination (Jones et al., 2019) but lower bleaching tolerance (Darling et al., 2012). Sediment clearing mechanisms are diverse and complex, as are morphologies designed to maximize photosynthesis (Stafford-Smith and Ormond, 1992; Guest et al., 2016b). Traits of corals surviving under low light, high sediment, and high temperature conditions along gradients of remote, urban, and developing coastlines is a critical area of research. A large portion of sites considered in our analysis did not record taxonomic data, and therefore, we were unable to evaluate taxonomic explanations for our results. Future studies should use innovative methods like rapid Structure-from-Motion photogrammetry to capture geospatial patterns in morphology over large areas, and over time (Fukunaga et al., 2020).

Another possible reason for urban refugia is that urban and rural sources of turbidity differ. The types of dissolved and particulate matter associated with urban areas may be, counterintuitively, less harmful to corals than rural runoff at high temperatures. Urban turbidity is often associated with point-sources of pollution like wastewater treatment plants and stormwater outfalls (Carlson et al., 2019), whereas rural areas are linked to more diffuse, nonpoint sources of pollution like unpaved roads (Ramos-Scharron and MacDonald, 2005), agricultural and golf course runoff (Delevaux et al., 2018), erosion from scalded hillsides (Bartley et al., 2014), soil disturbed by feral ungulates (Tribble et al., 2016), and new construction areas, i.e., future urban sprawl (Maragos, 1993). Mills and Sebens (1997) showed that corals reject “clean” sediments. “Natural’, terrigenous sources of turbidity may be most damaging to corals in Hawaiʻi, where lateritic, volcanic soils have small grain sizes (Storlazzi et al., 2015) that are energetically expensive for coral polyps to clear, reducing energy stores that are severely needed during thermal stress (Erftemeijer et al., 2012); however, these particles also take longer to settle out of suspension (Bartley et al., 2014). While soil runoff and fine siltation certainly occurs in urban areas, it might represent a more dominant size fraction in steep, erosive watersheds with few retaining structures, and where new construction is accelerating. Future research in Hawaiʻi should investigate sediment constituents and settling rates alongside bleaching severity to determine which sedimentation profiles are most associated with coral bleaching, mortality, and resilience.

Point sources of pollution, which concentrate in urban areas, are better regulated than nonpoint sources in the United States, the country of our study region (Humenik et al., 1987). Several case studies demonstrate the benefits of legally mandated, infrastructural upgrades on reef health in turbid urban waters like wastewater outfalls zones in Hawaiʻi. For example, after a municipal sewage system in "Māmala" Bay, Hawaiʻi received treatment upgrades in 1977, reefs in the bay showed no significant evidence of sewage impact (Grigg, 1995). Likewise, sewage upgrades in 1978 in Kāneʻohe Bay resulted in a doubling of coral cover from 12% in 1971 to 26% in 1983 (Hunter and Evans, 1995). If rural areas are more difficult to regulate, or last in line for infrastructural improvements, they may drive a higher toxic load on reefs from facilities like leaking cesspools/underground storage tanks (Yoshioka et al., 2016; Abaya et al., 2018). Finally, our results align with a recent meta-analysis on bleaching susceptibility in human-dominated and isolated reefs (e.g., atolls), which provide limited evidence for “remote refugia” but indicate some resilience in human-dominated habitats (Baumann et al., 2022). In particular, our results suggest that regions that are not yet developed but show symptoms of new development (non-urban turbidity) may be particularly susceptible to coral bleaching.

Our results do not imply that urban turbidity creates a ‘safe haven’ for coral reefs. We did not evaluate coral disease, diversity, fecundity, size classes, calcification rate, or other elements of coral health in our analysis, and it is possible that corals in silty, urban areas are stress-tolerant but sparse, homogenous, or otherwise represent ecosystem phase shifts. However, our results support past case studies (e.g., Guest et al., 2016b) by indicating that urbanization can counterbalance thermal stress during bleaching events. We further add that, in rural areas, turbidity does not counteract heatwaves, and runoff from new development, agriculture, unpaved roads, and nonpoint sources of pollution may still pose a critical threat. In Hawaiʻi, rural areas may represent a combination of unchecked earthmoving, new dredging, and aging infrastructure, which require urgent intervention.

A recent study indicated that mortality is highest when thermal stress and dredging conditions (high suspended sediment, low light) occur together (Fisher et al., 2019). While corals bleached less under low/moderate light reduction, severe light limitation acted synergistically with thermal stress to increase mortality. We likewise find that turbidity only acts antagonistically against bleaching under certain conditions (low heat stress). More research is needed on the nuanced thresholds of light, irradiance, and heat that alter interactions between stressors for various coral species.

Our findings are inconsistent with one hypothesis of turbid refugia—that sediment shades corals from excessive sunlight (Anthony et al., 2007). Equivalent levels of rural and urban turbidity should have the same effect on coral bleaching if protection comes from light-shading alone. Instead, in our study, turbidity associated with different land uses affects bleaching differently. At low levels of SSTa, rural areas show a slightly steeper decline in bleaching severity than urban areas. Conversely, at high SSTa, an equivalent increase in turbidity provides less bleaching protection in rural areas compared to urban areas ( Figure 5 ). However, turbidity levels may be relatively moderate in the Hawaiian Islands compared to other regions where turbid refugia have been studied. For example, Morgan et al. (2016) found that turbidity in Paluma Shoals, Australia caused 100% attenuation of Photosynthetically Active Radiation (PAR) at 4 m water depth, relative to turbid reefs in Hawaiʻi like Puhi Bay where we previously measured ~72% attenuation of PAR at 5 m depth in preliminary experiments. Rural turbidity may still protect corals from excessive irradiance during bleaching events in regions where turbidity levels are high (van Woesik et al., 2012; Cacciapaglia and van Woesik, 2016; Morgan et al., 2017) though, as discussed, suspended sediment and excessively low light can also induce bleaching and mortality in corals (Fisher et al., 2019).

When estimating the impacts of SSTa alone, we found that SSTa was an extremely strong predictor of bleaching occurrence, but not severity. However, the interactions between SSTa and local environmental variables (turbidity, urban turbidity, and wave power) had a very strong effect on bleaching severity. Specifically, under increasing SSTa, turbidity exacerbated bleaching severity, but urban turbidity and wave power reduced bleaching severity ( Figure 4B ). Overall, our results show that temperature, expectedly, is the most important factor determining whether corals bleach, but SSTa interactions with local environmental variables determine whether bleaching is minor or severe. Often, studies do not differentiate between bleaching occurrence and severity, and our results highlight additional insights gained through this approach. Irradiance showed no significant relationship to coral bleaching occurrence or severity, perhaps because our remotely sensed irradiance metric captured light on the ocean surface but not at depth, and corals in our analysis experienced a wide range of in-water light-scattering conditions. Indeed, since turbidity is a light-based metric, turbidity likely captured some variability in irradiance at depth. In addition, the effect of SSTa, irradiance, and other environmental variables may have been moderated by the lower resolution of those data (e.g., 5 km for SSTa, 4 km for irradiance, 500 meters for total effluent) relative to turbidity; nonetheless, these were the best data available.

Our results highlight the importance of SSTa-environment interactions, which is consistent with growing literature on the impact of local drivers on coral reef health (Graham et al., 2015; Cinner et al., 2016; Ellis et al., 2019; Donovan et al., 2021). Local variables (e.g., wave action, upwelling, internal waves, tides, mangrove refugia, and other features) alter coral exposure to thermal stress and ability to cope with and recover from bleaching (West and Salm, 2003). In two global analyses using Reef Check data, interactions between SSTa and local conditions including turbidity, wave power, and macroalgae were found to mediate bleaching severity (Sully and van Woesik, 2020; Donovan et al., 2021). In our analysis, wave power negatively influenced SSTa effects on bleaching severity, perhaps because waves increase mass transfer of oxygen radicals from coral tissues (Nakamura and Van Woesik, 2001), thereby reducing host stress. Waves may also have caused ocean cooling at finer scales than the thermal signal detected by satellites. Live coral cover was linked to increased bleaching occurrence but decreased bleaching severity. One possible explanation is that high coral cover increases the probability of observing any (even minor) bleaching on a survey, but also indicates areas where corals survived prior bleaching in Hawaiʻi (high bleaching resilience). Overall, interactions between local environment and regional heatwaves is a critical area of ongoing research, and advances in remote sensing can enable a granular, scalable understanding of these dynamics.

Including island as a fixed effect increased model fit, and in particular, three islands (Maui, Molokaʻi, and Oʻahu) showed increased bleaching severity relative to Hawaiʻi. Molokaʻi is largely undeveloped but exhibited the highest turbidity of any island ( Table 1 ), likely because numerous surveys were located on the heavily sedimented south shore of the island where terrigenous sediment is deposited from several drainage basins (Ogston et al., 2004; Jokiel et al., 2014). Oʻahu, in contrast, is highly urbanized with variable turbidity and low coral cover proportional to substrate area (Asner et al., 2020b), though Oʻahu surveys displayed a mean 43% live coral cover and surveys may therefore have focused on higher-coral areas. Maui represents a mixture of urban and rural areas experiencing sewage discharge and bottom sediment resuspension in urban centers (Laws et al., 2004) and increasing development and agriculture in remote regions (Oleson et al., 2017). Given this complex array of factors, it is difficult to interpret island-based results without additional investigation, but it is possible that site selection may have played a role in island influence (e.g., excluding the north and west sides of Molokaʻi).

Finally, it should be noted that including ‘organization name’ as a proxy for survey method substantially increased model fit in all cases, pointing to the need to control for survey biases in coral research. Nonetheless, when using broad metrics like percent bleaching, past research has shown that coral survey techniques yield comparable results (Jokiel et al., 2015).

This study is an early step in using high-resolution satellite monitoring to characterize possible turbid refugia. We considered broad categories of turbid landscapes (urban and rural), which should be refined in the future to include a full spectrum of land-use classes such as agriculture, forestland, grassland/scrub, suburban and residential areas, and others. In addition, recent research points to the importance of temporality in turbidity regimes (persistent, fluctuating, and transitional turbidity) in dictating coral health (Zweifler et al., 2021). Planet Dove satellites, which offer a daily revisit frequency, can help characterize differences in the timing of turbidity, and identify which temporal regimes are linked to bleaching protection or harm. Importantly, the methods we highlight here may be highly useful to local managers. High-resolution satellite data can help coastal decision-makers identify specific plumes, resuspension areas, and other sources of turbidity to address runoff on a case-by-case basis. For example, where heavy plumes are linked to bleaching refugia, managers may prioritize coral protection by addressing persistent stressors like pathogens and disease, while frequently monitoring for algal overgrowth. Turbid areas associated with coral decline, on the other hand, should be prioritized for complete remediation. High-resolution satellite data can also be used to help refine existing models for predictive purposes, and scaled to local suspended sediment concentrations to monitor sediment against known coral health thresholds.