Unusual but consistent latitudinal patterns in macroalgal habitats and their invertebrate communities across two countries

The physical characteristics of biogenic habitats and environmental conditions are important determinants of biodiversity, yet their relative importance can change across spatial scales. We aimed to understand how relationships between the physical characteristics of macroalgal habitats and their invertebrate communities varied across spatial scales and whether general ecological patterns occurred across two countries.

The physical characteristics of habitats can significantly change localized environments, which may, in turn, alter associated community structure (Heck & Orth, 1980;Jones, Lawton, & Shachak, 1997). For example, at local scales, species abundances and diversity may correlate with a range of physical characteristics of habitats including patch size and shape (Airoldi, 2003;Andrén, 1994;Bruno & Kennedy, 2000), vegetation size (Kelaher, 2003a), density (Gribben & Wright, 2014;Heck & Orth, 1980;Nilsson, 1979), structural complexity (Matias, Underwood, Hochuli, & Coleman, 2010;Stein, Gerstner, & Kreft, 2014;Taniguchi, Nakano, & Tokeshi, 2003), frond length and surface area (Stelling-Wood, Gribben, & Poore, 2020) and multivariate structural gradients (Tews et al., 2004). In some cases, communities may vary due to the composition or diversity of habitats (Tews et al., 2004), whereas in others they may vary due to changes in the characteristics of individual biogenic habitats (e.g., density or frond length; Lawton, 1987) or both (Stelling-Wood et al., 2020). Interspecific variation in the physical characteristics of morphologically distinct biogenic habitats can have important consequences for associated communities (e.g., Angelini et al., 2011). Much less, however, is known about how variation in physical characteristics within individual habitat-forming species or complexes of morphologically similar biogenic habitats influence associated communities (but see Badano & Cavieres, 2006;Kelaher, 2002) even though variation in the physical characteristics of individual biogenic habitats may have similar strong consequences for associated communities. In support of this, at one site, Stelling-Wood et al. (2020) found that intraspecific variation in morphological traits was more important than species identity in predicting epifaunal abundances. Improving our understanding of the spatial relationships between the physical characteristics of biogenic habitats and associated communities is critical to developing effective conservation management strategies (Byers et al., 2006) and for predicting how biodiversity may respond to global environmental change (Wernberg, Thomsen, Tuya, & Kendrick, 2011). Some studies have addressed these relationships at local scales (e.g., Airoldi, 2003;Kelaher, Underwood, & Chapman, 2003;Palomo, People, Chapman, & Underwood, 2007); however, as biogenic habitats and their associated communities often co-occur across broad geographic ranges, there is a need to examine how these relationships change over large spatial scales, considering changes in the physical characteristics within (Fowler-Walker, Connell, & Gillanders, 2005a;Ralph, Morrison, & Addison, 1998;Rice, Kenchington, & Chapman, 1985) and among habitats (e.g., Messier, McGill, & Lechowicz, 2010). At present, we do not have a strong understanding of how these relationships vary across biogeographic gradients (but see Heck & Orth, 1980;Stein et al., 2014).
There is a range of environmental conditions-operating across multiple spatial scales-known to influence community structure. Latitude, and its common covariate temperature, is often related to the composition of communities (Cruz-Motta et al., 2010;Hillebrand, 2004;Schemske, Mittelbach, Cornell, Sobel, & Roy, 2009), as are moisture (Meynard et al., 2013), tidal height (Underwood & Chapman, 1998) and site exposure (Blanchette, 1997;Blanchette, Thornber, & Gaines, 2000). Anthropogenic factors such as pollution (Terlizzi, Scuderi, Fraschetti, & Anderson, 2005) and human population (Bloch & Klingbeil, 2016) have also been associated with Main conclusions: Our results support other studies showing that large-scale patterns can emerge from systems where there is high local-scale variability. The results show that communities within macroalgal habitats respond to both the physical characteristics of the habitat and external environmental conditions (e.g., temperature), suggesting that local-scale environmental factors, including anthropogenic stressors, may modulate environmental gradients over larger scales.

