The role of herbivores in shaping subtropical coral communities in warming oceans

Tropicalization is rapidly restructuring subtropical marine communities. A key driver for tropicalization is changes in herbivory pressure that are linked with degrading ecosystem stability. Consequently, subtropical algal beds are being displaced by climate-mediated colonisation of coral communities. This process is thought to be aided by the elevated herbivory resulting from tropicalization, but the relative contribution to herbivory by different taxa is not fully understood. Evaluating herbivory pressure and its effect on coral cover and rugosity across a subtropical latitudinal gradient will help predict how these processes may change with further tropicalization and ocean warming. Herbivory pressure exerted by fishes and urchins across this subtropical latitudinal gradient remains unquantified. Using in-situ feeding observations, we quantify fish and urchin herbivory pressure at seven sites across non-accreting coral communities, and warmer accreting coral reefs in southern Japan. We then relate herbivory pressure to respective fish and urchin community structure and coral cover and rugosity. Urchin herbivory is greater on non-accreting coral communities than on true coral accreting reefs; a result which is reversed for fish herbivory. Overall, herbivory pressure is greater on accreting coral reefs than on coral non-accreting communities, but is dependent on reef characteristics as community structures differ more strongly among reefs than between regions. These factors are linked to coral cover and rugosity that differ between reefs, but not between climatic regions, further emphasising the influence of local factors on the benthic cover and the associated fish and urchin community, and thus herbivory pressure. Our findings provide a foundation for understanding how non-accreting coral communities may respond to ongoing tropicalization, given the fish and invertebrate herbivores they host.


