Introduction

Marine ecosystems provide a wide range of ecological, economical and socially valuable resources to people, through services such as tourism, food supply and cultural heritage (Beaumont et al. 2007; Barbier et al. 2011; Hoegh-Guldberg et al. 2015). As a consequence of this high cultural and economic importance, human populations in coastal regions are typically very large, with approximately 41% of the world’s global population and roughly 60% of the world’s ‘megacities’ found within 100 km of the coast (Martínez et al. 2007). The combination of high population densities and easy access to markets or trade routes often leads to intense pressure on surrounding marine ecosystems via chronic overfishing, alongside land-based pollution and progressive coastal urban development (Worm et al. 2006; Martínez et al. 2007; Cinner et al. 2018).

The size-selective nature of many fishing methods often targets adult life stages and larger-bodied species, which can cause pronounced changes in predator–prey dynamics (Bascompte et al. 2005; Robinson et al. 2017). Within an overexploited or unmanaged fishery, where extraction exceeds replenishment, the loss of high numbers of large and adult fish is capable of producing destabilising trophic cascades (Pace et al. 1999). Such effects can cause a resultant shift in size spectra (i.e. the distribution of community-wide body sizes) towards a dominance of smaller fish with r-selected life-history traits (Jennings and Blanchard 2004; Robinson et al. 2017). Simultaneous reductions in functionally important reef-associated groups, such as the herbivorous scrapers, grazers and browsers, can cause additional indirect changes to the benthos because these species tend to increase reef resilience and can help to prevent post-disturbance algal phase shifts (Hughes et al. 2007; Graham et al. 2013, 2017; Hempson et al. 2018).

In response to these potential broad changes in fish community demographic structure, there has been a recent global push towards increasing the coverage of large marine protected areas (MPAs) or reserves. This push is further linked to meeting international Aichi and UN Sustainable Development Goal targets to effectively protect at least 10% of coastal and marine areas and the biodiversity that they contain by 2020 (Thomas et al. 2014). However, small-scale MPAs are also important for meeting this target and can have significant benefits for the conservation of fisheries and habitat if designed and managed well (Lester et al. 2009; Thomas et al. 2014). Despite the Philippines being considered to hold the highest global marine biodiversity, only 2.7–3.4% of the coral reefs are currently protected (Weeks et al. 2010), and small-scale MPAs are the most popular format for local management in this region. Philippine MPAs have a median size of only 0.15 km2 (Carpenter and Springer 2005; Weeks et al. 2010), 95% of which are community-managed (i.e. designated under local or municipal-level ordinances) and are therefore typically without government funding (Weeks et al. 2010).

While country-wide coral cover in the Philippines has increased overall between 1981 and 2010 (Magdaong et al. 2014), most small-scale fisheries are declining and are labelled as unsustainable (Muallil et al. 2014a, b). These changes in Philippine-wide reef fisheries and coral cover take place in the context of intense localised fishing, pollution and development pressures (Burke et al. 2012; Selgrath et al. 2018), alongside a range of large-scale environmental impacts, including coral bleaching and frequent typhoon, tsunami and earthquake damage (Shimozono et al. 2015; Anticamara and Go 2017; Jamero et al. 2017; EM-DAT 2020). The major 1997–1998 and 2010 bleaching events in particular caused widespread coral mortality and the loss of much fish habitat across the reefs of the Philippines (Goreau et al. 2000; Magdaong et al. 2014).

We analysed reef monitoring data collected between 2000 and 2011 in the central region of the Philippines to assess the effect of a network of fully protected small-scale MPAs on long-term reef community dynamics following major environmental disturbance and intense local resource extraction. We describe the rate and direction of change in family-level fish abundance, size and trophic group composition following disturbance as well as changes in coral cover and benthic composition. We ask whether the presence of such small-scale local MPA management in this area was sufficient to temper the expected declines from intensive fishing and disturbance, and whether certain fish sizes or trophic groups have been more or less successful.

Materials and methods

Study area and site information

The Danajon Bank reef system is located within the central Visayas region of the Philippines and covers a total area of 271.7 km2 encompassing the Cebu, Bohol and Leyte provinces (Fig. 1).

Fig. 1
figure 1