K E Y W O R D S
biodiversity, biogeography, community structure, habitat, invertebrates, latitudinal gradient, macroalgae, physical characteristics changes in community structure. Given that biogenic habitats can facilitate associated communities by dampening particular environmental stressors at small spatial scales (Bruno & Kennedy, 2000;Dijkstra, Boudreau, & Dionne, 2012;Wright & Gribben, 2017), then it could be hypothesized that the physical characteristics of habitats may be particularly important predictors of their associated communities as they can weaken gradients in external abiotic conditions at biogeographic scales (McAfee, Cole, & Bishop, 2016;Silliman et al., 2011). However, some abiotic stressors (e.g., pollution) may not be moderated by habitats and as such may either affect habitats and communities to similar degrees or may only affect either the community or the habitat.
Macroalgae create some of the most conspicuous and ecologically important habitats on rocky reefs in temperate regions around the world (Bellgrove, McKenzie, Cameron, & Pocklington, 2017;Schiel & Foster, 1986). To investigate the relationships between the physical characteristics of macroalgae and the invertebrate communities within them, we conducted biogeographic surveys of intertidal macroalgae on the temperate east coasts of Australia and New Zealand. Surveys were conducted across multiple spatial scales including between countries, sites within countries, and macroalgae within sites. We related observed spatial patterns in macroalgal physical characteristics and their invertebrate communities to 26 abiotic variables (natural and anthropogenic) sourced for each site. First, we tested the hypothesis that the physical characteristics (individual and patch level) of each macroalga were related to latitudinal gradients. Second, we hypothesized that invertebrate communities would vary in their structure and composition with latitudinal gradients. Third, we hypothesized that patterns in macroalgal physical characteristics and their invertebrate communities would relate to similar external environmental variables across the study area, but that the invertebrate communities would mostly correlate with changes in macroalgal physical characteristics. We did not sample across the total distribution of each macroalgal habitat but selected a study area where the habitats co-occurred across a similarly large latitudinal range in both countries. Biogeographic studies often use univariate measures of species richness or abundance to describe changes in community composition. In this study, we used multivariate analyses of invertebrate community structure and composition across a broad range of taxa, but with a lower taxonomic resolution as this approach can be effective at detecting large-scale diversity patterns even when there is high local-scale variability (Anderson, Connell, et al., 2005). F I G U R E 1 Study area including 18 sites along the east coast of Australia (10 sites) and New Zealand (8 sites). Symbols show the macroalgal habitats that were sampled at each site: Hormosira banksii (•), Coralline (▲), Sargassum (*) and Cystophora (×). Coordinate system: GCS_WGS_1984 2 | ME THODS

| Study sites
We sampled a single rock platform (25-75 m long) at each of 18 sites across the temperate east coasts of Australia and New Zealand ( Figure 1). In Australia, we surveyed 10 sites from Bonny Hills in Northern NSW to Eaglehawk Neck in Tasmania, ranging across more than 1,300 km (linear distance; Figure 1). In New Zealand, we surveyed eight sites from Leigh (northern New Zealand) to Shag Point, Otago (southern New Zealand) ranging across more than 1,000 km ( Figure 1). We selected Bonny Hills as the upper latitudinal limit of the study, as this coincides with the transition from temperate to sub-tropical climate, based on Köppen climate classes (Bureau of Meteorology, 2014). Within countries, sites were at least 10 km apart; however, sites were generally over 100 km linear distance from each other ( Figure 1). The latitudinal gradients sampled in both countries overlapped by 6.73 decimal degrees (Figure 1)