Introduction
As sea surface temperatures rise and marine heatwaves increase in frequency and intensity (IPCC 2019; Robinson et al. 2019;Misra et al. 2021), shallow marine ecosystems are being placed under stress (Dalton et al. 2020;Graba-Landry et al. 2020). Worldwide ocean warming is now facilitating novel phase shifts in coastal marine systems, whereby tropical species are extending their poleward range into temperate regions through tropicalization (Wernberg et al. 2013;Vergés et al. 2014a, b;Kumagai et al. 2018;Donelson et al. 2019). Range contractions are slower than range expansion, resulting in dynamic transition zones where range-shifted tropical species overlap in geographic distribution with native communities (Vergés et al. 2014a;Donelson et al. 2019). Biogeographical transition zones are characterised by entirely new species assemblages and altered interspecific interactions (Vergés et al. 2014a;Donelson et al. 2019). These climate-driven changes in ecological communities can lead to altered biotic interactions, which have profound consequences for ecosystem functions within novel phase shifts (Vergés et al. 2014a;Pages et al. 2018;Kingsbury et al. 2020). For example, algal forests provide crucial ecosystem services such as habitat provision, coastal defences, primary and secondary productivity and act as nurseries for a range of marine species (Steneck et al. 2002;Smale et al. 2013). New biotic interactions as tropical herbivores expand their range into temperate regions have already resulted in largescale losses of ecologically productive algal forests (Ling et al. 2009;Vergés et al. 2014aVergés et al. , 2016. Deforested algae are subsequently replaced by non-accreting coral communities, yet the indirect trade-offs for community composition and successive biotic interactions from this phase-shift are not well understood (Vergés et al. 2014a(Vergés et al. , 2019Donelson et al. 2019). Therefore, a comprehensive understanding of species interactions within transition zones is needed to determine how these interactions may change under tropicalization.
Herbivory is a key interspecific interaction which underpins ecosystem stability (Bakker et al. 2016). Herbivory also reinforces ecosystem resilience and allows reef ecosystems to maintain algae-free states (Hughes et al. 2007;Bennett et al. 2015). Alterations to herbivory pressure can cause drastic changes in community structure (Bianchi et al. 2014;Ling et al. 2018;Vergés et al. 2016). For instance, a release from herbivory pressure can lead to a significant increase in algal cover (Steneck et al. 2017), whilst elevated herbivory pressure can cause complete deforestation of macrophyte-rich ecosystems (Wernberg et al. 2013;Vergés et al. 2014bVergés et al. , 2016Kumagai et al. 2018).
Environmental perturbations accompanied by changes in herbivory pressure can lead to irreversible ecosystem phase shifts (Vergés et al. 2014a;Dakos et al. 2015). For example, the 2010 marine heatwave in Western Australia destroyed extensive kelp forest areas and many temperate benthic invertebrates were replaced by subtropical and tropical species, resulting in a drastic ecosystem phaseshift (Wernberg et al. 2016). Similarly, the climate-driven influx of tropical herbivores due to ocean warming in southern Japan has deforested some kelp ecosystems, followed by the establishment of zooxanthellate coral communities (Nagai et al. 2011;Vergés et al. 2014aVergés et al. , 2016Wernberg et al. 2016). In the context of ocean warming, these and future changes in herbivore distribution and relative abundance will play a crucial part in shaping shallow marine communities. Understanding the mechanisms behind herbivore-mediated phase shifts in tropicalising reefs is vital to determine the impacts of future community composition and ecosystem functions (Vergés et al. 2014a(Vergés et al. , 2016. Uncertainties lie around how changes in marine herbivory pressure will interact with other climate-mediated stressors to facilitate algal deforestation and the tropicalization of temperate communities (Vergés et al. 2016;Donelson et al. 2019).
In conservation planning, ignoring biological interactions such as herbivory can lead to the development of inadequate management strategies (Kavousi 2019). Understanding which areas may exhibit changes in herbivory pressure could provide a basis for more effective conservation. For example, foundation species such as zooxanthellate corals rely on high herbivory levels (Vergés et al. 2014a). Ignoring herbivory pressure when identifying potential coral refugia can compromise their effectiveness (Vergés et al. 2014a;Kavousi 2019). Similarly, regulating fishing practices which alter herbivory pressure through the removal of primary consumers can reduce the stress on already threatened coral reefs (Edwards et al. 2014). To understand how herbivory pressure may change with ocean warming, current levels of herbivory need to be quantified across the tropics and subtropics.
Herbivorous fishes regulate algal cover, enabling zooxanthellate coral recruitment and growth via competitive release (Jouffray et al. 2015). Different functional groups of fish herbivores target different benthic types (Jouffray et al. 2015). Consequently, a high biomass and diversity of herbivorous reef fishes is associated with high coral reef rugosity and cover in the tropics (Chong-Seng et al. 2012;Graham and Nash 2013;Russ et al. 2021). Fish herbivory pressure on tropical reefs changes based on feeding preferences, that in turn change with benthos type (Hamilton et al. 2014). Further, feeding rates can serve as a measure of ecosystem function (Lefcheck et al. 2019; but see Gouhier and Pillai, 2019 for issues associated with this work). Similarly, phase shifts between coral and algae-dominated communities have been largely associated with changes in fish herbivory pressure (Molina-Hernandez et al. 2018;Martinez-Rendis et al. 2020). Although species richness and abundance of fish herbivores decrease with latitude (Floeter et al. 2005;Steneck et al. 2017), there has been almost no quantification of the herbivory pressure across a subtropical latitudinal gradient (Bennett and Bellwood 2011;Longo et al. 2014Longo et al. , 2019. Despite being a critical component in shallow marine ecosystems, invertebrate herbivores are relatively not well understood compared to herbivorous fishes (Francis et al. 2019;Mudge et al. 2019). Sea urchins are keystone herbivores in shallow marine ecosystems (Steneck 2020) and can significantly reduce algal cover and density (Westbrook et al. 2015;Ishikawa et al. 2016). Changes in sea urchin density have been associated with phase shifts in tropical (Hughes et al. 1987;Idjadi et al. 2010) and temperate regions (Ishikawa et al. 2016). However, whether there are differences in urchin herbivory pressure within the transition zones between true accreting coral reefs and non-accreting coral communities that grow on hard substrates in higher latitudes (Buddemeier and Smith 1999;Beger et al. 2014) remains unclear. Although urchins and fish herbivore assemblages are closely associated with present algae communities (Easton et al. 2018), the collective and cumulative herbivory pressure of both taxa found in transition zones has not yet been assessed. This knowledge gap hampers our understanding of the relative herbivory pressure exerted by both fishes and urchins and their relationship with coral reef and coral community characteristics.
Here, we quantify the relationship between herbivory pressure and reef characteristics on true, accreting coral reefs and non-accreting coral communities in higher latitudes of Japan. We first compare fish and urchin community composition across accreting true coral reefs in Okinawa and non-accreting coral communities in Kochi. Next, we quantify and compare the herbivory pressure exerted by fishes and Diadema urchins across sites. Finally, we assess the relationship between herbivory pressure and accreting coral reef and non-accreting coral community characteristics, and the relative influence of fish and Diadema herbivory in such interactions. To our knowledge, this is the first study to quantify and contrast herbivory pressure exerted by fishes and sea urchins across accreting coral reefs and non-accreting coral communities. Furthermore, we use a holistic approach to understand the relationship between fish and urchin community structure, their relative herbivory pressure, and coral reef characteristics.

Site description
Our fieldwork took place between June and July 2019, in two latitudinal regions in Japan; higher-latitude Kochi and lower latitude Okinawa. Okinawa's coral reefs are true accreting reefs (Hibino and van Woesik 2000), whereas with ongoing tropicalization, some reefs within higher-latitude Kochi Prefecture have transformed from kelp-dominated to supporting non-accreting coral communities growing on hard substrates (Nakamura et al. 2013;Vergés et al. 2014a). Both the main Japanese islands and Okinawa are influenced by the Kuroshio Current running northwards from the tropics, facilitating the tropicalization of surrounding reefs (Kumagai et al. 2018). Three non-accreting coral communities (Nishidomari Sheltered, Nishidomari Exposed, and Komame) in Kochi (> 8 m), and four warmer, accreting shallow reef sites (< 12 m) (Miyagi Channel, Kin Bay, Hamasango Point, and Hamatsu) in Okinawa were selected for the study ( Fig. 1; Appendix S1, Table 1). Moreover, we classified the sites into "exposed" and "sheltered" based on the coral reef or coral community exposure to waves. Exposed sites include Nishidomari, Miyagi Channel, Kin Bay and Hamatsu. Sheltered sites include another Nishidomari location, Komame and Hamasango Point.

Community surveys
We recorded abundance of all fishes and urchins via SCUBA visual census following the Reef Life Survey Methods Manual (https:// reefl ifesu rvey. com/ wp-conte nt/ uploa ds/ 2019/ 02/ NEW-Metho ds-Manual_ 150815. pdf), with modifications. We conducted three replicate 10 m video belt transects at each site to film reef fish within ~ 5 m either side of the transect and cryptic fish within ~ 2 m on either side of the transect. The same GoPro Hero 2 was used to film the fish to ensure the same area size was captured in each transect. We used the video footage to count and identify all fishes to the lowest taxonomic resolution using Coral Reef Fishes (Lieske and Myers 2002) and Japanese Saltwater Fish Photo Search (Yoshino and Seno 2008) guidebooks. We used 10 m transects (shorter than the standard visual census method) as we were unable to identify coral patches exceeding 10 m in length at most non-accreting coral sites. We used 10 m transects across all sites to ensure standardised sampling, and we generated species accumulation curves to account for missed species in our results interpretation (Appendix S1, Fig. 1, Table 2). We assigned each fish to a trophic group based on their feeding habits (Appendix S1, Table 3) inferred from FishBase (Froese and Pauly 2000) and field observations (Beger, unpublished data). We photographed every urchin within 1 m distance of the transect, identifying individuals to species level using SeaLifeBase (Palomares and Pauly 2020). As the reef sites at Kochi and Okinawa are not as structurally complex as typical tropical coral reefs, the crevices in which Diadema urchins rest during the day were easily accessible, allowing for urchin identification.
We captured the benthic cover by placing and filming a 50 × 50 cm 2 quadrat every 1 m along each 10 m transect. We identified all zooxanthellate corals to genus level using Coral Finder 3.0 (Kelley 2016) and Corals of the World (Veron et al. 2016). Morphological diversity of corals provides the basis for the structural complexity of the habitat (Denis et al. 2017) and is a key driver of coral rugosity (Darling et al. 2017). We recorded coral growth form (branched, columnar, digitate, encrusting, foliose, solitary, laminar and massive) and used its growth diversity (Shannon's H) as a proxy for coral rugosity where greater diversity of growth forms corresponded to greater rugosity. Using Coral Point Count with Excel (CPCe) (Kohler and Gill 2006), we extracted benthic cover data by assigning 30 points per quadrat to benthic cover categories. The benthic cover categories assigned were: epilithic algal matrix (EAM), macrophytes, coralline algae, Alcyonacea (soft coral), zooxanthellate corals, other cnidarians, Echinodermata, Porifera, other invertebrates, substrata, and the tape, wand, and shadow categories as required by CPCe. We filmed all video footage using GoPro Hero 2 and GoPro Hero 5 cameras and took photographs using an Olympus Tough TG5 with underwater housing.

Quantifying fish herbivory
To quantify fish herbivory, we considered all herbivorous fishes, as well as omnivorous and planktivorous fishes that were observed feeding on algae. We initially aimed to estimate fish bite rates per 1 m 2 of area. Initially, we set up a remote GoPro camera facing a 1 m 2 quadrat placed haphazardly over a reef substrate. However, we found that the quadrat altered fish behaviour, as fish actively avoided feeding within the quadrat and instead appeared to choose to feed around the quadrat boundaries. To overcome this (Sheltered); N2 = Nishidomari (Exposed); K = Komame (Sheltered); MC = Miyagi Channel (Exposed), KB = Kin Bay (Exposed); HP = Hamasango Point (Sheltered); and H = Hamatsu (Exposed) (Appendix S1; Table 1). Shapefile curtesy of https:// www. igism ap. com challenge, feeding fishes were instead filmed by a diver for seven minutes at a time as another diver conducted the community surveys. After the seven minutes of filming, the diver approached a different group of feeding fish and filmed it for seven minutes. We conducted all feeding observations within proximity of the transect but not on the transect itself to avoid interfering with the visual census of the fishes. Across three separate data collection dives on each site, a total of six to eleven videos of feeding fishes were filmed per site. The number of videos filmed per site was dependent on the number of feeding groups of fishes present within the sample site. The seven minute videos consisted of a two minute acclimation period during which no data on feeding rates was collected. We removed from the analyses any fishes which remained in the video frame for less than five seconds, and those species which appeared less than three times. We identified and observed all fishes for the entirety of each video. We recorded the time spent in the frame by each fish and counted the number of bites taken by each fish. We defined a bite as a clear contact between the fish's mouth and the EAM substrate or macroalgae. We divided the number of bites taken by each fish by the time of observation (seconds) and multiplied it by 30 to yield a bite rate r in bites/30 s. As the bite rates of the same species did not vary between sites (Appendix S1, Table 5), we calculated an average bite rate for each species across all the sites (Appendix S1, Table 6). For each transect, we estimated the fish biomass m of each species by combining the fish counts and the relative length-weight ratios available from FishBase (Froese and Pauly 2000; Appendix S1, Table 4). In the case where length-weight data for a species was not available, the average length-weight ratio across the genus was used. As fish sizes cannot be determined accurately using GroPro video footage, and local fish size information is limited, we used maximum length (cm) for all the observed fish species. We then multiplied the average bite rate r for each species by their relative biomass m on each 10 m belt transect. This calculation yielded the mass standardised herbivory pressure P for each species with three repeats per site (one per 10 m belt transect). Using fish biomass rather than abundance to quantify herbivory pressure is important as fish size affects its carbon intake with each bite (Hoey and Bellwood 2009;Lefcheck et al. 2019) however, we acknowledge that using maximum fish lengths likely overestimates herbivory across all sites. We calculated the total herbivory pressure exerted by fish P from all algae-feeding species per transect as:

Quantifying urchin herbivory
We estimated urchin feeding rates using Diadema as the model genus. Diadema urchins are important reef grazers (Nozawa et al. 2020) and are more sensitive to increases in water temperature than other echinoids (Sherman 2015), making them more susceptible to ocean warming. Furthermore, as a non-burrowing urchin, Diadema can be moved with minimal force in comparison to other abundant echinoids. We estimated herbivory rates by collecting a randomly selected Diadema from the vicinity of the transect but not from the transect itself and placed it in a 50 × 50 × 50 cm mesh cage. The mesh size was approx. 1 cm 2 , allowing water flow but preventing other herbivores from entering the cage. We cleared the substrate off any algae and secured the cages using pins or string to prevent them from being uplifted. We weighed down the edges of the cages with rocks to prevent the urchins from escaping. We chose Sargassum as the test food macroalgae as it is widespread across Pacific Japanese coasts (Haraguchi and Sekida 2008), it has been sucessfully used as a food source in Diadema laboratory studies (Coppard and Campbell 2007), and it is consumed by Diadema on shallow, Japanese reefs (Abe et al. 2008;Fujita et al. 2013). We collected Sargassum fronds from the littoral zone in Kochi and Okinawa before each dive. Laboratory studies report Diadema grazing rates on Sargassum to be < 30 g/24 h (Coppard and Campbell 2007). To avoid food limitation, we used 50 g of wet Sargassum in each cage with a single urchin. To weigh out 50 g of Sargassum, we blotted it dry and weighed out 50 g until its mass no longer changed due to water evaporation or dripping. After 24 h, we removed the algae from each cage and reweighed it. We used the change in mass to calculate the grazing rate in grams of wet algae ingested per 24 h (Wg/24 h) for each urchin. We completed at least two experimental periods on each site and moved the cages to a different location within the site between experiments. Depending on the Diadema presence on the site, we used between two and four test cages for each experimental period. In addition to the test cages, we used one urchinfree control cage for each experimental period to account for any changes in Sargassum mass not associated with urchin herbivory. The control cages also contained 50 g of wet Sargassum, which was reweighed after 24 h. Any changes in the Sargassum masses in the control cage were used to adjust the final Sargassum mass values for the test cages within the same experimental period.

Comparing the combined fish and urchin herbivory pressure across sites
We compared the overall herbivory pressure exerted by fishes and Diadema on each site by converting the herbivory pressure into comparable units of grams of dry algae ingested per 24 h (Dg/24 h) using the equations from Van Rooij et al. (1998) and Sissini et al. (2017) (Appendix S1, Table 6). For each independent transect, we multiplied the daily species-specific dry algae intake by the abundance of each species. We summed the dry algae intake for each species, yielding three repeats for each group per site.

Data analyses
We conducted all data analyses in R Studio v.3.6.1 (RStudio Team 2019). We carried out a non-metric multidimensional scaling (NMDS) analysis with Bray-Curtis distances (Faith et al. 1987), a square root transformation and a Wisconsin's double standardization to assess the fish and urchin community structure across all the sites (Gardener 2014). We fitted the benthic cover categories onto the ordination as environmental vectors. We added the abundances of fishes and urchins onto the ordination as additional vectors and assessed the significance of all vectors using 999 permutations. We used an Analysis of Similarities (ANOSIM) with Bray-Curtis distances and 999 permutations, to determine whether there is a significant difference in species composition between accreting coral reefs and non-accreting coral communities, between the sites, and between exposed and sheltered sites (Gardener 2014) (Appendix S1, Table 1).
We conducted an Analysis of Variance (ANOVA) to compare the coral cover and coral rugosity across the sites and regions (Dytham 2011). We tested the data for normality of residuals, homogeneity of variance and homoscedasticity of residuals through a visual inspection of diagnostic plots and Shapiro-Wilk tests (Beckerman et al. 2017). We applied a square root transformation to the coral cover data to meet the assumptions of ANOVA.
To compare the number of bites taken by different fish species, we used Quasi-Poisson Generalized Linear Model (GLM) with sites as a fixed factor and a log offset for the length of fish observation. We chose this model as the variance of the data amongst the groups was heterogeneous and a Poisson model exhibited large overdispersion (Wedderburn 1974). We applied a square root transformation to the total fish herbivory pressure to better match the normal distribution. We compared the total fish herbivory pressure between sites and regions using Welch's ANOVA (Moder 2010) with a Games-Howell Post-Hoc (Ruxton and Beauchamp 2008) as the assumption of equal variances was violated. Using ANOVA, we compared the urchin grazing rates. We applied a Poisson GLM and a Quasi-Poisson GLM to assess the differences in Diadema herbivory pressure between the sites, and between the accreting coral reefs and non-accreting coral communities, respectively (Wedderburn 1974).
After a square root transformation, we compared the combined urchin and fish algae intake across sites and regions using Welch's ANOVA to account for heterogeneity of variances (Moder 2010). We conducted a Games-Howell Post-Hoc for pairwise comparisons of combined algae intake across the sites (Ruxton and Beauchamp 2008). We tested for a relationship between taxa-specific herbivory pressure, combined algae intake, coral cover and rugosity using Spearman's Rank correlation (Gardener 2014).

Community structure
Across the seven sites, we recorded 69 species of fish and five species of urchins. During the feeding observations, an additional 23 species of fish were observed outside of the transects. These fish species were not included in the community structure analysis to maintain equal area sampling across sites. The optimal NMDS (with a stress level of 0.13) revealed clear clustering of the accreting coral reefs and nonaccreting communities (Fig. 2). NMDS ordination showed community structure of accreting coral reef sites was associated with greater soft coral cover and scleractinian coral rugosity, but not scleractinian coral cover. This was further apparent with ANOSIM with Bray-Curtis Dissimilarity which revealed significant differences in fish and urchin species composition between accreting coral reefs and nonaccreting coral communities (r 2 = 0.368, P = 0.002), between the sites (r 2 = 0.766, P = 0.001), and between exposed and sheltered sites (r 2 = 0.190, P = 0.026; Fig. 2, Table 1). Fourteen fish and two urchin species significantly contributed to the differences in species composition across the sites (Appendix S1, Table 8). Thirteen fish and four urchin species were found in both latitudinal regions (Appendix S1, Table 9). EAM, coral cover, coral rugosity and soft coral cover were significantly correlated with the fish and urchin community structure (Appendix S1, Table 10).

Fish herbivory
We collected a total of 667 feeding observations of 25 species. We included 14 fish species and one Scarus sp. that were observed feeding on algae and were present on the transects in the herbivory analyses. Species which were not detected on the transects would have a total biomass of zero, leading to the herbivory pressure exerted by this species also equalling zero. This means our transects did not identify at least 11 algae consuming fish species across all sites. We observed Acanthurus maculiceps and Siganus unimaculatus during the site surveys but not during the feeding observations. We obtained the mean bite rates of S. unimaculatus from the literature (Nanami 2018). There is no published literature on the bite rates of A. maculiceps and thus for the purpose of this study, we inferred the mean bite rates from other observed acanthurids. The mean number of bites for each species is listed in Appendix S1, Table 6. Herbivory pressure (√P) varied significantly between sites (F = 10.503, P = 0.007; Fig. 4a) with Kin Bay (Okinawa) experiencing fish herbivory pressure (130.0 ± 24.7) greater than Komame (23.7 ± 22.1, P = 0.04) and Nishidomari (Exposed) (Kochi) (6.46 ± 5.8, P = 0.04). The herbivory pressure exerted by fish was significantly greater on the accreting coral reefs (164.0 ± 100) than on the non-accreting coral communities (38.2 ± 46.6; F = 14.703, P = 0.001; Fig. 4b; Appendix S1, Table 12).

Urchin herbivory
We conducted a total of 29 urchin grazing experiments and 9 controls across all non-accreting coral community sites and the accreting coral reef site of Kin Bay. There were no Diadema available for the experiment at the rest of accreting coral reef sites. There was no difference in Diadema grazing rates across the sites (F (2, 25) = 1.643, P = 0.205; Appendix S1, Fig. 2). This allowed for the Diadema abundance to be used as proxy for herbivory pressure. Komame had Fig. 2 A three-dimensional ordination of the fish and urchin community structure from all 21 transects across the seven studied sites: non-accreting coral communities in Kochi: N1 = Nishidomari (Sheltered), N2 = Nishidomari (Exposed), K = Komame; accreting coral reefs in Okinawa: MC = Miyagi Channel, KB = Kin Bay, HP = Hamasango Point, H = Hamatsu. Community structure differs significantly between accreting coral reefs (orange triangles) and non-accreting coral communities (blue circles) (r 2 = 0.368, P = 0.002). Community structure also differs between the individual sites (r 2 = 0.766, P = 0.001). The arrows show the direction in which the benthic cover vectors drive the distribution of the sites within the ordination (Appendix S1, a significantly greater Diadema abundance (and therefore herbivory pressure) than Nishidomari (Sheltered) (Z = 3.506, P < 0.001), Nishidomari (Exposed) (Z = − 4.110, P < 0.001), Kin Bay (Z = 3.506, P = 0.005), or Hamatsu (Z = 3.506, P = 0.005; Fig. 4c; Appendix S1, Table 12). There was no Diadema present at Miyagi Channel and Hamasango Point, hence urchin herbivory experiments could not be conducted. Herbivory pressure exerted by Diadema was significantly greater on the non-accreting coral communities than on accreting coral reefs ( Fig. 4d; Appendix S1, Table 12; T = − 2.356, P = 0.029).

Combined urchin and fish herbivory
The combined urchin and fish herbivory calculated as the total daily algae intake (g/24 h) differed significantly between sites (F = 9.428, P = 0.008). A Games-Howell Post-Hoc pairwise comparison of sites revealed that only Fig. 3 Differences in coral cover where A coral cover varies significantly between sites (F (6, 14) = 19.46, P < 0.001; Appendix S1, the herbivores at Nishidomari (Sheltered) and Nishidomari (Exposed) differed significantly in the total daily algae intake ( Fig. 5a; Appendix S1, Table 12; T = 7.420, P = 0.020). Comparing the total daily algae intake between the two climatic regions revealed that herbivory was significantly greater on accreting coral reefs than on the non-accreting coral communities ( Fig. 5b; Appendix S1, Table 12; F = 6.651, P = 0.018).

Relationship between herbivory and coral reef characteristics
There was no relationship between fish herbivory pressure and coral cover on either the accreting coral reefs (P = 0.130) or the non-accreting coral communities (P = 0.380; Fig. 6a). Similarly, we detected no relationship between fish herbivory pressure and coral rugosity on accreting coral reefs (P = 0.150) or on the non-accreting coral communities (P = 0.680; Fig. 6b). There was no relationship between Diadema abundance and coral cover on the non-accreting coral communities (P = 0.170) or accreting coral reefs (P = 0.420; Fig. 6c); nor between Diadema abundance and coral rugosity in the non-accreting coral communities (P = 0.850) or the accreting coral reefs (P = 0.690; Fig. 6d; Appendix S1, Table 13). There was also no relationship between total urchin abundance (all urchin species identified) and coral cover on the accreting coral reefs (P = 0.420) and the non-accreting coral communities (P = 0.170; Fig. 6e). Total urchin abundance was weakly negatively correlated with coral rugosity on the accreting coral reefs (r s = − 0.61, P = 0.037) but not the non-accreting coral communities (P = 0.220; Fig. 6f). Similarly, combined daily algae intake showed no correlation with coral cover on the accreting coral reefs (P = 0.670), or on the non-accreting coral communities (P = 0.230; Fig. 6g). However, coral rugosity was positively correlated with daily algae intake on the nonaccreting coral communities (r 2 = 0.78, P = 0.013). There was no relationship between the daily algae intake and coral rugosity on the accreting coral reefs (P = 0.200; Fig. 6h).

Discussion
Our study provides the first urchin and fish herbivory pressure quantification across accreting coral reefs and higherlatitude non-accreting coral communities, and the complex relationship between herbivory pressure and reef characteristics. Overall, the herbivory pressure exerted by fishes was significantly higher on accreting coral reefs than on the non-accreting coral communities, whilst Diadema urchins exerted higher pressure on the non-accreting coral communities (Fig. 4). However, this result is likely driven by the higher abundance of Diadema at Komame, compared to at all other sites. This suggests site-specific factors influence herbivory pressure exhibited by Diadema urchins more than latitudinal location. Due to the logistics associated with quantifying herbivory in situ, only Diadema urchins were considered in the urchin herbivory analyses. Other urchins observed across the sites may also exert herbivory pressure, but due to their behaviour, were not used in this study. For example, Echinometra mathaei is known to feed predominantly on algae, but transplantation into the cages utilized in our experiments was not possible as it is a burrowing urchin (Khamala 1971). Consequently, we likely underestimated the herbivory pressure exerted by urchins. The use of Diadema urchins also likely underestimated urchin herbivory on exposed sites as Diadema urchins are thought to dominate sheltered reefs, whilst burrowing urchins, such as Echinometra, dominate exposed reefs (Bronstein and Loya 2014).
Whilst some evidence suggests that urchin feeding rates increase with temperature (Pages et al. 2018) we found no evidence for latitudinal effects on Diadema herbivory (Appendix S1, Fig. 2). When accounting for all species, overall urchin abundance was not correlated with coral cover on accreting coral reefs or on non-accreting coral communities (Fig. 6e). However, there was a weak negative correlation between urchin abundance and coral rugosity on accreting coral reefs (Fig. 6f). This is interesting as although urchins can be important herbivores, they are also important bioeroders (Alvarado et al. 2016); herbivory and bioerosion co-occur when herbivores remove algae from the reef and consequently erode calcium carbonate (CaCO 3 ) (Peyrot-Clausade et al. 2000;Perry et al. 2012) from calcified structures. Bioerosion rates, like herbivory, largely underpin coral reef stability (Barkley et al. 2015). Bioerosion is critical for healthy reef function as it plays a role in calcium carbonate cycling, allowing particles from dead coral to be reused during bioaccretion (Tribollet and Golubic 2011). However, excessive bioerosion can lead to a switch from net calcium carbonate accretion to net erosion (i.e. negative carbonate budget) (Alvarado et al. 2016). As such, increased bioerosion can be related to lowered coral rugosity. It is possible that on accreting coral reefs, which already exhibit high herbivory pressure, urchins are relatively redundant herbivores and predominantly result in greater bioerosion of those reefs and reduce the coral rugosity (Alvarado et al. 2016). Together, the combined herbivory pressure exerted by fishes and urchins was significantly greater on accreting Fig. 4 A Herbivory pressure (P) exerted by fishes varies significantly across the sites (F = 10.503, P = 0.007). B P is significantly greater for accreting coral reefs in Okinawa (blue) than on the non-accreting coral communities in Kochi (orange) (38.2 ± SD 46.6, F = 14.703, P = 0.001). C Diadema abundance (and thus herbivory pressure) is greater at Komame than any other site (Z = 3.506, P < 0.001). D Diadema abundance is significantly greater on the non-accreting coral communities in Kochi than on accreting coral reefs in Okinawa (T = − 2.356, P = 0.029) coral reefs that on non-accreting coral communities (Fig. 5). Additionally, we highlight the importance of considering urchins in herbivory analyses, as we observed that the relationships between herbivory pressure and coral reef characteristics vary significantly between the latitudinal regions. There was a positive link between the overall combined herbivory pressure and coral rugosity on the non-accreting coral communities, but not on the accreting coral reefs, however this relationship was not detected when fish and urchins were considered separately (Fig. 6). This result was further reflected in the differences of the fish and urchin communities that are present in each region (Fig. 2).
There was a clear difference in urchin and fish community structure between accreting coral reefs and non-accreting coral communities (Fig. 2). However, tropical reef species such as Chaetodon lunulatus, C. melannotus, Chromis margaritifer, Dascyllus trimaculatus, Neoglyphidodon melas, Pomacentrus coelestis and Thalassoma lutescens were also observed on higher-latitude non-accreting coral communities. The presence of tropical species in non-accreting coral communities contributes to the existing evidence of tropicalization of this area of the Japanese Archipelago (Nakamura et al. 2013). Furthermore, with ongoing tropicalization and in particular the expansion of tropical species into higher latitudes, the fish communities hosted by non-accreting coral communities in Kochi may in the future resemble those currently hosted by accreting coral reefs in Okinawa. Such community shifts on non-accreting coral communities could modify the patterns of herbivory currently observed. As tropical herbivores expand into higher latitudes, they are likely to increase the herbivory pressure in their extended ranges (Vergés et al. 2014a). It is therefore possible that with ongoing tropicalization, herbivory pressure will increase on non-accreting coral communities in Kochi, however, this will likely be dependent on local factors affecting the community structure.
Fish and urchin community structure differed more strongly between sites than between latitudinal regions ( Table 1), suggesting that factors other than latitude play an important role in which fish and urchin species are present locally. Here, we found that exposure had a weak effect on the community structure of fish and urchins, but we also found significant differences in herbivory pressure between two sites that mostly differed in wave exposure (Nishidomari "Sheltered" and "Exposed"). This result suggests local micro-habitats such as exposure (amongst other factors) Fig. 5 A Nishidomari (Sheltered) exhibited significantly greater algae consumption (g/24) by fishes and Diadema than Nishidomari (Exposed) (Games-Howell Post-Hoc; T = 7.420, P = 0.02). There was no significant difference in algae consumption between other sites. B Total daily algae consumption is significantly greater on accreting coral reefs in Okinawa (orange) than on non-accreting coral communities in Kochi (blue) (Welch's ANOVA; F = 6.651, P = 0.018) affect the community structure, which in turn affects the strength of herbivory pressure. Exposure to waves is known to influence community structure (Easton et al. 2018;Taylor et al. 2018). Life history traits further affect the structure and dynamics of coral associated communities (De Roos et al. 2003), which could explain the influence of reef exposure on the fish and urchin community structure. As exposure was determined to have only a weak effect on the community structure (Table 1), it is likely that a combination of factors including latitude, exposure and other site-specific Fig. 6 The detected relationships differed between accreting coral reefs in Okinawa (orange) and non-accreting coral communities in Kochi (blue). Spearman's rank correlation detected a negative relationship between urchin abundance and coral rugosity in Okinawa (r 2 = − 0.61; P = 0.037); and a positive relationship between coral rugosity (H) and total algae consumption in Kochi (r 2 = 0.78; P = 0.013) features drive the composition of fishes and urchins. For example, the fish and urchin community structure were strikingly different at Komame than those at any other site. This was also reflected by the high abundance of Diadema urchins at this site. Unlike the other sites, Komame is heavily influenced by freshwater run-off, which may affect the reef characteristics as well as fish and urchin community structures. Although reef proximity to freshwater discharge is associated with lower coral diversity in the tropics (West and Van Woesik 2001), the non-accreting coral community at Komame exhibited the greatest coral cover of all sites (Fig. 3a, Table 2). However, this site was largely homogenous, dominated by large foliose Pavona corals with some massive Porites colonies. This low diversity of corals and their growth forms is reflected by the site's relatively low coral rugosity (Fig. 3c, Table 2), and is likely driving the observed simplified community structure of the fishes and urchins (Gratwicke and Speight 2005).
It is apparent that coral cover and EAM cover are associated with opposing fish and urchin communities but not on a subtropical latitudinal gradient (Fig. 2). EAM is the low algal turf combined with detritus which creates the primary grazing surface for reef fishes (Wilson et al. 2003). Algal turf consists of early successional stages of macroalgae and small filamentous algal species (Graham and Nash 2013), and in the absence of herbivory, it can completely overgrow corals, leading to a phase-shift from a coral-to an algaedominated state (Jouffray et al. 2015). However, in this study, coral cover was not related to the herbivory pressure exerted by fishes, Diadema, or the combined pressure of both groups (Fig. 6). This is interesting, as typically a high biomass and diversity of herbivorous reef fishes is associated with high coral cover in the tropics (Chong-Seng et al. 2012;Graham and Nash 2013). Moreover, coral rugosity, which is more consistently associated with fish species richness and biomass than coral cover (Gratwicke and Speight 2005;Hall and Kingsford 2021), was not associated with fish or Diadema herbivory pressure.
These findings pose several questions. Firstly, are nonaccreting coral communities and warmer subtropical accreting coral reefs influenced by herbivory in the same way as tropical coral reefs, or is the biomass of herbivores in those regions not high enough to influence coral cover or rugosity? To answer this question, a comprehensive assessment of herbivory pressure across a larger latitudinal gradient using the same methodology is required. For example, Bennet and Bellwood (2011) observed a poleward decline in herbivory pressure on the tropical Great Barrier Reef (GBR), but their study includes different taxa and uses different methods than our work.
Secondly, as tropicalization continues, will herbivore biomass become sufficient to shape subtropical non-accreting coral communities and accreting coral reefs? Or do herbivory pressure and herbivore biomass not necessarily relate to coral reef characteristics in the same way despite being closely associated? For example, high coral rugosity provides a structurally complex habitat for reef organisms, providing prey species with refuge, but also may restrict predator detection through physical obstacles (Graham and Nash 2013;Gauff et al. 2018). If a reef exhibits high coral rugosity but insufficient algal cover (low food availability), the trade-off between predation risk and food reward may be too high, leading to low herbivory pressure at a site with high rugosity (Holbrook and Schmitt 1988). On the other hand, if coral rugosity is high, but food quality and availability are also high, we may expect to observe greater herbivory. The matter is further complicated by the regional differences in species behaviour. For example, the regional difference in fish herbivory pressure on the GBR is not a result of differences in fish biomass but changed fish behaviour: the fish herbivores were less likely to consume transplanted Sargassum in the southern region of the GBR (Bennet and Bellwood 2011). We did not detect differences in fish feeding behaviour between our studied regions, as the biting rate of the same species did not vary between the sites (Appendix S1, Table 5). These contrasting observations further highlight the need for a comprehensive assessment of herbivory pressure across a larger latitudinal gradient, using the same methodology.
Urchins are evidently important herbivores in nonaccreting coral communities (Ishikawa et al. 2016;Steneck 2020). Herbivory exerted by Diadema urchins was significantly greater in the non-accreting coral communities than on Okinawan accreting coral reefs ( Fig. 4d; Appendix S1, Table 12). This pattern was reversed for fish herbivory ( Fig. 4b; Appendix S1, Table 12). When combined, overall herbivory pressure estimated as daily algae intake was significantly greater on accreting coral reefs ( Fig. 5b; Appendix S1, Table 12). This finding suggests that the majority of the herbivory pressure is driven by fishes, but as combined herbivory pressure only includes herbivory estimates for Diadema urchins, we are likely underestimating the overall herbivory pressure on non-accreting coral communities. Nevertheless, our findings support the paradigm that lower latitudes experience greater herbivory pressure (measured as herbivorous fish richness and biomass) than coral habitats at higher latitudes (Floeter et al. 2005;Steneck et al. 2017).
As coral cover did not differ significantly between the accreting coral reefs and non-accreting coral communities, but did differ between sites (Fig. 3a, b; Appendix S1, Table 10), it is likely that local factors influence coral cover more than latitude. On the other hand, coral rugosity was significantly greater on accreting coral reefs than on the non-accreting coral communities (Fig. 3d). Furthermore, combined herbivory pressure was strongly associated with coral rugosity on non-accreting coral communities, but not on accreting coral reefs (Fig. 6h). This relationship may be due to new herbivorous species increasing the herbivory pressure on non-accreting coral communities, reducing the competitive pressure on coral by maintaining low algal cover, enabling coral recruitment and consequently leading to greater coral rugosity and more suitable habitat for reef fish in a positive-feedback effect (Denis et al. 2013;Jouffray et al. 2015;Steneck et al. 2017).

Conclusions
There are clear signatures of herbivory pressure influencing coral reef characteristics, however this relationship is highly complex and to accurately disentangle it we must consider local factors influencing community structure and thus, we must better understand species interactions. With continuing ocean warming (IPCC 2019), further tropicalization of non-accreting coral communities is likely. This tropicalization may result in an elevated herbivory pressure on nonaccreting coral communities, associated with the increase in fish herbivory. It is therefore possible that non-accreting coral communities may experience further algal loss, with the potential recruitment of new coral communities. However, as it is evident that other local factors affect the fish and urchin community structures on reefs, the effect of increased herbivory pressure will vary between sites. Moreover, as urchins are also present on accreting coral reefs, it is difficult to speculate whether they will remain at similar densities under future conditions on non-accreting coral communities. If upon tropicalization of the subtropics urchin density remains unaltered, it is possible that coral growth may be hindered by high levels of bioerosion. Bioerosion is therefore another biological process which requires further research along a subtropical environmental gradient.
Our study provides a better understanding of herbivory pressure across a subtropical latitudinal gradient, providing a basis for future research and understanding of how herbivory may change in the future. As herbivory pressure underpins ecosystem stability (Barkley et al. 2015;Bakker et al. 2016), and thus influences the likelihood of phase shifts, herbivory and local factors influencing it should be considered in future conservation planning. Firstly, subtropical areas of conservation concern which currently exhibit signs of tropicalization can be assessed on their vulnerability to phase shifts underpinned by elevated herbivory pressure. This assessment can inform management practices aiming to minimise the risk of a phase-shift. On the other hand, areas of currently lower ecological importance which are experiencing tropicalization could be considered as potential coral refugia. Herbivory pressure should therefore be considered in future coral reef management. Overall, the findings of this study can enhance identification of sites which may be vulnerable to phase shifts associated with changes in herbivory, and further the current knowledge of the relationship between herbivory pressure and coral reef characteristics.
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:// creat iveco mmons. org/ licen ses/ by/4. 0/.