Fishery catch is affected by geographic expansion, fishing down food webs and climate change in Aotearoa, New Zealand

Introduction

Shifts in marine species’ geographic distribution with ocean warming are increasingly well documented (Cheung, Watson & Pauly, 2013; Pinsky, Selden & Kitchel, 2020; Chaudhary et al., 2021). Further shifts may alter community structure (García Molinos et al., 2016; Poloczanska et al., 2016), food web interactions (Tekwa, Watson & Pinsky, 2022), and fisheries catch (Cheung et al., 2010; Cheung, Watson & Pauly, 2013). Ocean warming can also impact the fitness of marine ectotherms, particularly large, active predatory species (Forster, Hirst & Atkinson, 2012; van Rijn et al., 2017; Lavin et al., 2022). Concurrently, historical fishing pressure has resulted in deleterious impacts on marine communities and ecosystems (Howarth et al., 2014). Prolonged, intensive fishing results in reductions of population biomass (Worm et al., 2009), as well as the maximum size of harvested species, as fisheries preferentially target the largest species and individuals within a population (Jackson et al., 2001) due to their higher market price (Tsikliras & Polymeros, 2014). Once the largest species, often of the highest trophic levels and long-lived, have been depleted, fishing pressure has been observed to shift to smaller, shorter-lived species of lower trophic levels—termed ‘fishing down marine food webs’ (Pauly et al., 1998; Pauly & Palomares, 2005). The occurrence of fishing down can have deleterious consequences to the transfer of energy through marine food webs, thus impacting marine biodiversity and ecosystem function though the depletion of high trophic level predators (Pauly et al., 1998; Baum & Worm, 2009; Boyce et al., 2015). Thus, both fishing impacts and climate-change related outcomes to marine populations should be considered in concert in order to inform fisheries and conservation management (Lynam et al., 2017).

One method to infer changes in the geographic distribution of marine fish and invertebrate populations includes the Mean Temperature of the Catch index (MTC, Cheung, Watson & Pauly, 2013). This index represents the catch-weighted mean temperature preference of species in the fisheries catch in an area, and has been applied to a variety of marine ecosystems globally (see Leitão et al., 2018 and Dimarchopoulou et al., 2021b and references therein). Overall, an increase in MTC over time indicates an ocean-warming induced increase in the proportion of species of warmer thermal affinity relatively to those of colder thermal affinity recorded in the catch (Cheung, Watson & Pauly, 2013; Tsikliras & Stergiou, 2014; Liang, Xian & Pauly, 2018). Similarly, one method to infer the occurrence of fishing down includes the Mean Trophic Level index (MTL, also known as the Marine Trophic Index) (Pauly et al., 1998; Pauly & Watson, 2005), i.e., the catch-weighted mean trophic level of fisheries within an area. Also utilized extensively (see Liang & Pauly, 2017, and references therein), a decrease in MTL indicates the transition of catch from large, long-lived, higher trophic level species to smaller, short-lived, lower trophic level species, as a result of overexploitation and thus fishing down (Pauly et al., 1998). Conversely, increasing trends in MTL indicate the increasing catch of higher trophic level species, as has been observed during the historical expansion of fisheries (Leitão, 2015).

Yet, for both of these well-utilized indices, apparent trends in catch-weighted values can be masked by a variety of factors. One masking factor includes the ‘skipper effect,’ in that skippers will continuously and preferentially target certain species, usually of high trophic levels, within a multispecies catch (Pinnegar et al., 2002; Liang & Pauly, 2020). This may also arise when catch quotas prevent shifting catch from one species to another. Another masking effect includes the taxonomic resolution at which analyses are performed, as trends in catch-weighted values may not be apparent when analyzed at coarse taxonomic resolutions (Pauly & Palomares, 2005; Liang & Pauly, 2017). Lastly, when such indices are calculated across large geographic areas (e.g., exclusive economic zones, EEZs) and time scales, the gradual geographic and depth expansion of fisheries to include new habitats and species (i.e., deep-water species of colder thermal affinity, or higher trophic level species) can mask trends of fishing down or MTC warming (Bhathal & Pauly, 2008; Kleisner, Mansour & Pauly, 2014).

