Anthropogenic pressure is increasing the variety and frequency of environmental disturbance events, limiting recovery and leading to long-term declines in wild plant and animal populations. Coral reefs and associated fish assemblages are inherently dynamic due to their susceptibility to a host of disturbances, but regional-scale nuances in the drivers of long-term change frequently remain poorly resolved. Here, we examine the effects of multiple potential drivers of change in coral reef fish assemblages across 4 inshore regions of the Great Barrier Reef Marine Park (GBRMP), Australia, over 12-14 years (2007-2021). Each region had a unique disturbance history, in conjunction with and long-term changes in physical habitat variables. Phases of recovery were apparent in the years between disturbance events at all locations, but these were not long enough to prevent substantial declines in reef fish density (by 33-72%) and species richness (by 41-75%) throughout the study period. The main drivers of change in fish assemblages varied among regions, however the most rapid changes followed cyclone and flood events. Limited recovery periods resulted in temporal shifts in fish species composition from typically coral-associated to algae-associated. Most trophic groups declined in density except farmers, grazers, omnivores, and parrotfish. No-take marine reserves (NTMRs) had small and inconsistent effects on total fish assemblages, but delivered benefits for fishery-targeted piscivores. Our findings suggest that coral reef responses to local stressors and cumulative escalating climate change impacts are highly variable at regional scales., and that small NTMRs are unlikely to mitigate the impacts of increasingly frequent and climatic disturbances. Nearshore coral reefs worldwide are high-value habitats that are either already degraded or vulnerable to degradation and the loss of important fish groups. Global efforts to reduce greenhouse gas emissions must be coupled with effective local management that can support the functioning and adaptive capacity of coral reefs.

The data were collected at four inshore island locations of the Great Barrier Reef, located 10 - 30 km from the mainland coast and spanning 4.5 degrees of latitude, from 18.603S to 23.19S. Standardised underwater visual census protocols were used to survey benthic and fish assemblages at long-term monitoring sites on fringing reefs of the Palm Islands (30 sites), Magnetic Island (8 sites), Whitsunday Islands (42 sites) and Keppel Islands (20 sites) (Figure 1) between 4 and 8 times during the period 2007 - 2021. The Palm and Whitsunday Islands were surveyed in 2007, 2008, 2009, 2012, 2014, 2016, and 2018, with an additional survey in 2017 at the Whitsunday Islands, the Keppel Islands were surveyed in 2007, 2008, 2009, 2011, 2013, 2015, 2017 and 2021, and Magnetic Island reefs were surveyed just 4 times due to weather constraints, in 2007, 2012, 2016 and 2019 (Table S1). Within each island group, monitoring sites were evenly distributed among reefs that are open to fishing (General Use and Conservation Park Zones) and no-take marine reserves (NTMRs) that were closed to fishing in either 1987 or 2004.

At each of the 100 sites, five 50 m transects were deployed on the reef slope along a single depth contour between 4 m and 12 m, depending on the reef slope depth and topography at each site. Fish and benthic surveys were conducted by trained and experienced observers on SCUBA, and all species of diurnal, non-cryptic reef fish were recorded. Large-bodied, mobile fishes were surveyed using a transect width of 6 m (i.e. 300 m² survey area) by two divers swimming side by side, with a third diver laying out the transect tape behind them. Small-bodied fishes (Pomacentrids and small Labrids) were surveyed by one diver during the return swim along each transect, using a transect width of 2 m (i.e. 100 m² survey area). The same three observers conducted all fish surveys throughout the monitoring period (DHW, DMC, RDE). All recorded fish species were assigned to trophic groups, and counts were converted to density (individuals per 1000 m2) for all analyses except the Generalized Linear Mixed Model, where individuals per 300 m² was used to satisfy the requirement of integers for the preferred negative binomial distribution.

Benthic communities were surveyed using a standard point-intercept survey method by one diver during the return swim along each transect. A single benthic point sample was recorded at every 1 m graduation mark along each transect tape (i.e. 50 samples per transect). Benthic biota and substrata were classified into the following categories: live and dead hard coral with further subdivision into morphological categories (branching, tabular, digitate, solitary, massive, foliose, encrusting), soft coral, sponge, clams (Tridacna spp.), other invertebrates (such as ascidians and anemones), macro-algae, coral reef pavement (covered in turf algae), rock, rubble and sand. Additionally, for the live hard coral categories (branching, tabular, and digitate), each colony was further classified as either Acropora or ‘other’. The structural complexity of the reef habitat at each site was estimated using a simple method that gave a rank (1-5) to both the angle of the reef slope and the rugosity of the benthos for each ten-meter section of each transect.