The location and characteristics of the eight marine protected area (MPA) sites along the Danajon Bank reef system within the Philippines. Site characteristics include measures of management and compliance based on (Alcala et al. 2008; Yasué et al. 2012). *The Coastal Conservation & Education Foundation (CCEF) management scoring system ranges from 0 (lowest) to 40 (highest). Refer to (White et al. 2006) for further details of score calculation and grading. In 2015 the Philippines changed their system for rating the effectiveness of MPAs to the ‘Management Effectiveness Assessment Tool’ (MEAT scores available at: https://database.mpasupportnetwork.org), however, we use the previous Coastal Conservation & Education Foundation (CCEF) system here for consistency over a longer period of assessment relevant to this study timeline. The islands of Bohol and Cebu (around the population hub of Cebu city) are additionally shown for geographical context

Full size image

Assessments were conducted inside and outside of eight community-led MPAs, within five municipalities of Bohol (Fig. 1). The sites are all small (ranging between 10 and 50 ha) and range in their age and level of enforcement, but are all fully no-take, with the oldest site being designated in 1995 (although the earliest baseline survey data begins in 1999).

These sites are all considered ‘well enforced’ according to the Coastal Conservation and Education Foundation (CCEF) scoring system, with scores ranging from 33 to 40 (of a 40 maximum), and were chosen as they held complete information on the level and duration of MPA management, have complete and consistently recorded replicate datasets using a standardised protocol (Samoilys et al. 2007) and are spread across the whole Danajon Bank region of Bohol.

Unprotected areas across this region are subject to fishing pressure from multiple gear types (including blast fishing, fish corrals, gleaning, hook-and-line, nets, poison fishing, skin diving and traps), which target a mixed-species fishery, primarily for local subsistence use (Selgrath et al. 2018). Alongside localised fishing pressures, there is extensive pollution, sedimentation and coastal development (primarily from the population hub of Cebu, ~ 20 km to the North West of northern Bohol.

Survey data

Biannual transect surveys were conducted along the Danajon Bank in the dry (March–May) and wet seasons (September–December) from 2000 to 2011 at eight sites. All surveys were conducted on reef flats at depths of 4–10 m, along 50 m fixed transects in a hierarchical nested design with four–five replicates per site, per management type (inside or outside the MPA), per season (wet or dry), per survey year. This gave a range of 17–51 replicates per site per management type for benthic transects and 20–63 replicates per site per management type for fish surveys, depending on length of enforcement at the sites (8–12 years) and the ability to survey during wet season (Appendix 1). For each transect we assessed benthic composition, family-level fish abundance and fish total length. Samoilys et al. (2007) showed wet and dry season surveys in this region do not contribute significantly to fish variation among sites; therefore, data were pooled by year for all fish and benthic surveys to increase within-replicate sample sizes.

Fish were assessed in belt surveys, recording abundance and estimated size (in cm) of all diurnally active non-cryptic reef-associated fish families within a 5 × 5 m viewing window. Size estimation and family ID were conducted visually by trained observers. Fish length estimates were validated against in situ bars of known length during training and then periodically tested through field surveys. Carangidae (jacks/trevallys) were excluded from the analysis of biomass across time as large mobile piscivores such as these are often strongly over-estimated in localised visual censuses and typically over-inflate and skew results (Williams et al. 2015; Robinson et al. 2017).

Following surveys, all observed families were classified into groups (Appendix 2) according to their broad commercial food and aquarium trade value (low, medium or high), feeding type (corallivore, detritivore, herbivore, invertivore, piscivore, planktivore, scavenger or spongivore) and their degree of reef association (low, medium or high) using data from Fishbase (Froese and Pauly 2018) and related literature (Christie et al. 2006; Samoilys et al. 2007; McClanahan 2014). Dominant trophic feeding strategy was chosen for each fish family if multiple strategies existed across species or life stage.

Fish biomass estimates were calculated for each family using the standard weight–length relationship (Froese 2006), where W = a Lb, with the values of constants a and b based on species-specific calculations obtained from Fishbase (Froese et al. 2014; Froese and Pauly 2018) using the Rfishbase() web-interface package in R. Given that the fisheries data were collected at the family taxonomic level, we pooled size-based growth values for all species known to occur in this region (Bayley et al. 2019) and then used the mean values for the family averaged groups.

The composition of benthic substrate cover was assessed using Line Intercept Transects (English et al. 1997), with two–three transect replicates per site, per management type (inside or outside the MPA), per season (wet or dry), per survey year. Broad substrate groups were recorded, including seven coral growth morphologies (branching, encrusting, foliose, massive, mushroom, submassive and tabular).

Data analysis

A mixed-effects linear regression, fit by restricted maximum likelihood, was used to assess the significance of MPA management and length of enforcement in explaining variations in total fish biomass. ‘MPA management’ (In/Out) and ‘length of enforcement’ (years) were set as fixed factors, and ‘MPA site’ (n = 8) was set as a random effect with random slope and intercept included (following log likelihood ratio test assessment). A log transformation was applied to biomass to achieve a normal distribution of model residuals, and the plotted residuals were checked for homoscedasticity prior to using the results of the model. Following the overall analysis, the data were then split into eight trophic groups, with individual models run for each group (with site remaining as a random factor with random slopes).

To test for overall changes in size spectra through level of management and time, linear regressions were calculated for (log-transformed) fish density (250 m−2 transects) inside and outside protected sites, for the period 2000 to 2011. Data were pooled across eight sites and all families, and these dates were used in order to cover the largest time span where all sites were monitored consistently with adequate sample replication. Fish counts were separated into 10 cm (total length) size bins, analysing each separately. The 40 cm + size class includes all size classes 40–100 cm, due to limited abundances within these classes. The linear regressions were conducted using the R function lm (). Differences in fish total length frequency distributions between MPA management (In/Out) were tested for the eight pooled sites at 0, 5 and 10 years since MPA designation using the permutational ‘sm.density.compare’ function in the R package ‘sm’ (Bowman and Azzalini 2018).

Further, mixed-effects linear regressions were used to assess the significance of MPA management and length of enforcement in explaining variations in live coral cover through time. ‘MPA management’ and ‘length of enforcement’ were again set as fixed factors, and ‘MPA site’ (n = 8) was set as a random effect with random slope and intercept included (following log likelihood ratio test assessment). A square root transformation was applied to percentage cover values to achieve a normal distribution of model residuals, and the plotted residuals were checked for homoscedasticity prior to utilising the results of the model. While the error distribution of percentage cover data is technically binomial in nature, in practice the range of cover values did not exceed 0.2–0.8, residuals fitted a normal Gaussian distribution, and model predictions were bounded within the expected 0–100 range, so a simpler linear model was chosen over a generalised linear model. Site-level differences were investigated using the function Lmlist(), dividing the pooled data into eight groups, based on the random factor of site.

Linear mixed-effects models were generated using the lmer() function in the R package lme4. Significance (p) values for the selected models were calculated using Kenward–Roger standard errors and degrees of freedom. Confidence intervals are calculated at the 95% level.

Community analysis

To test for differences in coral trait composition between sites (n = 8), MPA management (n = 2) and years since designation (n = 16), we used permutational multivariate analysis of variance (PERMANOVA) over 9999 permutations on a log-transformed Bray–Curtis dissimilarity matrix using the R ‘vegan’ package functions adonis() and vegdist() (Oksanen et al. 2017). Principal co-ordinate analysis was then conducted on benthic community structure and coral growth morphology including years since designation, for all sites, with MPA management type and site name as factors.

Forest plots were used to show effect sizes for differences in fish community through time and between management. Plots were based on natural log-transformed response ratios (LnRR) of fish family mean total abundances, pooled for all MPAs (n = 8). The LnRR is a robust measure of relative differences within a community and can be widely compared with larger meta-analyses (Goldberg et al. 1999; Cote et al. 2001). Samples inside and outside of protection were analysed separately, with zero and single individual counts removed. We further analysed the responses to protection for the year 2011 alone, following 10–15 years of enforcement. Forest plots were created using the escalc () function within the R Metafor package.

All analyses were carried out using R: version 3.4.1 (R Core Team 2016).

Results

Effect of enforcement and time on fish communities

Fish biomass

Logged average fish biomass (grams per 250 m2) from N = 570 surveys across 8 sites did not significantly change over the study period inside of MPAs (slope estimate = − 0.03, t value = − 0.72, CI = − 0.10, 0.05), from an initial intercept estimate of 6.43 kg biomass/250 m2 (t value = 28.86, CI = 4.19, 14.19 kg/250 m2). Fish biomass did, however, decrease significantly over this period outside of MPAs (slope estimate = − 0.12, t value = − 4.17, CI = − 0.15, − 0.07) from an initial intercept estimate of 6.43 kg biomass/250 m2 (t value = − 1.35, CI = 3.50, 11.85 kg biomass/250 m2). Total fish biomass did not differ significantly inside and outside of MPAs at the time of designation.

The majority of the variance (44% of the 55% explained) within the mixed-effects model was among-site variation (Appendix 3), indicating an array of responses to protection, but an overall trend towards losses outside of MPAs. The Danajon Bank region-wide change in biomass outside of MPAs equates to a loss of ~ 802.9 g fish biomass per transect, or ~ 3.2 g m−2 yr−1 over the time period 2000–2011 (Fig. 2).

Fig. 2
figure 2

Logged fish biomass (g m−2 yr−1) recorded inside (red line) and outside (blue dashed line) of marine protected area (MPA) protection within the Danajon Bank, Philippines. Linear models plotted with biannual time since designation, for each site (left) and for pooled sites (right). 95% CI of models shown in grey

Full size image

Fish size distributions through time

The size of reef fish in the Danajon Bank has shifted towards proportionally higher abundances of smaller-bodied (0–10 cm) individuals through time (88% of total abundance) in unprotected reef areas (Fig. 3a; Appendix 4). Inside MPAs the size range is more evenly spread with 67% of fish being 5–25 cm long, but still right skewed. In the first year of designation, fish size kernal density distributions were not significantly different inside and outside of MPAs (p = 0.32), but became significantly more right skewed outside of MPAs in year 5 (p < 0.001) and year 10 (p < 0.001).

Fig. 3
figure 3

a Kernel density distribution estimates for fish size (Total Length) counts for periods 0, 5 and 10 years since marine protected area (MPA) designation, both inside (red) and outside (blue) of MPA management for eight sites surveyed along the Danajon Bank, Philippines (grey denotes overlap); and b linear regressions for (logged) fish count density (250 m−2) inside (red line) and outside (blue dashed line) of sites, against years since designation of the MPA site, within the same eight MPAs. Fish are separated into five 10 cm (total length) size bins and pooled for all size classes. The 40 cm + size class includes all size classes 40–100 cm. Note data points < 0 years since designation represent pre-designation monitoring

Full size image

Following separation of fish density data into five binned size classes (Fig. 3b), inside of MPAs the ≥ 20 to < 30 cm size class category showed significant decline (slope estimate = − 0.08, t value = − 2.46, CI = − 0.14, − 0.02) and the ≥ 30 to < 40 cm size class category showed a similar strong decline (slope estimate = − 0.08, t value = − 1.65, CI = − 0.14, − 0.01). However, outside of MPAs, the rate of change in fish density through time was shown to decline significantly faster than inside MPAs within the 0–10 cm (slope estimate = − 0.10, t value = − 4.16, CI = − 0.14, − 0.07), 10–20 cm (slope estimate = − 0.06, t value = − 2.61, CI = − 0.10, − 0.02), and 20–30 cm size classes (slope estimate = − 0.13, t value = − 1.76, CI = − 0.20, − 0.09), with the greatest overall rate of decline seen in the 20–30 cm size class fish outside of protection (Fig. 3b; Appendix 5). Fish density was significantly reduced outside of protection in the 20–30 (estimate = 6.42 fish/250 m2, t value = − 1.85, CI = 2.32, 18.17) and 30–40 cm size classes (estimate = 2.18 fish/250 m2, t value = −1.79, CI = 0.58, 8.33) and significantly elevated density found outside in the 0-10 cm size class (estimate = 295.89 fish/250 m2, t value = 3.89, CI = 138.38, 626.41). While no significant decline was observed for the > 40 + size class, the analysis of this class was limited by lack of data, (having to group five separate 10 cm size classes together), due to generally extremely low or zero abundance of large individuals across sites.

Relative fish abundances by management type

After 8–15 years of protection across sites, significant differences in relative mean abundances are seen inside versus outside of MPAs in Carangidae (jacks), Chaetodontidae (butterflyfish), Haemulidae (grunts), Holocentridae (soldierfish) and Serranidae (grouper), which were all more abundant inside managed areas (Fig. 4). All these families are large-bodied commercial valuable food fish, or in the case of the Chaetodontidae are valuable for the aquarium trade and sensitive to disturbance, being typically corallivores (Appendix 2). The only family with significantly higher abundance outside of MPAs was Zanclidae (Moorish idols).

Fig. 4
figure 4

The first sentence of the caption should read: Forest plots of the natural log-transformed response ratio (LnRR) of fish family-level mean abundances inside and outside of eight Philippine marine protected areas (MPAs), in 2011. Effect size from zero illustrated with squares, with lines illustrating 95% CI and significance indicated with (*). Mean and CI values also printed right

Full size image

Fish trophic composition through time

Following separation of the fish family biomass into trophic groups, declines through time were seen both inside and outside of MPAs for the piscivores, scavengers and invertivores (Fig. 5). Declines in fish biomass through time were most pronounced in the scavengers (i.e. the Haemulidae (grunts), Labridae (wrasse) and Mullidae (goatfish)), with significant declines both inside (slope estimate = − 0.10, t value = − 2.19, CI = − 0.20, − 0.01) and outside of MPA protection (slope estimate = − 0.18, t value = − 2.40, CI = − 0.23, − 0.09), but with significantly steeper declines outside of protection (Appendix 6). The piscivores (i.e. the Carangidae (jacks), Lethrinidae (emperors), Serranidae (groupers), etc.) similarly showed significant declines both inside (slope estimate = − 0.18, t value = − 4.07, CI = − 0.26, − 0.09) and outside of protection (slope estimate = − 0.21, t value = -0.83, CI = − 0.29, − 0.12), but with no difference in the rate of decline between management. Invertivores (i.e. the Gobiidae (gobies), Balistidae (triggerfish) and Monacanthidae (filefish) etc.) again saw significant declines in biomass inside (slope estimate = − 0.21, t value = − 0.83, CI = − 0.29, − 0.12) and outside of protection (slope estimate = − 0.31, t value = − 1.05, CI = − 0.48, − 0.14), but no significant difference in effect between management types. There was significantly reduced biomass outside of protection for corallivores (estimate = 0.41 kg/250 m2, − value = − 3.30, CI = 0.09, 1.88 kg/250 m2), herbivores (estimate = 2.64 kg/250 m2, t value = − 4.20, CI = 0.46, 15.37 kg/250 m2), piscivores (estimate = 7.91 kg/250 m2, t value = − 2.52, CI = 0.001, 0.002 kg/250 m2) and scavengers (estimate = 9.41 kg/250 m2, t value = − 1.93, CI = 2.12, 40.54 kg/250 m2). All other groups saw no significant declines or differences through time following MPA designation.

Fig. 5
figure 5

Log-transformed fish biomass (g/250 m2) inside (red line) and outside (blue dashed line) of marine protected areas (MPAs) through time, for eight family-level trophic groups in the Danjon Bank. Data pooled across all eight MPA locations and SE shown

Full size image

Effect of enforcement and time on benthos

Coral cover through time

Linear mixed-effects models showed no significant individual or interaction effects between live coral cover, management and length of time since designation for pooled sites, or for any other substrate category (Fig. 6; Appendix 7). Average percentage live coral cover from N = 472 surveys across 8 sites, did not significantly change over the study period inside of MPAs (slope estimate = − 0.04, t value = − 0.83, CI = − 0.12, 0.05), from a widely spread initial intercept estimate of 27.35% cover (t value = 10.28, CI = 17.89, 38.69% cover). Cover did not differ significantly outside of MPAs and similarly saw no decline through time (slope estimate = − 0.05, t value = − 1.13, CI = -0.12, 0.01), from an initial intercept estimate of 26.32% cover (t value = 9.87, CI = 13.54, 43.3% cover). Average coral cover did not differ significantly inside and outside of MPAs at the time of designation.

Fig. 6
figure 6

Percentage live coral cover recorded inside (red line) and outside (blue dashed line) of marine protected area (MPA) management within the Danajon Bank, Philippines. Percent cover recorded over 250 m2 transects. Linear models plotted with annual time since designation, for each site (left), and linear mixed-effects model outputs from all sites (right). 95% CI of models shown in grey

Full size image

The majority of the variance in coral cover (45% of the 46% explained) within the model was explained by site-level variation, but while there is an indication of declining cover from negative slopes for each management type, no significant trends were seen across the Danajon Bank during the observation period (Fig. 6).

On a site level, only Bilangbilangan and Pandanon changed significantly through time (Appendix 8). Within Bilangbilangan the management type caused significant differences in cover through time with cover outside MPAs declining (slope estimate = − 0.262, SE = 0.091, t1,440 = − 2.895, p = 0.004) and cover inside MPAs increasing (slope estimate = 0.134, SE = 0.057, t1,440 = 2.350, p = 0.019). For Pandanon a significant negative effect was observed inside MPAs (slope estimate = − 0.336, SE = 0.070, t1,440 = − 4.803, p < 0.001), while the slope outside the MPA was slightly (and non-significantly) less negative (slope estimate = 0.141, SE = 0.121, t1,440 = − 1.160, p = 0.247).

Benthic community composition

Separation of benthic groups into eight substrate classes allowed dissimilarities in community composition to be explored through time and between site and management factors (Appendix 7). Principal component analysis showed sample stations were clustered by site, showing distinct local variation in substrate composition, but showed no strong differences between management type (inside or outside) for any site (PCA1 = 32.17%, PCA2 = 18.20%, Appendix 9A).

PERMANOVA showed no interaction effect between time and management level (F18, 143 = 0.409, p = 1, R2 = 0.054) for the substrate groups and no effect of management (F1, 143 = 1.301, p = 0.223, R2 = 0.010). However, there was a weak interaction effect between time and site (F46, 143 = 1.795, p = 0.75, R2 = 0.182), with both time (F18, 143 = 1.524, p = 0.001, R2 = 0.095) and site (F7, 143 = 19.658, p < 0.001, R2 = 0.475) signalling a shift in community similarity towards sites with greater proportions of rubble, dead coral, sponge and sand.

A further PCA focussing on the coral growth morphological difference alone showed that the sample station coral growth forms were again clustered by site (PCA1 = 23.87%, PCA2 = 16.03%, Appendix 9B), but PERMANOVA showed no interaction effect between time and management level (F18, 143 = 0.553, p = 0.998, R2 = 0.074) for the coral trait composition, and no interaction effect between time and site (F46, 143 = 1.059, p = 0.324, R2 = 0.195). However, there was a strong effect of site (F7, 143 = 16.027, p < 0.001, R2 = 0.424) causing differences in coral community traits, and also a smaller effect of management type within each site (F8, 143 = 3.052, p < 0.001, R2 = 0.092).

Discussion

Broad-scale community shifts

The Danajon Bank reef system is exposed to intense localised pressures, primarily from heavy resource extraction, land-based pollution and the use of highly destructive fishing techniques, coupled with occasional large-scale natural disturbance events. Our study finds that these pressures have resulted in sustained declines in unmanaged fish populations. While small locally managed MPAs are maintaining fish biomass and a range of commercially exploited fish families within their boundaries, there are still significant reductions in total biomass and abundance of fish even inside of protection at individual sites, particularly of the piscivores/scavengers and mid-large (30–40 cm) fish size groups. There is typically very low or no large (> 40 cm) fish present anywhere throughout the reef. This study region appears to additionally still be recovering from the 1997–1998 large-scale coral bleaching event and frequent typhoons, with mean coral cover levels inside and outside of protection ~ 10% lower than country-wide averages and > 20% lower than pre-bleaching peak cover at the end of the study period (Magdaong et al. 2014).

Mean total fish biomass and abundance declined significantly in unprotected areas across the Danajon Bank region over the 12-year study period. Declines in fish density were most pronounced for larger fish (20–40 cm), both inside and outside of MPAs. These declines in commercially valuable size classes both inside and outside of protected areas could indicate potential infractions into the protected sites, such as at the Batasan MPA where the site guardhouse was noted to have been damaged, or at Pandanon where evidence of heavy blast fishing has been observed both inside and outside the MPA (Bayley et al. 2019). While small MPAs can be effective for the protection of mobile species (Claudet et al. 2010), the observed declines of larger fish in this study could also indicate an important limitation to these very small areas, which are typically less than the 2 km minimum diameter suggested for even partial protection of most reef fishes (Krueck et al. 2017). The declines indicate the need for larger reserves which fully capture species’ home ranges and for wider fisheries management in this area in order to fully protect stocks.

The strong positive response in abundance of large, slow-growing species, such as groupers, jacks and grunts inside the MPAs (along with butterflyfish and soldierfish), was surprising given the greater relative abundance outside of faster-growing (but also commercially exploited) species such as the parrotfish and rabbitfish. Home ranges of butterflyfish, soldierfish and grunts are typically very small (< 1 km), with groupers more variable in range, but with many species (e.g. Cephalopholis spp.) also showing high site fidelity, (Green et al. 2015). Rabbitfish and parrotfish meanwhile tend to typically be more wide ranging in their movement (> 2 km up to tens of kilometres) and so are less likely to benefit from the small MPAs studied here, which are typically ≤ 1 km in width (Green et al. 2015). Given the generally large home ranges of the jacks and parrotfish, and the tendency of some of these species to school/aggregate in large groups, it is possible that the biomass surveys of these family groups have been skewed/biased by their behaviour, as discussed by Richards et al. (2011). Parrotfish can also show a counter-intuitive increase in density in response to reduced live coral cover (Russ et al. 2015). Finally, the Philippines has a widespread trade in ‘live reef food fish’, with two minor hubs around Bohol (Mamauag 2004). This trade focuses on 22 fish species (including grouper, wrasse, parrotfish, snapper and rabbitfish), with the most targeted fish in this trade tending to be grouper (Mamauag 2004; Nañola et al. 2011), indicating the grouper family would be expected to respond well to protection.

Outside of protected areas, all size groups from 0 to 30 cm declined significantly over the period of study, with larger-bodied fish (20–30 cm) declining fastest of all size groups, indicating targeted removal of these more valuable larger fish (Mamauag 2004). The shift of the fish community size–frequency distribution outside of protection towards a greater proportional dominance of small-bodied fish could also potentially be a result of shifting fish community tropho-dynamics arising from reduced predation from the large piscivores/scavenger fish targeted throughout the reef (Wilson et al. 2008; Ritchie and Johnson 2009; Mumby et al. 2012).

Despite the observed changes in fish community structure, there was surprisingly no evidence of corresponding significant declines in live coral cover over this period. Likewise, there was no significant shift in the composition of the major substrate categories within each site, even between protected and unprotected areas, similar to observations by Graham et al. (2020) in the Seychelles over broadly the same period. Despite a weak (non-significant) shift towards communities with less branching, table, foliose and mushroom growth forms, no overall aggregate changes were observed in benthic form, cover or composition through time from the low initial mean coral cover (following coral bleaching mortality). Any differences in cover were instead largely driven through site-level heterogeneity. The fishing methods in this localised region therefore appear to typically not be permanently damaging to the benthos, unless intensive localised blast fishing is employed, such as at the site ‘Pandanon’ within this study (Selgrath et al. 2018; Bayley et al. 2019).

It is noted that the lack of any significant change in benthic composition could also be in part due to the survey method used missing species-level changes, because only broad substrate types and growth forms were recorded (Bertrand et al. 2006; Murphy and Jenkins 2010). If we are to assume the lack of any observed significant benthic community change is real, we can reasonably conclude that the size-selective reduction in fish populations would be largely due to focussed fishing pressure of medium to large-bodied commercially valuable fisheries, rather than due to any significant overall declines in benthic habitat availability (Pauly 1998; Robinson et al. 2017; Cinner et al. 2018). The declines through time in small-bodied fish could then be due to a trophic cascade effect brought about by the loss of large piscivores and scavengers from fishing/poaching (Casey et al. 2017; Robinson et al. 2017), resulting in meso-predator release and greater predation of low trophic groups (Ritchie and Johnson 2009; Hempson et al. 2018). Alongside this potential effect, it is important to note that this particular region experiences some of the most intense exploitation levels in the country, with significantly lower species richness for both target and non-target species (Nañola et al. 2011), indicating over-exploitation and some indiscriminate fishing methods still being used, such as blast fishing, and cyanide fishing for the live fish trade.

An alternative explanation, given no pre-bleaching coral cover baseline, for the observed shifts in fish size–frequency distribution in this region is due to a lag effect of the 1998 global bleaching episode, which saw Philippines-wide reductions in coral cover, including around Bohol (Magdaong et al. 2014). The subsequent loss of fish habitat and food before focussed monitoring began for this study (Wilson et al. 2006; Robinson et al. 2018) would likely have affected both the abundance and diversity of the fish populations (Syms and Jones 2000; Darling et al. 2017). While we do not have the appropriate baseline data available in this study preceding disturbance, this explanation matches well to other study observations and modelled predictions of fish biomass loss following habitat degradation (Graham et al. 2007; Rogers et al. 2018a, b).

The Danajon Bank fishery is therefore likely being impacted by both selective fishing pressure through time and by the long-term impacts of a severe bleaching mortality event. Nonetheless, the use of MPAs as a management tool is slowing this reef degradation effect overall. Despite MPAs dampening the decline in fish stocks, the current extent of overall biomass reduction and continued decline in the larger more fecund fish seen throughout the Danajon Bank in unmanaged areas are likely to increasingly hinder these fisheries’ future sustainability (Pauly et al. 2002; Hixon et al. 2014).

Management implications

Changes in coral cover through time at the site level were only significant at Bilangbilangan (seeing a significant increase in cover inside the MPA) and Pandanon (seeing a significant decrease in cover inside the MPA). These observations match well to other studies of this reef system which found Bilangbilangan to be relatively well managed, whereas Pandanon was seen to have regular systematic blast fishing occurring (even inside the MPA limits at points in the reserve’s history), causing loss of structure and a hampered recovery due to the resultant mobile rubble substrate (Samoilys et al. 2007; Bayley et al. 2019).

Despite some signs of infractions of the management areas, the CCEF management effectiveness scores around all sites are still rated as high (Alcala et al. 2008). The higher rate of decline in fish biomass around the central sites of the municipality of Getafe relative to other sites is likely to be in part due to their close proximity to the population hub of Cebu city and its ports, which is the second most populous area of the Philippines. Additionally, these sites’ geographical proximity and accessibility to other fishers from the nearby island villages of Cebu such as Mactan and Olango, as well as the highly populated areas of Bohol, are likely leading to encroachment and heavier resource extraction in these nearby waters. As previously noted, Nañola et al. (2011) found this central Visayas region to have the highest density of municipal fishers in the Philippines, and that this intense local fishing pressure has caused widespread species richness depletion, particularly in the Bohol area (Lavides et al. 2009; Nañola et al. 2011).

The reduction in overall fish stocks (both inside and outside of sites) matches well with larger studies showing that as population size and accessibility to reefs, termed ‘gravity’ of human impacts, increase, the biomass of fish seen in these areas, even inside MPAs, will typically be lower (Heenan et al. 2016; Cinner et al. 2018). It is also important to note that while current enforcement and management is good within these sites, initial and ongoing management effectiveness is highly dependent on technical or financial support from individual organisations or Municipal Local Government Units (MLGUs). Within the Bohol municipality of Getafe (which encompasses the four central sites, 3–6), the MLGU is typically lacking of funds, facilities and incentives, leading to elevated levels of destructive fishing in unmanaged areas around MPAs, as well as higher possibility of occasional incursions inside protected sites. Such capacity shortfalls are well documented to limit MPA effectiveness, and therefore, regional investment in such capacity is essential to improve management performance of these areas (Gill et al. 2017).

As many of the sites in this study have declining fish biomass even inside the MPA limits, losses may also be a result of edge effects from the MPA’s small coverage, lack of connectivity or an inadequate patch size to viably maintain the fish community (Joint Nature Conservation Committee 2010; Edgar et al. 2014; Green et al. 2015; Roberts et al. 2017). The MPAs in this study (and across much of the Philippines) are typically much smaller than many fishes’ home ranges or dispersal distance (Weeks et al. 2010; Green et al. 2015; Abesamis et al. 2017). While the cumulative coverage of designated MPAs (with varying levels of protection) intersecting with the Danajon Bank reef is now ~ 12.3%, the median MPA size is 0.25 km2 (UNEP-WCMC 2020). Increasing the number of these small sites, their proximity to one another and their size would likely help slow the overall decline of fish populations in this region by better protecting the more mobile and highly exploited large adult fish (Gaines et al. 2010; Abesamis et al. 2017; Krueck et al. 2017). These changes typically enhance the observed localised benefits to biomass and function and improve spillover benefit and connectivity (Abesamis and Russ 2005; Botsford et al. 2009; Lester et al. 2009). However, the enlargement of protection needed for the full benefit of these MPAs to be realised is likely in the order of 30% total coverage of the whole reef (Gaines et al. 2010), or an increase in individual MPA sizes to at least larger than the home range and mean dispersal distance of targeted species, to allow site self-replenishment (Botsford et al. 2001, 2009; Green et al. 2015). Any such increases in MPA coverage would of course require increasing management and enforcement capacity which is currently severely limited (with most MPAs typically voluntarily managed by low-income communities). These increases in size and capacity may be particularly effective in existing sites such as Pandanon, Jandayan Norte and Handumon, which are heavily exploited and closest to dense populations of people and fishery markets/ports (i.e. close to Cebu city and Northern Bohol) and have shown continuing rapid present-day declines in fish biomass or benthic cover (Cinner et al. 2018; Bayley et al. 2019). Without these improvements, it is unlikely any local benefits will be sustained by the small sites, and instead population gains will leak out or be encroached upon by the unmanaged and exploited surrounding region.