Fisheries in Aotearoa (henceforth referred to as New Zealand) waters have experienced a significant transformation within the 20^th century through technological proliferation and government incentives. Following the establishment of their EEZ in 1978, New Zealand became the first country to implement a Quota Management System (QMS) in 1986, in recognition of declining stocks of important coastal fisheries (Johnson & Haworth, 2004). Since then, New Zealand has prioritized the transition towards ecosystem-based fisheries management in line with the targets outlined in the Convention on Biological Diversity (Cryer, Mace & Sullivan, 2016; Durante, Beentjes & Wing, 2020). Overall, the transition towards an ecosystem-based, multi-species management scheme requires in-depth understanding of the historical effects of fisheries (Pauly, 1995), as well as those effects presented by trends in ocean warming (Free et al., 2019) in order to inform management under multiple and interacting stressors in a changing ocean.

To investigate the historical impacts of fisheries expansion in New Zealand waters, Durante, Beentjes & Wing (2020) calculated the MTL using New Zealand fisheries catch data from 1930–2014. The authors found a positive trend in MTL coinciding with a period of geographic expansion of the fisheries, followed by a negative trend in MTL, thus displaying signals of fishing down. Yet, the authors did not correct the MTL trend for the masking effects of geographic expansion of fisheries across the time series (Durante, Beentjes & Wing, 2020). In the present study, we do so by applying the Fishing-in-Balance (FiB) index (Pauly et al., 1998; Pauly, Christensen & Walters, 2000; Pauly & Watson, 2005; Bhathal & Pauly, 2008), and the adapted Mean Trophic Level index (aMTL, also known as the region-based Marine Trophic Index, Kleisner, Mansour & Pauly, 2014; Liang & Pauly, 2017) to New Zealand fisheries catch data from 1950–2019. By accounting for the geographic expansion of fisheries in New Zealand waters, we investigated whether fishing down occurred concurrent with fishery expansion. Further, by separating catch data between species’ habitat and gear used, we also investigated whether signals of warming in the MTC were apparent across the study period in line with trends in ocean warming.

Materials and Methods

The study area was the New Zealand EEZ (Fig. 1). To quantify the rate of ocean warming, we extracted Sea Surface Temperature Anomaly (SSTA) data from NOAA’s Kaplan Extended SST V2 database (Kaplan et al., 1998). At 5° latitude by 5° longitude resolution, monthly values were averaged across years from 1950–2019. These SSTA data are relative to temperatures during the period from 1951–1980 (Kaplan et al., 1998), and SSTA was chosen as temperature anomalies are considered to present a more consistent indicator of climate-change related ocean warming across large areas relative to absolute Sea Surface Temperatures (Tsikliras & Stergiou, 2014; Tsikliras et al., 2015; Dimarchopoulou et al., 2021b). We then performed simple linear regression and segmented linear regression of SSTA on time (Years). Segmented linear regressions were performed in order to test whether trends were better explained by multiple regression lines (Dimarchopoulou et al., 2021b), as indicated by the presence of a significant breakpoint (Muggeo, 2016). Breakpoints were identified using the segmented package in R (Muggeo, 2008; Team R Development Core, 2022), whereby a non-zero difference in slope on either side of the breakpoint was tested using a p-score test from the segmented package (Muggeo, 2008).

Primary analyses were completed using fisheries catch data for the New Zealand EEZ extracted from the Sea Around Us database (Pauly, Zeller & Palomares, 2020). This database includes reconstructed catches of officially reported statistics from the Food and Agriculture Organization of the United Nations (FAO), plus estimated unreported catch, discards and by-catch from all fishing sectors operating in the area, from 1950–2019 (Simmons et al., 2016; Pauly, Zeller & Palomares, 2020, Fig. 2A). In order to complete MTC and MTL analyses, species traits’, including mean preferred temperature and trophic level were gathered from FishBase (Froese & Pauly, 2021) via the rfishbase package in R (Boettiger, Lang & Wainwright, 2012; Team R Development Core, 2022). The mean preferred temperature (°C) for each species is based on a compilation of modelled species’ distributions informed by occurrence data, and also adjusted for appropriate depth ranges (Cheung, Watson & Pauly, 2013; Froese & Pauly, 2021). A species’ trophic level (TL) is estimated via diet information in FishBase (Froese & Pauly, 2021), and is represented by a discrete value usually between 2.0 and 5.0 (Liang & Pauly, 2017). The trophic level represents the number of trophic links a species maintains from primary production (TL = 1), i.e., herbivorous fish have a trophic level = 2, piscivorous fish have a trophic level = 3, and so on, with the maximum trophic levels in marine systems exceeding a value of five (i.e., top predators) (Stergiou & Karpouzi, 2002; Dimarchopoulou et al., 2021a; Eddy et al., 2021).