We generated quantitative estimates of relative wave exposure at each monitoring site during each relevant cyclone. We used modelled wave height and direction data from NOAA WAVEWATCH III and the Commonwealth Scientific and Industrial Research Organisation (CSIRO) to identify which cyclones generated significant wave heights (Hs = average of top 1/3 of wave heights) of 3.5 m or more at each monitoring site. For each cyclone at each site, the distance to the nearest wave blocking obstacle was measured every 7.5 degrees around each site (fetch). These measures were weighted by the relative frequency at which cyclone-generated waves approached the site and their average magnitude. These distances were then summed and normalised to create a dimensionless index of relative cyclone wave exposure, as per previous studies.

Daily Moderate Resolution Imaging Spectroradiometer (MODIS) Level-0 data were acquired from the NASA Ocean Colour website (http://oceancolour.gsfc.nasa.gov) and converted into RGB colour images with a spatial resolution of 500 ×500 m using the SeaWiFS Data Analysis System (SeaDAS). The images were then (i) spectrally enhanced to transform them from RGB to the Hue-Saturation-Intensity (HSI) colour system and (ii) classified into three distinct water colour categories corresponding to the three optical water types (primary, secondary, and tertiary) commonly found in the GBR during the austral wet season.

We used the primary water type characterization to quantify frequency of exposure of the monitoring sites to highly turbid water from flood plumes and subsequent sediment re-suspension during the 2003-2017 Queensland summer wet seasons (December-April inclusive). The primary water type represents high turbidity, and high values of coloured dissolved organic matter (CDOM) and Total Suspended Sediment (TSS). TSS and Secchi Disc Depth (SDD) in the primary water type are typically around 18.3 ± 45.7 mg L-1  and 1.8 ± 1.8 m (mean ± 1SD), respectively. The primary water type is often associated with low salinity from flood plumes, but not always, as high turbidity can also reflect re-suspended sediment from wind and tides. We created twenty-two weekly composite images of daily images from December 1st to April 30th per wet season, to minimize the amount of area without data per image due to masking of clouds and sun glint. We assigned each weekly composite a presence/absence (0/1) value of primary water type in each pixel (500 × 500 m resolution).

Two measures of water quality were used: remotely sensed Chlorophyll-a, which provides an estimate of phytoplankton biomass and can act as a proxy for seawater nutrient concentrations, and Diffuse Kd490 (the Diffuse Attenuation Coefficient at 490 nm), which provides an estimate of turbidity. Chlorophyll-a and Kd490 composite monthly 4km data, collected using a Moderate Resolution Imaging Spectroradiometer (MODIS) satellite, from 2003-2017 were downloaded from the ERDDAP website. (Chlorophyll-a - https://coastwatch.pfeg.noaa.gov/erddap/griddap/erdMH1chlamday; Kd490 - https://coastwatch.pfeg.noaa.gov/erddap/griddap/erdMH1kd490mday). In-situ measurements of these variables are preferred as there is increased uncertainty in turbid waters, however in their absence, remotely sensed measurements can and have been used in a number of other studies. The Whitsunday Islands data for both Chlorophyll-a and Kd490 were anomalous, so they were excluded from the Whitsundays boosted regression tree analyses.

Degree Heating Weeks (DHW) values represent the accumulated thermal stress over the previous 12 weeks at a given pixel. DHW is calculated as the number of degrees above the coral bleaching threshold multiplied by the number of weeks that the elevated temperature persists. Coral bleaching is likely at 4 DHW, and this is routinely used to estimate thermal stress on coral reefs. Daily 5km data from 1998 to 2016 were provided by NOAA Coral Reef Watch (2018). The Maximum DHW reported between sequential surveys was used for each year, however if the period between surveys exceeded one year, the maximum DHW within the two previous years was used in the following year of the study.

Annual average SST and SST anomalies were calculated from multi-scale, ultra-high resolution (MUR), SST, and sea surface temperature anomaly (SSTA) data. Monthly 1km data from 2002-2017 were downloaded from the NOAA ERDDAP website (https://coastwatch.pfeg.noaa.gov/erddap/griddap/jplMURSST41mday.html and https://coastwatch.pfeg.noaa.gov/erddap/griddap/jplMURSST41anommday.html).

The data were analysed with generalised linear mixed models, non-metric multidimensional scalin,g and boosted regression trees.