Exploring mechanisms of spatial segregation between body size groups within ﬁsh populations under environmental change

Ample evidence has indicated shifts in distribution of ﬁsh populations in response to environmental stress. However, most studies focused at the whole population scale. This neglects the spatial dynamics between groups of diﬀerent body size (body size groups), that fundamentally shapes the spatial structure of a population. Here, we explored the mechanisms that modulate spatial dynamics of body size groups, and applied our analyses to three North Sea ﬁsh populations which experienced severe declines in biomass from 1977 to 2019: Atlantic cod (Gadus morhua), haddock (Melanogrammus aegleﬁnus), and whiting (Merlangius merlangius). All three populations exhibited strong declines in the overlapped area between body size groups in winter over 43 years, yet their mechanisms diﬀered. These declines were either due to (1) diﬀerent magnitudes of contraction of the distribution area of body size groups; and/or (2) diﬀerent speeds and directions of spatial shift among various body size groups, both increasing spatial segregation within populations. These patterns were either associated with ocean warming, and/or declining population biomass, and such associations often varied according to distinct body size groups. Our analytical approach provides a powerful tool for identifying vulnerable populations under environmental stress and can be generalized to study a variety of size/age structured populations at various ecosystem types.


Introduction
Many marine fish populations have undergone significant shifts in their spatial distributions over the past decades, largely related to ocean warming and declining population size (Perry et al. 2005, Sunday et al. 2015).Most of these studies focused at the whole population scale; however, several lines of evidence have suggested that the spatial shift varies in magnitude and direction for different body size groups within a population (hereafter, body size groups) (Bell et al. 2015, Barbeaux and Hollowed 2018, Frank et al. 2018, Yang et al. 2019, Li et al. 2022).For instance, the distribution of the middle size groups of some fishe populations in the Eastern Bering Sea shifted at a greater speed in warm seasons, compared to groups of smaller or larger body sizes (Barbeaux and Hollowed 2018).Another study across North Pacific, North Atlantic, and South Atlantic suggested that the distribution of large size groups within some fish populations shifted deeper, as a result of size-selective fishing at shallower water (Frank et al. 2018).These size-dependent shifts in distribution are likely to reduce overlapped areas between body size groups, that is, increase the spatial segregation within populations.However, temporal changes in spatial segregation (i.e., overlapped area) between body size groups have not been quantified for real-world populations, despite earlier efforts from theoretical approaches (Hughes and Grand 2000).
Changes in spatial segregation between body size groups of a population have various consequences on population dynamics.On one hand, a population with high spatial segregation between body size groups can reduce the stress from predation and competition.
On the other hand, a population with highly segregated size group is more vulnerable to local perturbations.These perturbations include size-selective fishing, size-selective predation, or unfavorable habitat conditions for certain body size groups (Hsieh et al. 2010b).
These perturbations can change the abundance of certain body size groups, which in turn alter the demographic structure and spatial structure of the population (Tao et al. 2021).
More generally, changes in the spatial structure of a marine population can influence life history and demographic variations, which potentially affect its resilience to perturbations (Ciannelli et al. 2013).
What are the potential mechanisms shaping spatial segregation between body size groups of a population?Within a population, the distribution area of each body size group, and the distance between their abundance-weighted centers of distribution area (hereafter, centers of abundance), determine the overlapped area between them.On the one hand, when the distribution areas of two body size groups contract, their overlapped area declines, provided that their centers of abundance are fixed.On the other hand, elongated distance between the centers of abundance reduces the overlapped area between body size groups, provided that their areas of distribution are fixed.
Ocean warming and population decline potentially impact the area of distribution and the center of abundance of body size groups (Barnett et al. 2017, Orio et al. 2017).
These impacts are likely size-specific.For example, earlier studies showed that ocean warming and fishing altered the abundance of body size groups at various extents (Barnett et al. 2017, Orio et al. 2017).Such size-specific changes in abundance could lead to differential changes in their area of distribution, based on abundance-distribution relationships and density-dependent habitat selection (MacCall 1990, Fisher and Frank 2004, Thorson et al. 2016).In addition, previous findings suggest that ocean warming and fishing contributed to size-specific shift in spatial distribution (Barbeaux andHollowed 2018, Frank et al. 2018).This is due to thermal tolerance, food requirements, spatial constraints, and mobility that vary with body sizes within a population (Dahlke et al. 2020, Ciannelli et al. 2022).Depending on the original positions of the center of abundance, the size-specific shift could increase the distance between their distribution.Linking body size-specific distribution response to ocean warming and population decline is key to understanding the mechanisms behind the changes in spatial segregation between body size groups of a population.
In this study, we quantified and explored the mechanisms of changes in spatial segregation over time between body size groups within fish populations.We asked the following question: did the overlapped area of body size groups within populations decline over time, and what are the mechanisms behind?We studied fish populations in the North Sea, a global warming hotspot that has experienced rising sea surface temperature over the past decades (Hobday and Pecl 2014).Particularly, we focused on those fish populations that are ecologically and economically important and experienced large geographical redistribution over the past century (Huserbråten et al. 2018), including Atlantic cod (Gadus morhua), haddock (Melanogrammus aeglefinus), and whiting (Merlangius merlangius).
The total biomass of these populations has declined since 1980s with slow recovery in recent years (Engelhard et al. 2014).Therefore, these populations are prone to distribution area contraction and fragmentation.We analyzed their spatial dynamics using 43-year (1977-2019) winter survey data.We hypothesized that the body size groups of these populations became spatially more segregated over time, which was associated with contracted distribution area of body size groups, and/or elongated distance between centers of abundance of these groups.In addition, these changes were caused by body size-specific responses to environmental stress, including ocean warming and population decline.