The number of species with available trait information contained in the catch data was 110. After identifying these species, we compared their Sea Around Us catch trends to officially reported statistics from the FAO (FAO, 2020). For several taxonomically-similar species, or species that were historically mislabeled, catch statistics have been reported under a group code or at higher taxonomic resolution (ex. Allocyttus niger and Pseudocyttus maculatus combined as ‘Oreo species,’ Simmons et al., 2016). In the present study, we restricted analyses to the species-level, as applications of the MTC and MTL at higher taxonomic coarseness have been demonstrated to mask the effects of fishing down trophic levels (Pauly & Palomares, 2005; Liang & Pauly, 2017). We therefore removed from analyses species that were previously reported in higher taxonomic groups (n = 7), were without available trait information (n = 19), and records labelled at a higher or an unidentified taxonomic resolution (n = 103) (Table S1).

Following this data treatment 102 species were included in analyses (Table 1). They represented 80% of the total species recorded and 48% of the total raw catch data contained in New Zealand’s Sea Around Us dataset (Fig. 2B). The majority of this analyzed catch data is of QMS species (96% of data), 67% of which is reported catch and 29% of which is unreported catch (i.e., reconstructed, Fig. 2B). We retained several ‘rare’ species (of low catch quantities) or bycatch-related species in our analyses. This was done as fish community composition can change over time, such as due to climate warming and shifts in species’ geographic distribution to new areas (Gordó-Vilaseca et al., 2023). We then calculated the MTC of New Zealand’s catch using the formula:

(1) $$MT{C}_{yr}={\displaystyle \frac{{\sum }_i^n{T}_i{C}_{i,yr}}{{\sum }_i^n{C}_{i,yr}}}$$where T_i is the preferred mean temperature of species i, C_i,yr is the catch of species i in year yr, and n equals to total number of species recorded (Cheung, Watson & Pauly, 2013). We then calculated the MTL of the same catch data, using the formula:

(2) $$MT{L}_{yr}={\displaystyle \frac{{\sum }_i^{m}T{L}_i{Y}_{ik}}{{\sum }_i^m{Y}_{ik}}}$$where TL_i is the trophic level of species i, Y_ik is the catch of species i in year k, and m is the number of species recorded (Pauly et al., 1998). For both MTC and MTL trends, we calculated simple and segmented linear regressions of each index on time (Year) across the study period. In order to relate MTC trends to trends in ocean warming, we computed Kendall’s τ correlation between MTC and SSTA using the stats package in R (Team R Development Core, 2022); with the null hypothesis that there is no correlation between the two variables. Thereafter, in order to first address the geographic expansion of fisheries and thus the masking of MTC and MTL trends, we calculated the Fishing-in-Balance (FiB) index. The FiB index is designed to identify surplus fisheries catch (i.e., from adjacent areas or stock assemblages) based on assumptions of energy transfer between trophic levels. The FiB was calculated using the formula:

(3) $$Fi{B}_k=lo{g}_{10}\left[Y_k{\left({\displaystyle \frac{1}{TE}}\right)}^{MT{L}_k}\right]-lo{g}_{10}\left[Y_0{\left({\displaystyle \frac{1}{TE}}\right)}^{MT{L}_0}\right]$$where Y is the catch and MTL is the mean trophic index value for year k, while Y₀ and MTL₀ are the catch and MTL value for the first year of data, and TE is the transfer efficiency between trophic levels (TE = 0.1), as estimated and assumed constant by Pauly & Christensen (1995) and used in Bhathal & Pauly (2008). Under this assumption of trophic transfer efficiency, a decrease in TL of 1 should yield a 10-fold increase in catch, or conversely, an increase in TL of 1 should yield a 10-fold decrease in catch (Kleisner, Mansour & Pauly, 2014; Liang & Pauly, 2017). As such, fisheries are considered to be “fishing in balance” (FiB = 0) when a decline in MTL value corresponds to the expected increase in catch, or vice versa, in line with assumptions of constant TE (Kleisner, Mansour & Pauly, 2014). However, when catch increases more than what is expected from decreases in MTL (FiB > 0), this implies that excess catch came from adjacent stocks and thus indicates the geographic expansion of fisheries (Kleisner, Mansour & Pauly, 2014; Liang & Pauly, 2017; Dimarchopoulou et al., 2021a). We also applied linear and segmented regression to the trend of FiB across the study period (Muggeo, 2008).