| Study organisms
We sampled three macroalgal species/complexes with distinctive physical characteristics in each country; two of which were shared between countries: Hormosira banksii (Turner) and red turfing algae (hereafter Coralline; see Figure S1). Hormosira banksii is a prostrate brown alga with beaded vesicles that are connected in chains 10-30 cm long. It is distributed in Tasmania and NSW in Australia and is widely distributed on both islands in New Zealand (Edgar, 2008).
Coralline included several morphologically similar species from the family Corallinaceae (e.g., Corallina officinalis, Amphiroa spp., Jania spp.). Species from the family Corallinaceae are widely distributed in temperate Australia and New Zealand (Atlas of living Australia website, 2018 at http://www.ala.org.au. Accessed 01 July 2018) and different species occurred interchangeably throughout the study area.
Grouping of Coralline species at the family level as a morphologically similar complex has also been done in other similar studies on habitatcommunity associations (Kelaher, 2002(Kelaher, , 2003a. The third macroalgal habitat sampled was Sargassum spp. (hereafter Sargassum) in Australia and Cystophora spp. (hereafter Cystophora) in New Zealand, which are closely related brown algae that occur in the lower intertidal zone (Edgar, 2008). These two habitats are both brown frondose, branching seaweeds, with receptacles either on branches in Sargassum or on vegetative fronds in Cystophora (Edgar, 2008). Cystophora was sampled at the genus level as multiple species occurred throughout the study area (e.g., Cystophora retroflexa, Cystophora scalari, and Cystophora torulosa) that would provide a broadly similar physical habitat structure compared to the other habitats, as fucoids with branching fronds.
For Sargassum, numerous, morphologically similar, species co-occur in Australia and accurate identification is difficult, being based on the seasonal size and shape of receptacles (reproductive structures at the end of the algal branches; Edgar, 2008); therefore, this habitat was categorized to genus. Sargassum are broadly distributed in Australia (though absent at some specific study sites, see below), and Cystophora are widely distributed in New Zealand (Atlas of living Australia website, 2018 at http://www.ala.org.au. Accessed 01 July 2018; Edgar, 2008). At each site, we sampled three macroalgal habitats except in: (a) Leigh and Picton in New Zealand where H. banksii F I G U R E 2 Survey design showing replicates at each scale of the study area. Each country (n = 2) included 8-10 sites sampled across the latitudinal gradient. Within each site, 2-3 macroalgal habitats were sampled with six replicate patches of each (see Figure 1 for the number of sites and habitats at each scale). Within each patch multiple habitat characteristics and invertebrate communities were measured. The spatial scales are shown at the base of the figure. Patterns were determined at large (latitude, country) and small (site, patch) spatial scales. All environmental covariates were sourced or collected at the site scale was absent, (b) Cook's Beach in New Zealand where Cystophora was absent, (c) Coles Bay in Australia where Coralline was absent, and (d) Eden and the two Tasmanian sites in Australia where Sargassum was absent ( Figure 1).

| Spatial patterns in macroalgal habitat characteristics
All macroalgal taxa were surveyed from January 2012 to early April 2012. Australian sites were sampled in a random order between January and April and New Zealand sites were sampled over a threeweek period in February. As ocean temperatures lag seasonally, the sampling period represented summer water temperatures. At each site, we sampled six replicate patches of each macroalgal habitat during low tide across the length of the rock platform ( Figure 2). The habitat patches selected occurred as discrete mono-specific patches with less than 10% of other habitat-forming organisms present.
Two patch-level characteristics were measured (patch area and percentage cover), plus two individual-level characteristics (frond length and biomass) of the macroalgae. Patch area was estimated by multiplying the longest and widest dimensions of each patch.
Frond length was determined from the average of 10 randomly selected fronds measured at the patch centre. Percentage cover of algae was approximated using a grid of regularly spaced points in a 25 × 25 cm quadrat. Macroalgal biomass was determined from two replicate core samples per patch. PVC cores (10 cm diameter) were collected near the centre of each patch, with algae scraped off at the rock surface with a paint scraper and placed into labelled plastic bags (Kelaher, Castilla, & Seed, 2004;Thrush et al., 2011). Biomass samples were rinsed over a 1 mm sieve to remove trapped sediment.
After excess water was drained, the algae were weighed in the field on digital scales (nearest 1 g). The two samples were then pooled to determine patch biomass. To ensure wet weight was an appropriate measure of biomass, samples of each macroalgal taxa were taken back to the laboratory and oven-dried at 60°C for 48 hr to determine dry weight (n = 12 cores/habitat). For each macroalga, wet and dry weights all were significantly correlated (Pearson's Correlation coefficient; r > 0.90).

| Spatial patterns in communities
Invertebrate communities in macroalgal patches were sampled using one of the PVC cores and collecting the invertebrates retained on the 1 mm sieve. To capture large or benthic invertebrates that may not be collected in the cores, a 25 × 25 cm quadrat with a 5 × 5 cm grid was used to survey larger, benthic invertebrates in each patch.
The survey was conducted by searching the fronds and substrate in each of the quadrat grid cells for macroinvertebrates (>2 cm).
All invertebrates from each patch (core + quadrat) were combined in a labelled plastic bag to capture one composite replicate patch ( Figure 2). The sample was later fixed in 7% formalin for a minimum of one week before being washed and transferred to 80% ethanol for preservation.
In the laboratory, all animals were identified and counted under a dissecting microscope. Molluscs were identified to family level and below (down to species), polychaetes to family level, crustaceans to order or suborder, echinoderms to class, Anthozoa to order, and foraminifera to phylum. The level of taxonomic identification related to the taxonomic group's dominance among samples and the condition of the samples required for fine-scale identification (e.g., although amphipods were a dominant group due to the high volume of collections and the time required to process the samples, some diagnostic features degraded after collection). It was also deemed more useful to include a large range of taxonomic groups identified to a coarse taxonomic level rather than a smaller range of taxa identified to species level in order to maximize chances of detecting habitat-community associations (Anderson, Connell, et al., 2005). Although our sampling methods may not capture all invertebrate taxa (e.g., barnacles, tube-dwelling polychaetes and colonial species such as sponges and bryozoans; Kelaher & Castilla, 2005), these did not appear common when sampling, most likely due to an absence of the bare rock they need for colonization (Edgar, 2008), and were excluded from the data set.

| Environmental variables
We sourced data on 40 environmental variables (13 of which were later excluded) related to coastal abiotic conditions (natural and anthropogenic) for each site from publicly available databases (Halpern et al., 2008), satellite images (e.g., NASA/NOAA; Meeus, 1991) and field observations (Table 1). All variables were sourced or determined at the site scale.

| Data analysis
All analyses were conducted using the PRIMER V7-PERMANOVA add on software (Anderson, Gorley, & Clarke, 2008;Clarke & Gorley, 2015). For all analyses on macroalgal habitat characteristics and environmental data, we used Euclidean distance matrices. Prior to the analyses, we normalized habitat characteristics and environmental variables as the data were measured in various units. For analyses on invertebrate communities, we used Bray-Curtis matrices. Community data were standardized by totals to reduce the influence of habitat patches with high biomass, and then square-root transformed to reduce the influence of a few highly abundant taxa on the results, while still preserving patterns of relative abundances among samples (Clarke, 1993;Clarke & Gorley, 2015). To visualize spatial patterns in habitat characteristics and communities at each site in each country, we produced principal coordinate analyses (PCO) plots. We did not evaluate patterns of spatial variation of environmental variables as they were used as predictor variables for the invertebrate and macroalgal habitat data.
We used RELATE analyses to test for linear latitudinal patterns in macroalgal habitats and communities, separately. RELATE is a nonparametric matrix correlation routine that measures how closely related two sets of multivariate data are by calculating a rank correlation coefficient (e.g., Spearman's) between paired samples in two similarity matrices (Clarke & Gorley, 2015). In this case, the biotic similarity matrices (habitats or communities) were compared to a "model matrix," which is an idealized distance matrix representing the structure of a hypothesized serial distribution representing a linear latitudinal gradient (Clarke & Gorley, 2015). This model determined Averages and standard deviations were estimated for the following: the full year during the day; the summer months (Dec-Feb) during the day a ; the full year during the night a ; and the summer months (Dec-Feb) during the night a (Meeus, 1991).
Chlorophyll-a (Chl-a; mg/m 3 ) Average and standard deviations of monthly values of the MODIS Aqua mission from July 2002 to December 2015. Averages and standard deviations were estimated for both the full year and for the summer months only (Dec-Feb; Meeus, 1991).
Chlorophyll-a anomalies a Number of events that surpassed 2 standard deviations of the average chlorophyll-a for a given year (Meeus, 1991).

Rainfall (mm)
Average and standard deviations of monthly accumulated rainfall from January 1979 through September 2015 obtained using TOVAS web-based application (Halpern et al., 2008).
Rainfall anomalies a Number of events that surpassed 2 standard deviations of the average rainfall for a given year.

Photoperiod a (index)
The difference between the sunrise and the sunset time, based on common astronomical formulae (Meeus, 1991 (Halpern et al., 2008).
Human population pressure (index) Estimated as the sum of total human population adjacent to the ocean within a 25 km radius. From LandScan 30 arc-second population data of 2005 (Halpern et al., 2008).

Ocean-based pollution (index)
Model based on combined commercial shipping traffic data and port data (Halpern et al., 2008).
Wave exposure Categorical variable with four levels indicating wave exposure at the site level. Categories were based on commonly used fetch measurements, with exposure defined by the openness of the site including the presence of offshore islands and protection provided by headlands (see Wernberg & Thomsen, 2005), though submerged reefs were not considered. a Variables that were excluded from the analysis: (i) the chlorophyll-a anomalies and (ii) rainfall anomalies indexes as they correlated with their respective standard deviation (>95%); (iii) photosynthetically active radiation (all variables-summer and annual averages and standard deviations) as they correlated with sea surface temperature annual day average (>95%); (iv) summer average and standard deviations for sea surface temperature during the day and night, as they correlated with the annual average of sea surface temperature day and night (>95%); (v) annual average and standard deviations of sea surface temperature during the night, which correlated with annual sea surface temperature during the day (95%); (vi) photoperiod as it correlated with annual average and standard deviations of sea surface temperature during the night.
whether sites located next to each other were more similar (in terms of the structure and composition of communities or habitat characteristics) than sites located further apart. Seriation models implicitly test for spatial autocorrelation in patterns of spatial distribution of the communities and habitat characteristics. As the macroalgal habitats had different physical characteristics (see Figure S1) and the communities had different taxonomic structure, we analysed each macroalgal habitat and community separately; however, our aim was to determine if the matrices for each habitat in both countries would correlate with the spatial models and with each other (i.e., habitats and their associated community). RELATE analyses compared the ranks of the habitat characteristic data matrix and the ranks of the community data matrix separately against each model. If patterns of spatial variation of the matrix were similar to the model, we would expect high and significant correlation values. RELATE analyses used Spearman rank correlations that tested for significance using permutations (999) of the original data to construct a null distribution.
Shade plots were used to visualize the relative abundances of taxa within each habitat across the latitudinal gradient in each country.
In each plot, taxa are clustered in a dendrogram based on an index of association (determined using a SIMPROF routine) with the dendrogram discriminating taxa that have correlated spatial distributions along the latitudinal gradient from those that do not (Clarke & Gorley, 2015).
For each macroalgal habitat, BEST analyses based on Spearman rank correlations were used to determine which environmental variables were most correlated with spatial patterns of habitat characteristics, or which variables and characteristics were most correlated with community structure. BEST is a nonparametric routine that finds the best match between two data matrices and searches over all possible subsets/combinations of variables to find the optimal model (Clarke & Gorley, 2015;Clarke & Warwick, 2001 of the community data are constrained by "optimal" combinations of the environmental data (Clarke & Gorley, 2015). Since environmental variables were estimated at the scale of sites and habitat and community were sampled at the patch level (different spatial scales), centroids per site were estimated for habitats and communities to run the BEST and LINKTREE correlation analyses.

| Spatial patterns in macroalgal habitat characteristics
Multivariate macroalgal habitat characteristics did not follow linear latitudinal gradients (i.e., the seriation model) for any of the habitats ( Figure 3; Table 2). However, the PCO plot suggested a nonlinear (C-shaped) gradient in the multivariate distribution in 5 of the 6 cases (except Coralline in Australia). To better describe the spatial pattern of macroalgal habitat characteristics (and test for the existence of a C-shaped pattern), RELATE was used with a cyclicity model, which tested whether sites at either end of the sampled distribution were generally more similar to each other than to adjacent sites, this showed a significant correlation for 5 or the 6 models (Table 2, Figure 3).

| Spatial patterns in communities
High abundances and diversity within invertebrate communities were observed throughout the study, with 10 s to 100 s of individuals collected in each habitat patch, and 77 taxonomic groups found across the study area. Taxon richness was not highly variable among habitats, with values per site ranging from 5 to 21 (mean 13, SD ± 5) in Hormosira banksii, 11-23 (mean 18, SD ± 4) in Coralline, 10-23 (mean 16, SD ± 5) in Sargassum, and from 6 to 18 (mean 11, SD ± 2) in Cystophora (see Figure S2). Communities in each habitat were generally distinct, except for H. banksii communities which had similarities with most other communities-except Coralline in Australia (see Figure S3). Abundances of invertebrates in H. banksii were the lowest and ranged from 62 to 628 per site (mean 224, SD ± 168), Coralline abundances were the highest and ranged from 86 to 8,048 per site (mean 1,486, SD ± 2,136), Sargassum ranged from 151 to 955 (mean 372, SD ± 269), and Cystophora from 412 to 3,040 (mean 1,475, SD ± 968) per site (see Figure S4). Though we expected differences in community structure between habitats because of their very different morphology, we were most interested in whether these distinct macroalgal habitats and communities would share similar spatial patterns and relationships with environmental drivers.
Consistent with analyses on macroalgal habitat characteristics, communities did not follow linear latitudinal gradients ( in Australia) appeared to have "C-shaped" latitudinal gradients in the PCO plots, similar to those described for the macroalgal habitat characteristics (Figure 3). Again, RELATE analysis confirmed that these communities were significantly correlated with the cyclical model ( Table 2), meaning that for these habitats' communities at the edges of the gradient were more similar to each other than to communities in the middle of the gradient. Although the remaining F I G U R E 3 Principal coordinate analysis (PCO) showing centroids of (a) habitat characteristics and (b) communities for each macroalgal habitat and site in Australia (AUS) and New Zealand (NZ). Vectors follow the latitudinal gradient to show the spatial relationships between sites. Symbols represent the three macroalgal habitats sampled in each country: Hormosira banksii (•), Coralline (▲), Sargassum (*) and Cystophora (×). Sites are numbered as shown in Figure 1  Note: The seriation model tested the prediction about linear latitudinal gradients, while the cyclicity model tested the significance of the "C-shaped" patterns (i.e., greatest similarities between southern and northern sites) that were apparent from the PCO results. The results show the Spearman rank correlations with the significance of the correlations (p values) shown in brackets. Significant results are highlighted in bold. Where there were two significant models for Hormosira banksii in NZ the highest correlation was considered to be the dominating pattern, accounting for some of the redundancy in the explicative models.

TA B L E 2 RELATE analyses showing
correlations of (a) habitat characteristics and (b) communities for macroalgal habitats sampled in Australia (AUS; n = 10 sites) and New Zealand (NZ; n = 8 sites)

habitats (Coralline in both countries and Cystophora in New Zealand)
were not correlated significantly with either a linear or cyclical distribution, they also appeared to have communities that were similar at the edges of the sampled distributions in the PCO plots, with sites at the north and south for both habitats generally clustering together, away from those in central sites (Figure 3). We do not present results from the third RELATE correlation as there were no significant correlations. Shade plots comparing relative abundances of taxa among sites and habitats showed that the dominant taxa in communities differed between macroalgal habitats. They also showed that some individual taxa exhibited spatial patterns across latitude that were consistent with the community patterns described above (see Figure S5) In New Zealand, Coralline communities had some taxa that were highly abundant in specific sites (e.g., the bivalve Kellia spp., and the gastropod family Eatoniellidae), as well as taxa that were common across sites (e.g., Nereididae and Amphipoda). The latitudinal patterns were not strong for these taxa, although they had smaller relative abundances at some central sites in Australia ( Figure S5c,d).
Sargassum communities also had some common taxa (e.g., Trochidae) across sites and some that were particularly abundant at specific sites (e.g., Amphipoda-which was also less abundant in the centre of the range; Figure S5e). Cystophora communities had large abundances of Eatoniellidae, Trochidae, Amphipoda and Isopoda. Again, the spatial patterns of individual taxa in Cystophora were not strong, but there were some large relative abundances of amphipods at the most southern and northern site ( Figure S5f).

| Correlations between communities, macroalgal habitat characteristics and environmental conditions
The communities, in contrast to macroalgal habitat characteristics, correlated significantly with environmental variables for all habitats in Australia and New Zealand, though the relationships between communities, macroalgal habitat characteristics and environmental conditions differed for each habitat (Table 3). Overall, human population pressure was an important component of the models and correlated with all communities. Sea surface temperature was also important for three of four communities. At least one habitat characteristic was included in the models for all communities. Biomass was included in the models for three of four communities, but other habitat characteristics were also important including frond length, canopy percentage cover and patch area. Some additional environmental or spatial variables were included in the models that were unique to individual communities, such as latitude, longitude, chlorophyll-a and ultraviolet light ( variables that were significant to multiple communities had contrasting spatial distribution patterns across latitude. Human population peaked towards the centre of the sampled range in Australia, though was relatively consistent in New Zealand ( Figure S6). Sea surface temperature, as is well known, decreased with increasing latitude towards the south of both countries. The habitat characteristic biomass was important for three of four communities, so we examined the univariate patterns in biomass across latitude for the relevant habitats. Coralline biomass peaked in the most northern site, though this appeared to be an outlier more than a trend and if removed there did not appear to be a latitudinal trend. Hormosira banksii biomass was smallest at the centre of the latitudinal range sampled. Sargassum biomass was also smallest at the centre of the sampled range-though it was sampled over the smallest range ( Figure 4).

The LINKTREE analyses enabled comparisons of invertebrate communities among macroalgal habitats in all sites in Australia and New
Zealand and indicated the relative importance of macroalgal habitat characteristics and external abiotic variables at different spatial scales.
These analyses highlighted the importance of sea surface temperature in explaining differences in communities between countries (split A in communities associated with the larger Cystophora were distinct from those in the physically lighter Coralline (split B, Figure 5), while H. banksii communities were divided between these two subgroups depending on macroalgal biomass (distinguishing sites on the north island from those on the south). Other macroalgal habitat characteristics were also important for further distinguishing communities in some sites for example per cent cover of algae related to Coralline and Sargassum communities in many Australian sites (split Q, Figure 5).

| D ISCUSS I ON
In this study, we investigated the role of environmental factors in shaping the physical characteristics of macroalgae and the invertebrate communities inhabiting them at continental scales. We found little support for our hypotheses that the physical characteristics of macroalgae and communities would follow latitudinal gradients.
Instead, patterns in both macroalgal habitat characteristics and associated communities generally shared unusual nonlinear patterns, where habitat characteristics and communities were more similar at the edges of the sampled distribution (i.e., north and south) compared to those at centrally located sites. Our models suggest that these patterns resulted from the influences of both small-and largescale environmental conditions and that communities within macroalgae were related to both macroalgal habitat characteristics and external environmental factors.
In our study, we wanted to test whether a suite of macrophyte characteristics would change across latitude, given that combinations of characteristics contribute to microhabitat environments (Airoldi, 2003;Heck & Orth, 1980;Kelaher, 2003aKelaher, , 2003bMatias et al., 2010;Nilsson, 1979;Taniguchi et al., 2003). We found unusual nonlinear gradients in suites of macroalgal habitat characteristics across latitude for five of the six macroalgae examined. The consistent pattern was observed across the macroalgal habitats, despite each of the taxa sampled having very different morphologies (e.g., turfing red algae versus. brown frondose algae; see Figure S1) and Australian macroalgae (Coralline and H. banksii) being slightly larger on average than the same taxa sampled in New Zealand.
Because the nonlinear "C-shaped" pattern in macroalgal habitat characteristics occurred in 5 of 6 tests, we expected there would be a clear large-scale factor driving the result; however, this was not the case. In fact, three of the models showed no significant variables and the others showed a combination of either localized environmental conditions (e.g., pollution and human population pressure) or largescale conditions (e.g., sea surface temperature). As the habitat characteristics did not correlate consistently with environmental variables, the nonlinear pattern is hard to reconcile. One explanation is that at the centre of the temperate ranges sampled, abiotic conditions related to latitude (e.g., temperature) may be generally favourable and macrophytes more strongly respond to local conditions (e.g., Maggi, Milazzo, Graziano, Chemello, & Benedetti-Cecchi, 2015). However, at the edge of the distributions sampled, large-scale environmental conditions may be the dominant driver of macrophyte characteristics. The result of this may mean that localized abiotic drivers at the centre of the range may mask macroecological patterns (e.g., Kerswell, 2006), which may only exert influence over communities at their edges. In fact, the peak hypothesis of range edges predicts that intraspecific traits change from the centre to the periphery of a species distribution due to systematic changes in biogeographic factors from the species core habitat (Gaston, Chown, & Evans, 2008;Lloyd, Murray, & Gribben, 2012). However, we did not cover the total distributions of each macroalgae in this study so we cannot generalize about their global distribution pattern.
The latitudinal diversity gradient is one of the most widely studied ecological patterns (Hillebrand, 2004;Kraft et al., 2011).
Similar to the habitat patterns, communities also generally displayed "C-shaped" biodiversity patterns. Interestingly, the nonlinear biodiversity pattern occurred despite each of the macroalgal habitats in each of the countries housing distinct communities (see Figure S3).
Examination of abundance patterns for individual taxa across latitude (see Figure S5) revealed that few taxa showed strong C-shaped patterns across latitude, indicating that the result was due primarily to the changing abundance and composition of the entire community. The distinction between communities in different macroalgal habitats and countries was expected because of their very different physical structures (see habitat-specific examples in: Edgar & Klumpp, 2003;Gemelli, Johnson, & Wright, 2019;Kelaher, 2003b).
For example, Coralline had large relative abundances of polychaetes, which was expected due to its dense turf-like structure (see Figure   S5c,d). Although the communities differed between habitats, the majority of the taxa sampled (75%) were shared between the two countries. Macroalgal characteristics can be an important driver of associated community structure and the similarities between macroalgal habitat and community patterns suggest that communities were responding to variation in the physical habitat characteristics.
Only H. banksii in New Zealand showed a traditional latitudinal gradient in community variation. The reason for this difference is not clear, though these communities had stronger associations with climatic variables including ultraviolet light and chlorophyll-a. Communities may have responded more strongly to these large-scale external abiotic variables as H. banksii occurs slightly higher on the shore than Cystophora and has less dense protective canopy cover than both Cystophora and Coralline.
When we investigated the environmental conditions that may be driving these patterns, we could not identify a consistent factor that correlated with all macroalgal habitats and communities. Although we predicted that communities would mostly correlate with habitat characteristics, we found that both habitat characteristics and external environmental variables were important. Of the characteristics that were significant, habitat biomass was the most common characteristic explaining differences in communities between sites within countries (see also Stelling-Wood et al., 2020). Sargassum also had a nonlinear pattern with biomass the smallest at the centre of the study area and largest at the north and south ( Figure 4). Biomass was not a significant variable for Cystophora communities. Canopy percentage cover was also an important factor structuring communities and other characteristics (e.g., patch size and frond length) also appeared within the habitat-specific models, suggesting the unique physical multivariate environment may promote distinctive communities by providing specific microhabitat conditions for the species they house (Angelini et al., 2011). The dominant factors such as habitat biomass and canopy percentage cover were expected to be an important factor in community structure as they are a proxy for habitat availability such as structural complexity and surface area (Kovalenko et al., 2012). In the exposed and highly variable system of rocky shores, space is a limiting factor and so the provision of structure for colonization is important to diversity (Matias et al., 2010;McGuinness & Underwood, 1986) as is their role in mediating environmental stressors (Dayton, 1972;Jones et al., 1994;Jurgens & Gaylord, 2017).
In the external abiotic environment, human population pressure was the only variable that correlated with all communities. Large human populations are known to have negative effects on rocky intertidal communities (Bloch & Klingbeil, 2016;Terlizzi et al., 2005;Thompson, Crowe, & Hawkins, 2002) and can cause reductions in canopy forming algae cover (Benedetti-Cecchi et al., 2001).
Abundance patterns of some of the dominant taxa showed that F I G U R E 4 Macroalgal habitat biomass (g) per site (mean ± SD) across the latitudinal gradient for (a) Hormosira banksii, (b) Coralline, (c) Sargassum in Australia (black circles) and New Zealand (white circles). Cystophora biomass is not shown as it was not a significant characteristic for those communities. Latitude is negative as sites are in the Southern Hemisphere. Sea surface temperature was the most important correlate with communities in the LINKTREE analysis at the country scale, with temperature strongly related to differences between the countries and also between mainland Australia and the island of Tasmania suggesting there may be a spatial hierarchy in its influence. Both large-and small-scale temperature gradients are frequently recognized as having effects on rocky shore biota (e.g., Cruz-Motta et al., 2010;Harley et al., 2012;Mabin, Gribben, Fischer, & Wright, 2013). Temperature influences community structure in both marine and terrestrial environments (Brown, 2014;Gaston, 2000). The differences in the invertebrate taxa sampled in each of the countries may have also influenced this distinction-though the majority of invertebrates sampled were common to both countries. The taxonomic resolution of invertebrate communities in our study may have also influenced the spatial patterns we identified, as coarser taxonomic identification (i.e., above species level) can have a higher likelihood of detecting large-scale F I G U R E 5 LINKTREE analysis showing the environmental and macroalgal habitat characteristics that correlated with spatial patterns in communities across countries, sites and habitats. Hormosira banksii = •, Coralline = ▲, Sargassum = *,and Cystophora = ×; Australian sites are the filled symbols and New Zealand sites are open. Sites are numbered as shown in Figure 1. The y-axis (A%) is the equi-spaced representation of the average between groups ranked Bray-Curtis dissimilarities, representing the magnitude of differences between the subset of samples. Continuous lines represent statistically significant splits, dotted lines represent non-significant splits or groupings. Letters on each split represent subgroups of the environmental or habitat variables that correlated with each split. Split A shows sea surface temperature (SST) differences across the study area (R = .48, B = 83%) left split < right split. Split B = habitat biomass within the lower SST group (R = .54, B = 61%) left split < right split. Split P = habitat biomass within the higher SST group (R = 0.55, B = 70%) left split > right split. Split Q = % cover of macroalgal habitat (R = 0.80, B = 76%) left split < right split. The remaining splits were not significant or were related to differences between individual sites diversity patterns (Anderson, Connell, et al., 2005). Nevertheless, we detected multi-scale variation in habitat characteristics and communities, indicating the resolution was not too coarse to detect a range of variation. Few studies have investigated the response of communities to large-scale differences in multivariate habitat characteristics both within and between biogenic habitats. In considering the importance of multivariate habitat characteristics to communities, the specific habitat characteristics that influenced communities differed between macroalgal habitats, but that across all habitats their characteristics generally played an important role in shaping communities. We also found that a combination of characteristics contributed to community structure patterns and so using a multivariate approach is important for observing associations at large scales. As the environmental conditions that related to community composition were idiosyncratic, our results highlight the importance of observing multiple scales of variation in biogeographic research (Fraschetti, Terlizzi, & Benedetti-Cecchi, 2005;Hewitt, Thrush, Dayton, & Bonsdorff, 2007). These results are important, as a lack of generality in biogeographic patterns is frequently cited as a reason for discounting the need for macro-  (Bearham, Vanderklift, & Gunson, 2013;Fowler-Walker et al., 2005a;Wernberg, Coleman, Fairhead, Miller, & Thomsen, 2003).
Although this research identified novel patterns in macroalgal habitat characteristics and associated communities, further research is needed to test the generality of these patterns, particularly in relation to the nonlinear (C-shaped) distribution patterns across latitude.
Understanding the relative contributions of large-and small-scale factors to community structure and how they interact with habitat characteristics to determine biodiversity patterns also needs more targeted research (Bruno et al., 2003;Bulleri et al., 2012;Maggi et al., 2015). Although we surveyed broad temperature gradients over 1,000 km in each country, we did not cover the entire distributions of the macroalgal taxa (this would be challenging as some of the genera sampled have global distributions This study demonstrated that multi-scale factors are important in determining patterns in community structure, including those operating across latitudinal and local scales. We also showed that it is important to consider variation in biogenic habitat characteristics when aiming to understand relationships between habitats and communities across large spatial scales (Gaston et al., 2008). In fact, our data suggest that habitat characteristics and the communities within biogenic habitats respond similarly to differing but potentially stressful environmental conditions. Though much research has focussed on identifying the most important spatial scale for driving biodiversity, this study and other research in this area (e.g., Kerswell, 2006;Meynard et al., 2013) show that a more nuanced approach is needed that considers that biodiversity patterns at local scales can be nested within a hierarchy of regional and global-scale processes (Gaston, 2000;Ricklefs, 2004).
Understanding these driving factors is important given the changes occurring in intertidal environments due to increasing anthropogenic stressors.

ACK N OWLED G EM ENTS
The authors would like to thank the Sydney Institute of Marine

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in Dryad at: https://doi.org/10.5061/dryad.51c59 zw53.