Fish populations and survey data
The North Sea is a European epicontinental sea connected to the Atlantic Ocean.
The north part of the North Sea is deeper, colder with higher salinity, while the south part is warmer, shallower with higher primary productivity.The North Sea has experienced rising sea temperature and intensive fishing activities over the past decades (Murgier et al. 2021).Fishing has been more intensive in the south part of the North Sea (Engelhard et al. 2014).
We focused on three fish populations in the North Sea: Atlantic cod, haddock and whiting.They belong to the Gadidae family and are demersal populations which live just above the bottom of the sea (for life histories of three populations see Table S1).They have spawning migration in winter (Tobin et al. 2010, González-Irusta andWright 2016).
Evidence have shown that North Sea Atlantic cod is a metapopulation composed of three subpopulations: South, Northwest, and Viking (ICES 2020).
We obtained the survey data of three target populations from the online database of the International Bottom Trawl Survey (IBTS) of International Council for Exploitation of the Sea (ICES) (https://data.ices.dk/).This survey follows a stratified sampling on survey rectangles of 1° longitude  0.5° latitude.The dataset is in the form of catch per unit effort (CPUE) per body size (in 10mm unit) for each rectangle and year-quarter.We extracted the winter data (January to February) between year 1977 and 2019 as our study period, because fishing gear was not standardized until 1977.We did not analyze the summer data, because the survey period is relatively short (starting from 1991), and that seasonal differences in the spatial structure is out of the scope of our study.

Define body size groups within populations
We examined the spatial dynamics at the body size level.We followed the most common approach for body size grouping through dividing a population into equal body size bins (Barbeaux and Hollowed 2018, Li et al. 2019, Yang et al. 2019).We first summed the CPUE for each body size bin (in 1mm unit) over time and survey rectangle, to derive body size distribution.As the distribution was right-skewed, we removed individuals below 5% and above 85% quantile to avoid extremely low abundance at both ends.Then, we divided the body size distribution into equal-interval body size groups.
We tested different body size group number (10, 15, and 20 groups) to see how it influenced the value of spatial dynamics.While higher size group number gave higher precision, the spatial dynamics did not differ with group number (Table S2).We thus reported the results with 20 body size groups in the main text.We did not use group number higher than 20, otherwise would leads to too few individuals for largest and smallest body size groups; this could raise uncertainty of the results.
Deriving fixed number of body size group for each population leads to wider body size bins for larger populations, and narrower body size bins for smaller populations.To confirm the temporal dynamics of spatial overlap within populations, we alternatively derived body size groups by using fixed bin width for all three populations (e.g., 5 cm).
We also examined the changes in the overlapped area over time between life stages within populations as a preliminary test.To do so, we grouped each population into juvenile and adult, based on the body size at 50% maturity (Table S1).
We did not analyze the spatial structure using age groups because existing agespecific data did not distinguish age groups older than six years.Thus, spatial dynamics calculated using this dataset would neglect the dynamics between older groups.In addition, body size interval differed from one age to another due to non-linear age-size relationships.
Because the results from age group or size group are not comparable, we reported only spatial structure between body size groups in this work.

Spatial structure indices
To explore the temporal changes in the spatial distribution of body size groups within populations, we calculated the following indices for each survey year: 1) area of distribution of each body size group, 2) center of abundance of each body size group, 3) overlapped area between pairs of body size groups, and 4) distance between centers of abundance of pairs of body size groups.There is a total of 190 (C 2 20 ) pairs of body size groups between 20 body size groups within a population.
The area of distribution of each body size group is the proportion of occupied area at any given year, over the maximum occupied area of the same body size group over the study period.This standardized measure accounts for variations in the occupied area between different body size groups.⁄ ).The concepts of distributional overlap has been used in inter-species co-occurrence at the community level (Griffith et al. 2018, Carroll et al. 2019), but not at the body size level.
The distance between centers of abundance is the longitudinal or latitudinal distance between centers of occupied area of a pair of body size groups.For body size group i and j at year t, the distance between centers of abundance in longitude is |  ,, - ,, |, while the distance between centers of abundance in latitude is Atlantic cod has three subpopulations in the North Sea (ICES 2020).Thus, we calculated the spatial structure of Atlantic cod at both the regional scale, as well as at the spatial scale concerning each subpopulation.