In order to quantify the geographic expansion of fisheries across the time series, we then computed the adapted Mean Trophic Level (aMTL) in line with Kleisner, Mansour & Pauly (2014) and Liang & Pauly (2017). The aMTL was developed under the assumption that geographic expansion occurs once a certain fished area or stock shows signs of depletion. Considering costs of fuel and equipment, it is assumed that most fishing activities are first concentrated within nearshore waters. Once the nearshore fisheries become overexploited, fishing effort then expands further away from shore and/or deeper (Kleisner, Mansour & Pauly, 2014; Liang & Pauly, 2017). Based on the FiB index (Eq. (3)), if fishing occurs “in balance”, where FiB = 0, with initial catch Y₀ and MTL₀, catch Y_k of year k can be computed:

(4) $$Y_k=Y_0\times {\left({\displaystyle \frac{1}{TE}}\right)}^{MT{L}_0-MT{L}_k}$$

As such, if $Y_k>Y_0\times {\left({\displaystyle \frac{1}{TE}}\right)}^{MT{L}_0-MT{L}_k}$, this indicates that a geographic expansion of the fisheries has occurred (Liang & Pauly, 2017), i.e., accessing previously unexploited assemblage of fish stocks in adjacent areas, habitats or depths (Kleisner, Mansour & Pauly, 2014). Upon identification of an expansion, a ‘node’ marks the year of occurrence, after which new aMTL values can be calculated for those newly identified stock assemblages (Kleisner, Mansour & Pauly, 2014). This identification of expansion operates on two main assumptions: that fishing does not cease in the original stock assemblages, and that fishing in the original stock assemblage is in balance (i.e., FiB = 0), or in decline (Kleisner, Mansour & Pauly, 2014).

Since the initial mean trophic level (MTL₀) in Eq. (4) may not be fully representative of all trophic levels present within the ecosystem, we corrected for this by assigning a range of possible trophic levels between the lowest TL (TL_lower) and highest TL (TL_upper) within the catch data (Liang & Pauly, 2017). Amongst J trophic levels, within the range [TL_lower, TL_upper], at each trophic level j, the catch potential p $Y_{kj}$ was calculated:

(5) $$\mathrm{p}{Y}_{kj}=Y_0\times {\left({\displaystyle \frac{1}{TE}}\right)}^{T{L}_j-MT{L}_k}$$while the corresponding maximum catch potential for year k is calculated:

(6) $$\mathrm{p}{Y}_{kj}={\sum }_{j=1}^J\left(\mathrm{p}{Y}_{kj}\times \mathrm{P}\mathrm{r}\left(\mathrm{T}{\mathrm{L}}_j\right)\right)$$with Pr( $\mathrm{T}{\mathrm{L}}_j$) equal to the probability that MTL₀ = TL_j (Liang & Pauly, 2017). By calculating p $Y_{kj}$, which is independent of MTL₀, we are able to estimate the maximum catch value that fisheries could achieve within a distinct stock assemblage, under the assumption of TE. When a reported catch Y_k > p $Y_k$, this indicates a geographic expansion has occurred, whereby year k is assigned node n_r, and r indicates the newly identified stock assemblage (Liang & Pauly, 2017). For any year and stock assemblage that follows n_r, catch and aMTL was calculated amongst distinct assemblages, with the maximum number of assemblages computed = 3, in line with the assumptions and caveats presented by Kleisner, Mansour & Pauly (2014). For a conceptual diagram of the workflow and calculations of the aMTL, see Liang & Pauly (2017).

To unmask the potential effects of the geographic expansion of fisheries on the MTC, we separated the fisheries catch based on species’ habitat classification (i.e., milieu) and also by the fishing gear used in order to re-calculate MTC trends for each group and gear type (Leitão et al., 2018). Since the overall MTC trend can be masked by the dominance of few (deep, cold-water) species, we separated catch by habitat associations and gear used, in an attempt to reduce the overall effect of dominant species on possible MTC trends of warming within different sectors of New Zealand fisheries. Species’ specific milieu identifications were done according to FishBase (Froese & Pauly, 2021, Table 1), while 23 distinct fishing gear classifications included in the catch data were aggregated into 12 general groups (Table S2). For the complete list of species caught by each gear type, see Table S3. We then calculate the MTC by gear, as the type of gear utilized will vary by the habitat exploited and thus determine the species composition of catch (Leitão et al., 2018). Similarly, we calculated MTC by milieu classification in order to group species based on similar environmental conditions and habitats utilized (Stephenson et al., 2020), and thus, distinguish between different habitats, geographic areas and species assemblages exploited throughout the period of fisheries expansion in the area. We also computed Kendall’s τ correlation between MTC values for each gear and milieu group and SSTA as well as the “full” MTC trend calculated for all analyzed species (n = 102).

