Introduction
Size structure plays an important role in maintaining reproductive potential and stability of fish population. For example, larger individuals tend to produce more and better eggs (Hsieh et al. 2010b, Hixon et al. 2014) and have a longer spawning season (Berkeley and Houde 1978, Lambert 1987, Hutchings and Myers 1993, Sogard et al.
2008); small and large individuals may spawn at different sites (Trippel et al. 1997, Lawson and Rose 2000, Vandeperre and Methven 2007). Such bet-hedging strategies provide resilience capacity for populations to sustain under unfavorable conditions (Hsieh et al. 2010b, Ripa et al. 2010, Schindler et al. 2010). Hence, investigating the change of size structure may provide insight of how resilient a fish population can be.
Several external forcings may alter the size structure of a fish population. The most well-known examples are fishing and temperature. Fishing represents size-selective removal of larger individuals that can truncate the size structure of a fish population (Barnett et al. , Bianchi et al. 2000, Berkeley et al. 2004, Ginter et al. 2015), which in turn may cause recruitment failure (Ottersen et al. 2006), reduce the reproductive outputs (Scott et al. 2006), and increase variability of fish populations (Hsieh et al.
2006b, Anderson et al. 2008, Rouyer et al. 2012). It may also lead to evolutionary consequences (Heino et al. 2015); for example, the genetic differences found in the populations of Atlantic cod (Gadus morhua) from Iceland is due to difference in depth-associated fishing mortality (Jakobsdóttir et al. 2011). As such, balanced
doi:10.6342/NTU201802401
!
!
!
26!
exploitation (Law et al. 2012) and harvest-slot-limit (Gwinn et al. 2015) have been proposed to prevent fishing-induced size truncation.
Apart from fishing, increasing sea water temperature caused by global warming may also lead to shrinking size structure of marine fish populations (Daufresne et al.
2009b, Cheung et al. 2013a). Elevated ambient water temperature directly influences fish metabolism at individual level (Gillooly et al. 2001a), which increases growth rate and causes earlier maturation at population level (Neuheimer and Grønkjær 2012). Also, temperature may indirectly influence the recruitment processes through trophic transfer (Sundby 2000) and thus change the size structure. Based on the match/mismatch hypothesis, the larvae survivorship relates to the match between the timing of larval feeding and the food production (Cushing 1990). For example, rising temperature since mid-1980s has modified the plankton ecosystem and reduced the survival of young cod in the North Sea (Beaugrand et al. 2003b). While the fishing and temperature effects have been well documented, their relative contributions on the size structure of fish populations remain poorly understood. This is a critical issue particularly for exploited stocks, because overfishing has been shown to enhance the sensitivity of fish abundance and distribution to climate (Hsieh et al. 2008b, Hidalgo et al. 2011, Rouyer et al. 2014);
nevertheless, whether such synergistic effect also occurs in size structure remains relatively unexplored.
Previous studies on quantifying fishing or temperature effects on the size structure of fish populations have been focused on univariate size-based indicators (SBI). Some studies used the upper 95-percentile of the length frequency (L95) to test fishing effect (Rochet et al. 2010, Brunel and Piet 2013), while the other used length class diversity to investigate the stability of population through time (Marteinsdottir and Thorarinsson 1998). Recently, European Commission Marine Strategy Framework Directive required
regional (or local) fishery reports to provide the information of basic SBIs (e.g. the mean length) in order to improve the management and maintain the sustainable development (Farmer 2012). However, it remains unclear whether these univariate SBIs could represent the entire size structure and the status of a population. It has been suggested that no single SBI can represent an effective overall indicator for external forces (Shin et al. 2005). Also, SBIs need to be selected carefully based on their implications. For example, L95 can only reflect the variety of large fish in fish population. The analysis with the North Sea cod, herring and plaice found L95 failed to reveal the effects of external forces on fish population, as it was rather insensitive in responding to fishing mortality (Brunel and Piet 2013).
To overcome the limitation of existing SBIs, it needs an alternative approach to (1) analyze complete information of size structure and (2) examine how external forces affect the size structure. Here, I employed the variation partitioning approach to conduct a size structure-based analysis that examines the variation of size class composition in response to external forcings. Variation partitioning can be best understood as a method for extending multivariate regression. In multivariate regression (y~x), y represents a univariate response vector and x represents multiple predictors, x1, x2, etc. (each is a vector) and possibly their interactions; the contribution of each predictor variable (xi) can be evaluated by partial R-square. Whereas in variation partitioning (Y~X), Y represents a multivariate matrix and X represents multiple predictors, X1, X2, etc. (each is a matrix); the contribution of each predictor matrix (Xi) and their interaction is also evaluated by partial R-square (Peres-Neto et al. 2006). Variation partitioning is commonly used in community ecology to examine the relationship between species composition (Y matrix) and various sets of explaining variables (e.g. 2 or 3 predictor matrices) (Peres-Neto et al. 2006). This method has also been extend to analyze
doi:10.6342/NTU201802401
!
!
!
28!
temporal and space-time variation of community composition data (Griffith and Peres-Neto 2006). Here, I borrow this concept to analyze temporal variation of size composition data in responding to fishing, temperature, and the interaction, with the simplification that fishing and temperature is just a vector. Specifically, for a given fish population, I apply the variation partitioning to quantify how the temporal variation of their size composition responds to fishing, temperature and their interaction (see Methods). The explained fraction of variation (partial R-square) by each factor then allows us comparing their relative contribution in affecting length composition through time.
Next, I perform a cross-stock meta-analysis linking the relative contribution of fishing or temperature (the output of variation partitioning as explained fraction of variation) to the life history traits of fishes, as well as long-term mean and variability of fishing or temperature (see Methods). This meta-analysis aims to examine which factor can explain the relative contribution of fishing, temperature and their interaction across stocks. This meta-analysis is motivated by the fact that life history traits are associated with the size structure of population (De Roos et al. 2003). I hypothesize that the large, slow growth, and late-matured species is more likely to be impacted by fishing in their size structure because size-selective removal (i.e. size truncation) may be more severe in these species (Rouyer et al. 2011) and their recovery will take longer time (Jennings et al. 1999, De Roos et al. 2003). I also hypothesize that small species is more vulnerable to temperature effects, because smaller species are more sensitive to temperature changes due to the constrains from metabolic allometries (Gillooly et al.
2001b).
Furthermore, I expect that fishing and temperature might exhibit interactive effects via multiple ways (Perry et al. 2010, Planque et al. 2010). For example, the long-term
fishing effect, such as long-term mean and variability of mortality ratio (fishing mortality divided by natural mortality, F/M) may affect the relative contribution of temperature in explaining temporal variation of size structure. Here, I standardize fishing mortality by natural mortality in order to have a fair cross-stock comparison.
Motivated by previous studies showing that fishing elevated sensitivity of exploited stocks to environmental changes (Hsieh et al. 2006b, Anderson et al. 2008, Rouyer et al.
2011), I hypothesize that the fish stocks experienced higher fishing pressure is more responsive to temperature effect in their size structure. I also hypothesize that habitat conditions, including mean and variability of temperature, affect the relative contribution of fishing effect. For instance, Wang et al (Wang et al. 2014) found that temperature affects the cod’s life history trait, making the cod population more vulnerable to fishing.
My objectives are, first, to apply variation partitioning to quantify how the variation of size structure responded to fishing, temperature and the interactive effects for 28 exploited stocks (Table 1) living in a wide range of habitats, including the west coast of US, Alaska, and North Sea (Figure S1). Secondly, I linked the fraction of explained variation by fishing (or temperature) to life history traits (including von Bertalanffy growth rate (K), length infinity (Linf), age at maturation (A50), length at maturation (L50)), as well as long-term mean and variability of fishing and habitat temperature conditions. Finally, to demonstrate the efficacy of our approach, I compared the performance of variation partitioning approach with the traditional univariate SBIs analyses. This comparison is straightforward, as both the univariate SBIs analysis and variation partitioning are computed using similar linear modeling of variance/covariance, with the difference only in the response variable- the response variable is a vector (y) in the univariate SBIs analysis whereas the response variable is a
doi:10.6342/NTU201802401
!
!
!
30!
matrix (Y) in the variation partitioning; two methods have the same explaining variables (i.e. fishing and temperature).
Material and methods
Size structure data of commercial stocks
I collected length frequency data from 28 exploited stocks, which contain temporal coverage over 20 years and annual fishing mortalities (or exploitation rates) estimated by stock assessment are available (Table 3-1, Table S1). These stocks came from 3 regions in the northern hemisphere (Figure S1): (1) the west coast of US (West US), which is part of Northeast Pacific; (2) Alaska, which separates into 3 fishing areas- Aleutian Islands (AI), Gulf of Alaska (GOA), and Bering Sea; and (3) North Sea.
Collectively, these 28 stocks in 3 regions were well studied, spanning distinct distribution of size structure, with a wide range of life-history traits and habitats, and therefore are representative of a compilation of global-scale fish stocks (Table S2).
I primarily use length frequency from the fisheries-independent surveys for each stock. These are the length frequency per size range of the given year as arranged in the stock assessment report. For Alaska region, the survey data includes bottom trawl survey in Aleutian Island (AI), East Bering Sea shelf (EBS) and Gulf of Alaska (GOA).
For the North Sea, I compiled the length frequency of ICES International Bottom Trawl Survey in the North Sea (NS-IBTS). I used the 1st quarter (Q1) in NS-IBTS for consistency because there was only annual Q1 survey prior to 1991. For the west coast of US, I used fisheries-dependent length compositions instead of bottom trawl survey because the data from fishery-independent surveys during 1980-1990 were limited.
Fishing and natural mortality
To quantify the fishing effect, I used time series of annual fishing mortality from the stock assessment reports (Supplementary Table S3-1). I focus on the single fishery whenever possible to minimize the uncertainty in fishing selectivity due to changing gears. I noted that analyses in the west coast of US use exploitation rate instead of fishing mortality in stock assessment. To make a fair comparison for meta-analysis, I transformed the exploitation rate into fishing mortality through the relationship between mortality and survival rate in fisheries (Ricker 1975). Here I first assume that these fisheries are type II fishery, in which the fishing and nature mortality operate concurrently. The exploitation rate (µ), fishing mortality (F), natural mortality (M), instantaneous total mortality rate (Z) and actual total mortality (A) have following relationships:
µ=F·A/Z (1) Z=F+M (2)
A=1-e (F+M) (3) the equation (1) can also be written as:
µ =F· (1-e (F+M))/(F+M)
With the known exploitation rate and natural mortality, this equation can be solved numerically and yields the fishing mortality. The natural mortality here is mostly compiled from the value of preferred model in the stock assessment reports (Supplementary Table S3-1, Table S3-2).
Temperature
In analysis, I primarily used empirical measurements of temperature along with the trawl surveys (see Supplementary Figure S3-2). The North Sea (53-59°N, 3°W-10°E) near-bottom temperature is the station observations of hydrochemical measurements
doi:10.6342/NTU201802401
!
!
!
32!
from ICES Oceanographic database (http://ocean.ices.dk/HydChem/HydChem.aspx) at Q1 (January-February), with coverage of almost entire North Sea. For the Alaska region, Aleutian Islands, East Bering Sea and Gulf of Alaska, the bottom trawl surveys provide bottom temperature measurements (data available at website of NOAA Alaska Fisheries Science Center, https://www.afsc.noaa.gov/RACE/groundfish/survey_data/data.htm). I took average for all the stations that the species occurred. Although many demersal species examined here have a pelagic phase at larval-juvenile stages and some species (e.g. Atlantic cod) forage on the whole water column, the preliminary analysis on surface (see Figure S3) and bottom temperature suggests that they are highly correlated (AI: 0.92, EBS: 0.72, GOA: 0.67, North Sea: 0.76). This suggests even the water is stratified, the interannual variability is very similar in temperate-sub arctic ocean.
Therefore, I used only the bottom temperature collected in trawl survey for these 4 areas in our analysis. For the west coast of US, I used the sea surface temperature from NCEP/NCAR Re-analysis monthly mean (data available at website of NOAA Earth System Research Laboratory http://www.esrl.noaa.gov/psd/data/timeseries/) and then took aerial average (31.4-48.6°N and 54.4-56.3°W) to represent the entire region.
Variation partitioning to quantify fishing and temperature effect on size structure
To determine how much of the temporal variation in size structure is explained by fishing and temperature effect, variation partitioning (aka. redundancy analysis) was used to decompose the total variation of length composition of each stock through years.
The response variable here is the matrix of size composition data through time. The explanatory variables were time series of annual fishing mortality and temperature. The unbiased estimation of adjusted R2 (accounting for sample size effect) provides a test with correct Type I error rate and good power for redundancy analysis (Peres-Neto et al.
2006). I reported the adjusted R2 for the pure effect of fishing, temperature and their interaction. Next, I perform permutation test (1000 times) to evaluate the significance of each fraction. Although the significance test cannot be done for the interactive component (Peres-Neto et al. 2006), I still report the interactive component for the sake of comparison.
To incorporate the lag effect of temperature, I additionally used 1 year- and 3 year-lag temperature as explanatory variables in variation partitioning. For the stocks in Alaska, sampling interval is 2- or 3-year, and thus I only considered 2 or 3-year lag. I did not investigate lagged effect of fishing; I assumed that fishing instantaneously affects adult population while temperature may affect the size structure through influencing the future recruitment at early life stages (e.g. egg and early larvae).
Because variation partitioning does not estimate log-likelihood in the procedure, the best model is selected according to the largest effective size (highest total adjusted R2 (see Table 2)). The fraction of variation explained by fishing and temperature from the best model was used in further analyses.
Variables affecting the relative contribution of fishing and temperature effect
To investigate what determines the relative contribution of fishing and temperature across stocks, I considered the following variables (see Supplementary Table S3-2): (1) life history traits: von Bertalanffy growth rate (K), length-at-infinity (Linf), age at maturation (A50) and length at maturation (L50); (2) indices of fishing effect: long-term mean and variability of mortality ratio (i.e. fishing mortality/natural mortality, F/M);
and (3) indices of temperature effect: long-term mean and variability of temperature. I used the mortality ratio (F/M) to reflect the fishing strength among different stocks.
Because fisheries management usually sets optimal fishing mortality (Fopt) proportional
doi:10.6342/NTU201802401
!
!
!
34!
to natural mortality (M) (Walters and Martell 2002), here I defined mortality ratio as fishing mortality normalized by natural mortality (F/M). When comparing multiple stocks with different natural mortalities, the mortality ratio can better reflect the impact of fishing pressure for cross-stock analyses.
To examine whether life history traits, long-term mean and variability of fishing and temperature have influenced the relative contribution of fishing and temperature effect, I first used simple univariate regression analysis to test each variable above. I then build a linear mixed effect model (LMM) with the same variable as fixed effect and habitat as random effect to check if the observed pattern remains.
Comparing the performance of univariate SBIs to variation partitioning
To demonstrate the efficiency of variation partitioning, I compared the P value of explained fraction in the variation partitioning to the P value in regression with univariate SBIs. The univariate SBIs considered in the analysis were: 95% percentile of length class (L95), mean length, Shannon diversity, Pielou’s evenness index, and skewness. After calculating the univariate SBIs, each SBI was fitted to a regression model where SBI ~ temperature + fishing + temperature*fishing. If the P value of variation partitioning were lower than that of multilinear regression, it would suggest that variation partitioning is more efficient than univariate SBIs in terms of rejecting the null hypothesis.
Computation and data availability
All computation was done with R version 3.3.3. Variation partitioning was carried out using R package vegan (Oksanen et al. 2017). The linear mixed effect model was done with R package lme4 (Bates et al. 2015). Further model testing and model
diagnostics were done with R package lmerTest (Kuznetsova et al. 2017). The original data (except few stocks via personal communications, see supplementary Table S1 in appendices for details), the R scripts to carry out all analyses, and the outputs are available at: https://zenodo.org/record/1211120.
!
Results and Discussion
My results of variation partitioning showed that the variance of size structure could be appreciably explained by fishing (on average of 10.9%), which is significantly higher than that of temperature and interaction (P=0.019 and P=0.00023 in ANOVA with paired t-test) (Figure 3-1). Specifically, 12 out of 28 stocks were significantly affected by fishing while 7 stocks were significantly associated with temperature (Table 3-2);
whereas, the interactive effect is small in most of stocks (Table 3-2).
I also found difference in fishing and temperature effect among regions (Figure 3-2) and habitat types (Figure 3-3). Fraction of variation explained by fishing is significantly different from temperature (P= 0.005) and interaction (P<0.001) in the west US. Both fishing and temperature effects are not significantly different from interaction in Alaska (P=0.74 and P=0.73) and North Sea (P=0.27 and P=1) (Figure 3-2). Among habitat types, fishing effect on size structure was significantly different from both temperature and interaction for demersal species (P=0.0014 and P <0.001 ANOVA with paired t-test) (Figure 3-3).
I employed the variation partition approach to analyze the time-series data of size structure of exploited stocks and found that fishing effect on size structure prevails (Table 3-2). Such commensurate result goes with the study that analyzed the stock biomass of 28 stocks in the Northeast Atlantics (Rouyer et al. 2014); that is, in a large-scale, fishing is the most critical factor. My analyses echo the increasing concern
doi:10.6342/NTU201802401
!
!
!
36!
over effects of fishing on size structure of exploited stocks (Blanchard et al. 2005, Planque et al. 2010, Cheung et al. 2013a, Ginter et al. 2015).
I note however, many of these studies applied univariate size-based indicators (SBIs) as proxies for the change of fish size structure. My approach incorporates the full size structure information in the analysis without assuming distribution of size data.
Comparison between variation partitioning and univariate SBIs suggests that variation partitioning is more efficient in rejecting the null hypothesis (p-value is smaller) than the univariate SBIs (Figure 3-5), particularly in detecting temperature effects (Probability of success=0.63, P=0.003) although only marginally in fishing (Probability of success=0.56, P=0.10). In other words, variation partitioning can be a useful complementary method to investigate the external forcings on size structure of fish populations, in addition to univariate SBIs.
My multi-stock meta-analytical framework also allows us to investigate how life history traits, and long-term mean and variability of fishing and temperature influence on the explained variation of size structure responding to fishing or temperature. The stocks experienced higher fishing variability (CV of mortality ratio) is more responsive to temperature effect in their size structure (Table 3-3, Figure 3-4d, P=0.038). This supports my hypothesis that fishing elevates the sensitivity of exploited stocks in responding to environmental changes.
Surprisingly, I found none of the life history traits is able to explain the fishing and temperature effect (Figure S4, S6). Also, the total adjusted R2 (Table 3-2) suggests that fishing plus temperature explained at most 44% of variation among 28 stocks. There may be other factors associated with body size, such as oxygen limitation of thermal tolerance (Pörtner and Knust 2007), affect the response to fishing/temperature effect in the size structure.
Through the application of variation partitioning, I had expected the efficacy to identify interactive effect of fishing and temperature on the size structure. However unexpectedly, my results indicate that the interaction effect is very weak (Figure 3-1).
This finding may superficially be interpreted as evidence for lack-of interaction of fishing and temperature on size structure because fishing effects have dominated.
However, I caution the interpretation of this finding, as variation partitioning is a linear variance decomposition method, which cannot account for nonlinear interactions.
While I demonstrate the efficacy of our size-structure based, meta-analytical framework to examine fishing and temperature effects on size composition of exploited stocks, I shall point out some caveats in our study. First, fishery-dependent data may lose some information due to the discard of small-size individuals. Second, there were fewer pelagic species than demersal and benthic species in our dataset, which might
While I demonstrate the efficacy of our size-structure based, meta-analytical framework to examine fishing and temperature effects on size composition of exploited stocks, I shall point out some caveats in our study. First, fishery-dependent data may lose some information due to the discard of small-size individuals. Second, there were fewer pelagic species than demersal and benthic species in our dataset, which might