Population biomass decline
We used the estimates of yearly total stock biomass from the ICES stock assessment (ICES 2016, 2018) as a proxy for population depletion level (Zhou et al. 2017).Total stock biomass showed a declining trend from 1977 to 2019 for all three populations (Fig. S1).

Ocean warming
We used sea bottom temperature as an indicator of ocean warming, because all three target populations are demersal species.We obtained the sea bottom temperature of the CTD stations across the North Sea region from the ICES online database.We obtained the yearly winter sea bottom temperature at the North Sea region by averaging the measurements from all CTD stations at each year.Sea bottom temperature in the North Sea exhibited a temporal increase from 1977 to 2019 (Fig. S1).

Statistical Analysis
The statistical models were constructed separately for each target population.We applied a four-step analysis are as follows: (1) We used linear mixed-effects models to examine the temporal trends in the overlapped area of pairs of body size groups.Overlapped area of 190 paired groups was included as the response variable (not averaged but as 190 measures).Overlapped area is countbased percentage data.Therefore, it was logit-transformed before model fitting for better homoscedasticity.Survey year was normalized and fitted as a fixed effect.The id of paired groups nested within the survey year was fitted as a random effect.This allows for random intercept and slope for each pair of body size group.We repeated the same analysis to examine the temporal trends in the distance between centers of abundance for pairs of body size groups, without transforming the response variable.Then, we repeated the analysis to examine the temporal trends in the area of distribution of body size groups.Area of distribution is continuous proportional data.Thus, it was logittransformed before fitting.The id of body size group nested within survey year was fitted as random effect.
(2) Then, we constructed a multiple regression model to test the relative importance of area of distribution and distance between centers of abundance on overlapped area.We regressed yearly mean of overlapped area (mean of 190 paired groups) against yearly mean of area of distribution (mean of 20 body size groups), and yearly mean of distance between centers of abundance in latitude and longitude (mean of 190 paired groups).
This resulted in 43 data points (43 years) in each model.To account for serial correlation in time series data, we included the temporal autocorrelation of one-step time lag (AR1).
For the initial model, we included an interaction term between area of distribution and the distance between centers of abundance.As none of the interaction term was significant for neither population, we removed the interaction term from the initial model.The final model wrote: Yearly mean of overlapped area across paired groups ~ 1 Yearly mean of area of distribution across body size groups + 2 Yearly mean of distance between centers of abundance in longitude across paired groups + 3 Yearly mean of distance between centers of abundance in latitude across paired groups + AR1, where  represents the fixed effects coefficients.All explanatory variables were normalized before fitting.All the explanatory variables had variance inflation factors < 6, suggesting no noticeable multicollinearity.We extracted the fixed effects coefficients with 95% confidence intervals to represent the relative importance of each explanatory variable.
(3) For each body size group, we evaluated the temporal trends in the area of distribution and center of abundance.To do so, we fitted a simple linear regression model for each body size group separately.We included the area of distribution (logit-transformed), or center of abundance of a body size group, as the response variable.We included survey year as the explanatory variable.We used the slope coefficient to indicate the rate of change in the area of distribution or center of abundance.Then, we examined how the rate of change varied with body size.To do so, we used nonparametric loess regression models.We included the rate of change in area of distribution, or center of abundance, as the response variable.We included body size group as a continuous explanatory variable.
(4) Finally, we examined whether the overlapped area was influenced by sea bottom temperature (Temperature) and total stock biomass (Biomass).In addition, we examined how the effects differed between body size groups within each fish population.We hypothesized that the overlapped area is shaped by the area of distribution, and center of abundance of each body size group.Therefore, we examined the effects of Temperature and Biomass on these two variables.Temperature and Biomass are highly colinear for three populations.Thus, we tested their effects using separate models.In We performed linear mixed-effects models using the function lmer from the lme4 package.
P-values were extracted using lmerTest package.We extracted Conditional R 2 (variance explained by both fixed and random effects) from the function r.squaredGLMM of MuMIn package.We performed the loess regression model with the geom_smooth function of ggolot2 package.We further used heatmaps to visualize the differences in the temporal trends of overlapped area between each pair of size groups.