Sensitivity analyses were completed in order to address the validity of the calculated metrics using the reconstructed Sea Around Us catch dataset. First, we computed the MTC, MTL, FiB and aMTL metrics using officially reported catch data from the FAO (2020). This dataset was composed of 42 species that were: both present in our calculations and also those of Durante, Beentjes & Wing (2020), had available trait information, and whose catch was reported at the species level (Table 1, Fig. 2C). Since gear information was not available in the FAO metadata, we thereafter calculated MTC by species’ milieu only. All species contained in this subset were QMS species, and represent 30% of the full catch data and 19% of all species/groups present in the FAO catch dataset. Second, we reduced the Sea Around Us dataset to those species contained in the FAO dataset (n = 42) and re-ran analyses. These 42 species represented 95% of the total analyzed Sea Around Us catch data (Fig. 2C).

Results

Across the study area and period, SSTA displayed a weak positive trend, increasing 0.04 °C per decade from 1950–2019 (Fig. 3A), with no significant breakpoint detected from segmented regression. When calculating the MTC of fisheries catch data, a segmented regression identified a significant breakpoint at 1986–whereby the MTC decreased 2.49 °C per decade from 1950–1986, while the MTC increased 0.44 °C per decade thereafter, i.e., from 1986–2019 (Table 2, Fig. 3B). Trends in MTL displayed a significant breakpoint at 1994, with the MTL increasing 0.15 per decade from 1950–1994, and decreasing from 1994–2019 (second slope = −0.03 per decade, Table 2, Fig. 3C). Results from the FiB index also reveal an overall positive trend, as the index increased 0.52 per decade until a detected breakpoint in 1999, while the subsequent slope decreased 0.11 per decade from 1999–2019 (Table 2, Fig. 3D). Thus, with modest levels of ocean warming across the study area, the MTC decreased sharply from 1950–1986, and has since gradually increased, whilst both the MTL and FiB indices increased sharply until the 1990s, and thereafter decreased gradually.

When calculating the aMTL, three distinct stock assemblages were identified when correcting for the geographic expansion of fisheries across New Zealand’s EEZ. In line with assumptions of the aMTL (Kleisner, Mansour & Pauly, 2014), the first stock assemblage includes the nearshore fisheries (Fig. 3E, green line), which maintained a stable aMTL from 1950–1964, thereafter decreasing until the 1990s. A second stock assemblage was identified in 1965 (Fig. 3E, blue line), indicating the first major geographic expansion of fisheries. In this case, fisheries started catching higher trophic level species compared to the initial nearshore assemblage, and aMTL also subsequently declined until the 1990s (Fig. 3E, blue line). The third stock assemblage was identified in 1969, which has since caught the highest trophic levels in New Zealand fisheries (aMTL values > 4) from 1969–2019 (Fig. 3E, red line).

Results from the aMTL confirmed the geographic expansion of New Zealand fisheries offshore throughout the study period. In turn, we separated fisheries catch data between species’ milieu and between fishing gear groups to re-calculate MTC trends. Results reveal a general pattern of MTC cooling for both bottom-associated species and fishing gears, with few groups displaying little to no MTC trend across the study period (Tables 3, 4, Figs. 4, 5). We do not report results for fishing gear groups ‘pots or traps’, ‘mixed gear’, ‘pelagic trawl’ or ‘pole and line’ as there was insufficient data to produce linear regression trends. The MTC for benthopelagic and demersal species displays strong reductions until 1979 and 1995, respectively (Table 3, Figs. 4C, 4D), with similar negative trends for bottom trawl (breakpoint: 1987), gillnet (breakpoint: 1999), longline (breakpoint: 1969), small scale (breakpoint: 1989) and unknown (breakpoint: 1979) fishing gear groups (Table 4, Fig. 5). Conversely, strong positive trends in MTC were observed for pelagic-oceanic species following 1965 (Table 3, Fig. 4F) and for purse seine fisheries (Fig. 5F) across the study period, as well as the ‘other’ fishing gear group after 1964 (Table 4, Fig. 5E). Overall, these results reflect the general trend of MTC cooling for bottom-associated species and fisheries, but also reveal warming MTC trends for pelagic and oceanic species across the study area.

Results from Kendall’s τ correlation found only the bathypelagic MTC trend to be significantly positively correlated with SSTA trends (Table 3), while the rest of the milieu and fishing gear MTCs were not correlated, including the primary MTC trend. Correlations with the primary dataset (n = 102) MTC trend were strongest for those milieu and gears that had strong MTC reductions, including benthopelagic and demersal species, as well as bottom trawl and small-scale gears (Table 4). These results show that milieu- or gear-specific MTC trends were not strongly associated with low levels of ocean warming across New Zealand waters, while the primary MTC trend is correlated with bottom associated species- and fishing gear trends of MTC cooling.

Sensitivity analyses calculating all of the aforementioned metrics using FAO catch data (n = 42) and a reduced Sea Around Us dataset (n = 42) reveal similar results from the primary Sea Around Us dataset (n = 102) analyzed (Tables S4, S5, Figs. S1–S4). This includes the identification of three stock assemblages via the aMTL, with the first two assemblages displaying fishing down trends, although the years of identified expansion into new stock assemblages are 6 and 8 years later, respectively, within the FAO dataset (Fig. S1D). Milieu- and gear-specific MTCs also reveal similar trends, with negative slopes in bottom-associated fisheries and species. For the pelagic-oceanic milieu MTC, FAO data displayed strong warming until 2009, thereafter decreasing (Fig. S2F), while the Sea Around Us data subset displayed a consistent positive slope (Fig. S4F). Although for purse seine fisheries, strong trends of MTC warming were not present within the Sea Around Us subset (Fig. S4L). Overall, sensitivity analyses support the validity of our findings from the primary Sea Around Us dataset.

Discussion

Our results have identified and unmasked aspects of MTC warming during the geographic expansion of New Zealand fisheries catch into offshore waters from 1950–2019. By separating fisheries catch between species’ habitat and by fishing gear, we observed warming trends of pelagic-oceanic species across the study period. Results from sensitivity analyses confirmed a warming trend of pelagic-oceanic species, although FAO data showed this trend to occur from 1950–2009, thereafter reversing (Fig. S2F). As well, sensitivity analyses using the reduced Sea Around Us dataset removed strong trends of warming in purse seine fisheries (Fig. S4L), and thus, must be considered when interpreting the primary results of purse seine MTC warming. Nonetheless, overall results reflect the increasing proportion of pelagic-oceanic species of warmer thermal affinity caught in New Zealand from 1950 until at least 2009. This trend was otherwise masked in the full MTC by the geographic expansion and deepening of fishing operations (Fig. 3B). Strong correlations between bottom-associated habitat- and gear-specific MTC trends and the full, primary MTC suggests that cooling resulted from deepening of New Zealand’s catch. This similarly masked result in New Zealand waters was reported by Cheung, Watson & Pauly (2013) during the initial application of the MTC.

Following the calculation of the primary, masked MTC, in a step-wise process we calculated the MTL of the same catch data. These results were in line with those reported by Durante, Beentjes & Wing (2020) who also analyzed New Zealand fisheries, displaying a strong positive trend of MTL values until the early 2000s, and thereafter decreasing. This increasing trend of MTL, that indicated harvesting of higher trophic level species, corresponds to the period of geographic expansion and deepening of fisheries into colder waters in the area (Leitão, 2015), as new habitats and species (i.e., of higher trophic levels) were gradually included into New Zealand fisheries (Durante, Beentjes & Wing, 2020). Thereafter, incorporating the MTL into FiB and aMTL indices, we identified two discrete periods of fisheries expansion, with the first occurring in 1965, and the second in 1969.