Results and Discussion
Temporal decline in overlapped area between body size groups For all three populations between 1977 and 2019, the overlapped area between pairs of body size groups declined; that is, the spatial segregation increased between 20 body size groups (Fig. 1).The declining trends were significant regardless of the number of size groups we defined for each population (from 10 to 20 size groups, see Table S2), or fixed size bin width (e.g., 5cm, Fig S2).In addition, the declining patterns were observed for each subpopulation of Atlantic cod (South, Northwest, and Viking) (Table S3), suggesting a universal declining spatial overlap for the Atlantic cod metapopulation.
For Atlantic cod, the decline in spatial overlap was strong between small groups, between large groups, and between small and large groups (Fig. S3).Supporting these results, we observed clear declines over time in the number of co-occupied survey rectangles between juvenile and adult stages (Fig. S4).In contrast, for haddock, the decline in spatial overlap occurred only between small size groups (Fig. S3).Similarly, whiting showed declining spatial overlap between smaller groups, but increasing spatial overlap between larger groups (Fig. S3).The lack of changes in the spatial overlap between small and large groups, for both haddock and whiting, explained why the changes in co-occupied survey rectangles between juvenile and adult stages are less drastic compared to Atlantic cod (Fig.

Contraction of the area of distribution of body size groups
One mechanism of spatial segregation between body size groups over time was related to the contraction of their distribution area, driven by rising sea temperature and/or population biomass decline.This mechanism was strongest in Atlantic cod and haddock.
The mechanism was weaker in whiting, which exhibited contraction of distribution area for smaller groups but expansion for larger groups over time.
Particularly, for Atlantic cod and haddock, the mean distribution area across body size groups declined over time (Fig. 2a-b).The mean distribution area was positively associated with the mean overlapped area across pairs of body size groups (Fig. 3a-b).In addition, total stock biomass positively contributed to the distribution area of each body size group (Fig. 4a-b).These results suggest that declining total stock biomass over time for these two populations (Fig. S1) led to contracted distribution area of body size groups, which in turn decreased their overlapped area.Particularly, during the latter years with lower stock biomass, larger size groups of Atlantic cod contracted their distribution areas in a greater rate than smaller groups (P < 0.0001 for an interactive term of total stock biomass  body size group, Fig. 4a, Table S4).This pattern implies a greater removal of larger groups under intensive fishing exploitation time period (Horwood et al. 2006, Hsieh et al. 2010a).The positive association between population biomass and the area of distribution of body size groups agrees with the positive relationship observed at the whole population level of many fish species, as a result of density-dependent habitat selection (Fretwell and Lucas 1970, MacCall 1990, Fisher and Frank 2004, Thorson et al. 2016).Our finding is also supported by earlier evidence that during low abundance years, the area of distribution of age-1 and age-2 North Sea cod contracted to less than half of that available, towards habitats that have near-optimal bottom temperatures (Blanchard et al. 2005).
Overfishing is a potential main reason for biomass decline and spatial segregation within the Atlantic cod population.However, we did not examine the direct impact of fishing activity (i.e., fishing mortality) on spatial dynamics of Atlantic cod.It is because Atlantic cod is categorized as overexploited species, and its biomass recovers very slowly even after relaxing the fishing pressure since 1990s (Köster et al. 2014).Thus, instantaneous fishing mortality measured at each year does not reflect the long-term impacts of fishing on the biomass and spatial structure of the population.Thus, in this study, we used estimated total stock biomass as the indicator of population depletion level (Froese et al. 2017) rather than fishing mortality, as a proxy to examine long-term fishing impacts on population spatial dynamics.
In contrast to Atlantic cod and haddock, whiting did not have a significant decline in the mean distribution area across body size groups (Fig. 2c).It was because larger groups expanded their distribution area while smaller groups contracted their distribution area over time (Fig. 4c).However, the mean area of distribution across body size groups was still positively related to their overlapped area (Fig. 3c).
In addition to the effect of population biomass decline, ocean warming also impacted the distribution area of body size groups, and the impacts varied among populations.For haddock, sea bottom temperature negatively explained the area of distribution of all body size groups (slope coefficient  standard error = -0.0120.004,P < 0.005, Fig. 4b, Table S4).That is, the rising temperature over the study period contributed to the contraction of the distribution area of all body size groups, which then reduced the overlapped area between them (Fig. 3b).In contrast, for whiting, rising sea temperature resulted in the contraction of the distribution area of smaller size groups, but expansion of distribution area of large sizes groups (P < 0.01 for the interactive term of sea bottom temperature  body size group, Fig. 4c).The differential responses between smaller and larger groups explains the lack of temporal patterns in the mean area of distribution across body size groups of whiting (Fig. 2c).In contrast to haddock and whiting, the distribution area of body size groups of Atlantic cod was determined by the population biomass but not by the sea bottom temperature (Fig 4a).
We speculate that the differences between haddock and whiting, in their distribution area response to ocean warming, may be due to their prey types.Haddock, regardless of body size, mainly feeds on benthic organisms which are spatially restricted under environmental changes (Schückel et al. 2010).In contrast, whiting is one of the top marine predators feeding on fishes, such as Norway pout, sandeel and sprat (Hislop et al. 1991).
These fish prey have higher dispersal potential than benthic organisms under environmental changes, and thus could lead to the expansion of distribution for adult whiting that followed their prey.This is supported by otolith microchemistry analysis, showing that adult whiting can travel long distances (>500 km) to faraway spawning areas (Tobin et al. 2010).Whereas, contrary to larger size whiting, the distribution area of small size whiting contracted over time (Fig 4c).These observations support the notion that larger groups of some fish populations can be resistant to adverse conditions related to warming, and could have better knowledge and higher mobility moving to the optimal foraging and spawning grounds (Hsieh et al. 2010a).
Haddock and whiting have shifted northward since 1977, and the shift of whiting was correlated with warming (Perry et al. 2005).If some fishes have shifted outside of the North Sea, then the population biomass within the North Sea may reduce, leading to contraction in the distribution area and then spatial segregation between body size groups.
Nevertheless, the spatial overlap indices used in our study is not sensitive to the spatial boundary of populations.This is because the indices are calculated based on the ratio of cooccupied area over occupied area by each body size group.Thus, these indices reveal the temporal variations in the degree of spatial overlap within the region analyzed in this study.

Distance increased between the centers of abundance between body size groups
In addition to the area of distribution, we hypothesized that the overlapped area between body size groups were negatively associated with their distance between the centers of abundance.In addition, such pattern was due to body size-specific shift in the centers of abundance, responding to rising sea temperature or population biomass decline.We found that whiting was the only population that exhibited this mechanism.In contrast, the spatial overlap within Atlantic cod and haddock was mainly determined by the area of distribution of body size groups.
Particularly, whiting showed temporal increases in the mean distance between centers of abundance across pairs of body size groups in longitude and latitude (Fig. 2 f, i).
In addition, the mean distance was negatively associated with the overlapped area across paired groups (Fig. 3c).These results suggest that an increase in the distance contributed to a decline in the overlap between body size groups.The center of abundance of larger whiting shifted westward, while smaller groups shifted eastward (P < 0.005 for an interactive term of sea bottom temperature  body size group id, Table S4, Fig. 4f).Therefore, depending on the original position of distributions, the body size-varying shift in the centers of abundance may have increased the distance between body size groups, hence reducing their overlapped area.
In contrast to whiting, for Atlantic cod and haddock, the distance between the centers of abundance across paired groups did not significantly explain the overlapped area (Fig. 3a-b).However, both populations showed an increase in the distance between centers of abundance (except for Atlantic cod at the longitudinal distance) (Fig 2d-e, 2g-h).These results suggest that the changes in the distance were too weak to influence the overlapped area between body size groups.Instead, the contraction of area of distribution of body size groups was the main driver for the spatial segregation under lower population biomass for Atlantic cod and haddock (Fig 3a-b).Interestingly, for haddock, the center of abundance of larger groups at latitudinal direction was more negatively associated with the sea bottom temperature, compared to smaller haddock (P < 0.01 for an interactive term of sea bottom temperature  body size group id, Fig. 4h, Table S4).Consequently, the centers of abundance of all size groups shifted northward in response to higher temperature, but larger size groups shifted faster than smaller size groups.Such different magnitudes of shift of the centers of abundance of body size groups in response to warming may have led to increased distance between their distributions for haddock.

Implications
While all the populations examined in this study demonstrated increased spatial segregation between body size groups over time, the underlying spatial dynamics of body size groups (i.e., area of distribution and center of abundance) and driving forces (i.e., ocean warming and population biomass decline) differed among the three studied populations (Fig. 5).These results have important implications for exploring the differences between populations in their physiological and biogeographic traits at the body size level.For example, body size groups within a population can exhibit different niches (e.g., thermal tolerance, food requirements (Ciannelli et al. 2013)).What drives different spatial responses among populations depends on the extent to which the niches of body size groups overlap.
For example, populations with stronger or weaker niche preferences between body size groups may respond differently to disturbances such as climatic or anthropogenic stress (Tao et al. 2021).
For Atlantic cod and haddock, the contraction of the distribution area of body size groups was the main driver for the spatial segregation among body size groups over time.
This finding has important implications to identify populations at risk of increased spatial segregation at body size group level.For example, both highly migratory pelagic predators (e.g., tuna, billfish) (Worm and Tittensor 2011) and species living in regional seas (e.g., Monterey Spanish mackerel at the coast of California (Collette and Russo 1984) and yellowtail flounder around Newfoundland (Brodie et al. 1998)) have shown contraction of their distribution area over the past decades.Although the contractions of distribution area were observed at the population level, these patterns may apply to the finer level of body size group.Furthermore, global projections estimated that the biomass of 77% of exploited fishes and invertebrates will decrease when high-temperature extreme will occur (Cheung et al. 2021).These pieces of evidence imply that many fish populations may have exhibited spatial segregation between body size groups, especially for those underwent reduced population biomass and contracted area of distribution, and for those living in climateunstable regions.Large distributional shift may also reduce population biomass and distributional area at the original habitats.For example, in the North Sea, nearly two-third of fish species have shifted northward or deeper between 1977 and 2001 (Perry et al. 2005).
Further investigations on these species which are "on the move" in the North Sea and beyond can help identify the state of the art of spatial dynamics within these populations, and to examine the spatial mechanisms and drivers for these vulnerable populations.These results are helpful to prioritize management and conservation efforts.
The ecological consequences of spatial segregation between body size groups of a population needs further investigation.While population growth may increase due to weakened cannibalism and competition, spatial segregation of size groups may increase the vulnerability of demographic structure to local perturbations.This merits future research to investigate the net effects of within-population spatial segregation on population dynamics and stability.

Conclusion
Recently, increasing evidence on aquatic and terrestrial populations has shown that the shift in spatial distribution varied between life stages or body size under environmental change (Bell et al. 2015, Máliš et al. 2016, Fei et al. 2017, Barbeaux and Hollowed 2018, Frank et al. 2018, Yang et al. 2019).However, it remains unclear to what extent different size groups within populations has segregated from each other over time.We develop a new analytical approach to deepen the understanding of spatial dynamics within populations under global environmental stress.This approach can be applied to populations at various terrestrial and aquatic ecosystems globally, to identify vulnerable populations under environmental stress.This approach also allows us to uncover the mechanisms of spatial segregation within populations, which have profound consequences in demographic connectivity and population stability.
addition, we included AR1 in the model to account for the temporal autocorrelation.The four full models were: i. Yearly area of distribution of each body size group ~ 1 body size group id  yearly Temperature + 2 CPUE + AR1, and ii.Yearly area of distribution of each body size group ~ 1 body size group id  Biomass + 2 CPUE + AR1, and iii.Yearly center of abundance of each body size group in longitude or latitude ~ 1 body size group id  Yearly Temperature + AR1, and iv.Yearly center of abundance of each body size group in longitude or latitude ~ 1 body size group id  Yearly Biomass + AR1, where log-transformed CPUE of each body size group was included as a covariate to account for abundance-distribution relationships.From each full model, we performed a backward stepwise model selection.We derived the most parsimonious model based on AIC and R 2 values.

Fig. 2 Fig. 3
Fig. 2 Area of distribution (a-c) and distance between centers of abundance in longitude (d-f)

Fig. 4
Fig. 4 Rate of change in the area of distribution and center of abundance between 1977 and

Fig. 5
Fig. 5 Conceptual diagram illustrating changes in overlapped area between body size groups Partial overlapped area is proportion of co-occupied area over the occupied area of each body size group of the pair and then taken average.For body size group i and j at time t, partial overlapped area is 0.5 × ( ,,,  , ⁄ + ,,,  , Thus, this measure allows us to directly comparing distribution area between different body size groups.The occupied area of a body size group at a given year is defined as the number of survey rectangles where the CPUE of this group is greater than zero.Therefore, the area of distribution of body size group i at year t is  ,  (  ) ⁄ , where  , is the number of rectangles with the non-zero CPUE of body size group i at year t, and  (  ) is the maximum number of rectangles of body size group i over the study period.Center of abundance is CPUE-weighted center of occupied area for each body size group.For body size group i at year t, the center of abundance in longitude is  ,, = ∑  ,,  =1 ×   ∑  ,,  =1 ⁄ , where   is the longitudinal center of rectangle r, and N is the number of survey rectangles where the CPUE of the whole population is greater than zero.Similarly, the center of abundance in latitude  ,, = ⁄ , where   is the latitudinal center of rectangle r.The overlapped area for a given pair of body size groups are indicated by union overlapped area, and partial overlapped area.Union overlapped area is the proportion of cooccupied area, over the area where either of the body size group occupies.For body size group i and j at year t, the union overlapped area is  ,,,  ,,, ⁄ , where  ,,, is the number of rectangles where both body size group i and j have CPUE greater than zero at year t, and  ,,, is the number of rectangles where either body size group i or j has CPUE greater than zero at year t.