The geographic expansion of New Zealand’s fisheries tracks the modernization of the country’s fisheries technology as well as policy. By the 1950s and 1960s, technologically-accessible nearshore fisheries were approaching saturation, and began to show signs of overexploitation (Stock Assemblage 1, Fig. 3E), coinciding with widescale deregulation measures for nearshore fisheries in 1963 (Clark, Major & Mollett, 1988; Pinkerton, 2017; Durante, Beentjes & Wing, 2020). Deregulation facilitated the rapid development of a domestic fisheries industry within the newly established 12 nautical mile Territorial Sea in 1965 (Bradstock, 1979; Clark, Major & Mollett, 1988; Fig. 1). This led to an expansion of fishing effort in shelf waters (Jackson et al., 2001; Kleisner, Mansour & Pauly, 2014), and coincided with our identification of the geographic expansion to include the second stock assemblage in 1965 (Fig. 3E). As we assume that this expansion occurred within shelf waters (sensu Kleisner, Mansour & Pauly, 2014), these two stock assemblages are parallel geographically (along New Zealand’s shelf) and thus maintain similar species assemblages (in contrast to the offshore Stock Assemblage 3). In turn, our results show similar trends of fishing down to have occurred in Stock Assemblages 1 and 2, which represents the full exploitation of species and thus the range of trophic levels available to fisheries within shelf waters during this period of fisheries expansion in New Zealand. A similar trend was observed between the first and second stock assemblage identified during the analysis of geographic expansion of fisheries in Indian waters (Kleisner, Mansour & Pauly, 2014).

As shelf-water fisheries continued to develop throughout the 1960s, high trophic level species were being depleted, thus Stock Assemblages 1 and 2 began to show signs of fishing down (Fig. 3E). This trend was revealed by the application of the FiB and aMTL indices, and was otherwise masked by the calculation of MTL of fisheries catch across the New Zealand EEZ (Fig. 3C, Durante, Beentjes & Wing, 2020). Similar results from Leitão (2015) reported signals of overfishing in nearshore stocks during a period of fisheries expansion, and thus support our results, as well as the assumptions of aMTL analyses, in that persistent fishing down in nearshore and shelf stock assemblages (i.e., Stock Assemblage 1 and 2) incentivizes expansion into previously unexploited stock assemblages (i.e., Stock Assemblage 3) (Kleisner, Mansour & Pauly, 2014). Concurrently, through further economic incentives including subsidies for new vessels, new processing plants, and an expanded list of exploitable species, New Zealand’s fisheries fleets expanded further offshore (Johnson & Haworth, 2004; Durante, Beentjes & Wing, 2020), coinciding with the expansion of fisheries into the third identified stock assemblage in 1969 (Fig. 3E). This third stock assemblage has maintained the highest mean trophic levels of catch across the time series (Fig. 3E), likely due to consistently high catches of the commercially-important hoki (Macruronus novaezelandiae, TL = 4.5, Table 1) since the early 1970s.

By completing sensitivity analyses with the reduced Sea Around Us catch data (n = 42) and the FAO catch data (n = 42), our results of fishing down within nearshore stock assemblages was confirmed throughout the period of fisheries expansion in New Zealand waters. The identification of node expansion into new stock assemblages occurred later within the aMTL from FAO data. This is expected, as FAO reported catch quantities are lower than those for the reconstructed Sea Around Us catch, including but not limited to the period from the 1950s to 1970s (Fig. 2C). Moreover, it is documented that the identification of node expansion can occur later if catch increases more gradually (Kleisner, Mansour & Pauly, 2014), as is the case for the analyzed FAO catch data.

By separating and re-calculating the MTC between species’ milieu and fishing gear type, MTC trends are in line with results reported via the aMTL. Several fisheries with bottom-associated gear or species, including deep-water species, displayed negative MTC trends from the mid- to late 1960s until the 1980s and 1990s (Figs. 4, 5). This cooling of the MTC coincides with the expansion of fishing operations to include Stock Assemblage 3 from 1969 onwards. As gear development promoted further expansion into deeper and cooler waters, new catches of barracouta (Thyrsites atun), southern blue whiting (Micromesistius australis) and hoki drove down MTC trends (Figs. 6A, 6B). By the 1980s, catches of hoki dominated not only gear- and habitat-specific fisheries catch, but also overall catch in New Zealand waters. This led to the levelling off of both hoki-dominated gear- and habitat-specific MTCs (e.g., Figs. 4C, 5A) and also the general trend of New Zealand’s MTC (Fig. 3B).

For pelagic-oceanic species and pelagic purse seine fisheries, trends of MTC warming were driven by the increasing catch of the subtropical skipjack tuna (Katsuwonus pelamis) and blue mackerel (Scomber australasicus), paired with the decreasing catch of the temperate southern bluefin tuna (Thunnus maccoyii, Figs. 6C, 6D). These results are in line with observations of pelagic species undergoing relatively rapid shifts in geographic distribution with oceanographic changes (Champion, Brodie & Coleman, 2021), albeit the changes in SSTA of the present study area are modest. Such changes in the composition of pelagic-oceanic species present in fisheries catch may be driven by the adult fishes’ ability to shift synchronously with their pelagic habitat to track local climate velocities (Pinsky, Selden & Kitchel, 2020; García Molinos et al., 2022). Moreover, these results suggest that the observed levels of ocean warming may favor the increased abundance of subtropical, smaller-bodied skipjack tuna and blue mackerel vs the temperate and larger-bodied southern bluefin tuna. In line with these results, and in recognition of the complexities of variable shifts in geographic distribution between different species (Dunn et al., 2022), we recommend that further research investigates the variability of shifts in geographic distribution of teleost fish species in line with ocean warming across New Zealand waters.

Our results of the positive trend in SSTA are in line with observations of the gradual warming of New Zealand waters. Overall, ocean warming within New Zealand’s EEZ has exhibited relatively weak trends and significant interannual variability due to the dynamic influences of the warm Subtropical Front (STF) from the north, north-west and the cold Subantarctic Front (SAF) from the south (Shears & Bowen, 2017; Sutton & Bowen, 2019). Management must consider the additive or synergistic effects of both fishing pressure as well as ocean-warming-induced changes in productivity and geographic distribution of fish populations. For example, fishing effort may be reduced on stocks that are overexploited and/or are shifting the trailing edge of their geographic distribution away from historically fished areas as they track their thermal affinity (Szuwalski & Hollowed, 2016; Gaines et al., 2018; Lamine et al., 2022). In the present study, this may apply to the cold-water southern bluefin tuna, which displays strong reductions in overall catch from the 1970s onward (Fig. 6C). This may also be relevant as other commercially-important species such as hoki are predicted to shift their geographic distribution polewards with future ocean warming (Dunn et al., 2022). Conversely, the leading edge of pelagic, warm-water skipjack tuna and blue mackerel may be penetrating deeper into New Zealand’s EEZ, thus facilitating higher catches. The joint application of the present indices can help inform such management decisions.

Interpretation of trends in the MTC should also consider historical fishing pressure and management within this particular system. Demersal or reef-associated species like the Australasian snapper Pagrus auratus were already displaying signs of overexploitation by the 1970s (Durante, Beentjes & Wing, 2020), as confirmed by our calculation of the nearshore and shelf aMTLs (Stock Assemblage 1 and 2, Fig. 3E). As a result, New Zealand implemented their quota management system (QMS) in 1986 (Clark, Major & Mollett, 1988), and since then, overall fisheries catches have decreased (Fig. 2A), and that of Australasian snapper has remained steady but well-below its maximum catch values seen in the 1960s and 1970s (Supplemental Files). Conversely, since the 1980s, skipjack tuna and blue mackerel have displayed oscillating but gradually increasing catches (Figs. 6C, 6D). This may be a result of the then newly-implemented QMS, the technological and geographic expansion of fisheries in New Zealand waters, as well as trends in ocean warming and the increased abundance of pelagic, warm-water fishes in the area. We therefore recommend that interpretation of the present MTC warming trend considers both the thermal affinity tracking of pelagic species with ocean warming, as well as fisheries management in New Zealand adapting to stock changes.

Conclusions

Overall, our results build on previous research of the MTL to correct for the geographic expansion of fisheries within the New Zealand EEZ from the 1950s onward. Using the FiB and aMTL indices, we identified trends of fishing down trophic levels in nearshore fishing stocks, a trend that was otherwise masked within the full calculation of New Zealand’s MTL (Durante, Beentjes & Wing, 2020). Thereafter, by separating catch between species’ milieu and fishing gear groups, increasing catches of subtropical skipjack tuna and blue mackerel revealed positive MTC trends for pelagic species and fishing gear. These trends were also otherwise masked by the offshore expansion and deepening (into colder depths) of fishing operations and catch, as species such as hoki came to dominate overall catches across the same time period. These results highlight the pervasiveness of fishing down food webs in New Zealand fisheries throughout a period of technological development and fisheries geographic expansion. These results also reveal how fisheries catch composition may change due to the responsiveness of pelagic, mobile species to track local climate velocities with ocean warming.

Supplemental Information