{'fraction_non_alphanumeric': 0.030747812251392205, 'fraction_numerical': 0.009029435163086715, 'mean_word_length': 4.968898385565052, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Regularity underlies erratic population abundances in marine ecosystems

The abundance of a species' population in an ecosystem is rarely stationary, often exhibiting large fluctuations over time. Using historical data on marine species, we show that the year-to-year fluctuations of population growth rate obey a well-defined double-exponential (Laplace) distribution. This striking regularity allows us to devise a stochastic model despite seemingly irregular variations in population abundances. The model identifies the effect of reduced growth at low population density as a key factor missed in current approaches of population variability analysis and without which extinction risks are severely underestimated. The model also allows us to separate the effect of demographic stochasticity and show that single-species growth rates are dominantly determined by stochasticity common to all species. This dominance---and the implications it has for interspecies correlations, including co-extinctions---emphasizes the need of ecosystem-level management approaches to reduce the extinction risk of the individual species themselves.


Introduction
Assessment of extinction risk and biodiversity loss is a central problem in ecology, which has direct implications for ecosystem management practices [1] and policies for the exploitation of natural resources [2]. It is estimated that currently up to 0.1% of all known species go extinct every year, which is over one thousand times above the background extinction rate observed in fossil records [3]. Whether caused by habitat degradation, interspecies competition, climate change, overexploitation, or the introduction of exotic species, the majority of such extinction events have not been anticipated. Existing inventories of the global conservation status, such as the IUCN Red List [4], are believed to include only a fraction of all endangered species and, even for those, significant uncertainty remains on their actual extinction risk. The central difficulty is that wild populations generally do not exist in a steady state from which minute deviations or trends could be detected. Instead, they tend to exhibit fluctuations over time [5][6][7][8][9].
Temporal variations in population abundance may be caused by species interactions, environmental changes, migration patterns, intrinsic nonlinearities, and/or human exploitation [10][11][12], and are sometimes sufficiently irregular to be regarded as stochastic. Such irregular time dependence, when combined with unavoidably imperfect sampling of the population, poses considerable challenges for population forecast. Significant previous research has focused on determining the frequency composition (so called "noise colour") of such fluctuations and their correlations with the environment [13,14]. Despite their immediate implications for sustainable ecosystem exploitation and management, much less understanding has been generated about the factors that influence the growth rate of individual species.
Perhaps not surprisingly, historically there has been a large number of both false positives and false negatives in the assignment of extinction status and extinction risk. For instance, the New Zealand's bird takahē, which was considered extinct by the end of the 19th century, was rediscovered in the wild 50 years later [15]. This is one of now many known examples of a Lazarus taxon [16], in which a sparse population passed undetected for an extended period of time. The passenger pigeon in North America, on the other hand, went from an abundant population to functional extinction in less than 20 years-which then led to actual extinction in the early 20th century [17]. Taken together, the picture that emerges is one in which the survival of a species appears to depend very subtly on the species' own abundance [18]. This picture is further complicated by the possibility of co-extinctions [19,20], whose actual role beyond directly dependent species (such as predator-prey and parasite-host) remains elusive [21][22][23]. Compared to the terrestrial case, extinctions of marine species caused by habitat loss or human exploitation have been rare, in part due to stabilizing effects such as changes in fish catchability at low population densities. Yet, the pace of marine defaunation is likely to accelerate dramatically as the strength and scope of human impact on marine ecosystems grow [24]. This underscores the need for a quantitative understanding of the growth dynamics of marine populations, including extinction risks.

Results
Here, we examine the dependence of growth rate on population abundance and stochastic factors, both when only the temporal variability of individual species is observed and when the dynamics of the entire ecosystem is taken into consideration. We base our analysis on data of the Large Marine Ecosystems (LME) portion of the Sea Around Us Project (http://www.seaaroundus.org) [25], which is a global-scale database on species abundance in 66 marine ecosystems. For each ecosystem, the data consist of the annual quantities (in tonnes) of its 12 most abundant species caught by fisheries over the period 1950 to 2006. Such landing data is a widely used tool for making inferences about marine populations [26][27][28][29], even if such use is occasionally controversial [30], as factors other than population abundance can influence commercial catches.
We will show, however, that our results hold true for marine stock assessments, which integrate research surveys and catch data with independent information (such as mortality rates and population size/age structure) to obtain a more accurate estimate of population abundance (supplementary material, Analysis of Stock Assessments). Therefore, we present our analysis based on the more widely-available landing data [25], followed by supplementary validation using stock assessment data.
For each species i in a given ecosystem, we assume that the annual abundance in year t can be approximated-up to a scaling factor-by the reported landing, denoted by x (i) t , although we will show that our results do not depend critically on this assumption (supplementary material, section 1, figures S1-S3). Our object of study is the year-to-year growth rate, To avoid ill-defined log functions associated with zero landings, we add 1 to each data point, so that the minimum value of x (i) t is 1. Without any further assumptions, the year-to-year change in the species' population abundance can always be written as where r t represents the growth rate in year t decomposed into the averagē r (i) and standard deviation σ (i) of the growth rate (calculated over the entire time series) and a time-dependent factor ξ the marine ecosystems we consider. In contrast with previous studies, here we will make no a priori assumptions on the growth rate, instead deriving its properties directly from the data. t , can be regarded as a random variable drawn from an appropriate distribution that is nearly stationary. The data show that, to an excellent approximation, this distribution is given by a double-exponential function also known as the Laplace distribution (figure 1B-C). This is itself an important, novel finding, which further allows us to devise a statistical model, as discussed below. We have verified that our addition of 1's to the zero-landing years does not affect this distribution, which incidentally also governs the fluctuations of growth rate for the population abundance of each ecosystem as a whole ( figure 1D). In addition, we have rigorously confirmed the fit of the individual species' growth fluctuations to the Laplace distribution versus the null hypothesis of a normal distribution, according to standard goodness-of-fit tests such as the Kolmogorov-Smirnov test and Akaike information criterion. These tests show that the former distribution is a significantly more plausible explanation than the latter distribution for the growth rates of a large majority of species (supplementary material, figure S4).
Having established that the normalised growth rate fluctuations obey a stationary stochastic process described by a parameter-free Laplace distribution, we can propose Eqs. 1-2 themselves as a stochastic model of population abundance. This model is expected to be appropriate away from the floor abundance x (i) t = 1. However, the growth rate at the floor is much smaller than at any other abundance, with probability 89% of being zero, as shown in the inset of figure 1C.
Naturally, because abundance is measured in integer units of landing tonnes, x indicates that the population is very low but not necessarily that it is zero. For this reason (and possibly due to migration and sample biases), the apparent "extinctions" considered here are local in space and in time. Indeed, the populations generally recover to detectable abundances in the systems under consideration. However, they do so more slowly than predicted by the overall growth rates. This reduced growth rate at low population density is consistent with the Allee effect (or depensatory population dynamics), which is a scenario previously observed for a number of species in diverse ecosystems [32]. Further analysis would be needed to conclusively address that particular issue in this case, which falls outside the scope of this paper. It should be noted, nevertheless, that the apparent reduced growth rate at low population abundances is counter to most existing models, including the classical Ricker model [33]. A recent study of exploited marine species stocks included in the RAM Legacy Stock Assessment Database [3] indicates that this missing element is in fact the probable explanation for the slow recovery of depleted stocks when compared with predictions from models commonly used in fisheries management [35].
We incorporate the floor effect into our model by taking the probability distribution of ξ (i) t to be and F 2 (ξ) = 0 otherwise. The parameter p (i) is a measure of the recovery probability once the species reaches the floor abundance, δ is the Dirac delta function, and ξ c (for simplicity in equation (3), we used the approximation c -is typically close to zero). Here, the average growth rater (i) and standard deviation σ (i) are estimated for x   figure 1D). The occurrence of pseudo-extinctions is systematically underestimated when stochasticity is dominantly species-specific, as illustrated in figure 4A for α = 0. In fact, the agreement of the model with empirical data is significantly better when stochasticity is taken to be dominantly common to all species (figure 4B), and the agreement becomes excellent for α = 0.8 ( figure 4A).
This corroborates the conclusion that variations in growth rate are largely synchronised within each ecosystem. We speculate that this synchronisation is partially rooted in external, environmental fluctuations, which are known to play a crucial role in the dynamics of terrestrial populations [11,37]. In marine systems, likely environmental fluctuations include the El Niño Southern Oscillation, the North Atlantic Oscillation, and riverine flood pulses. For fishes, the impacts of these fluctuations can be both direct and indirect-such as via externally-driven changes in the lower food web or via ecosystem regime shifts [38]. Other potential sources of synchronized growth fluctuations include interspecies interactions as well as correlations with human activity-for example shifts in overall fishing effort due to weather or changes in regulations. But regardless of the source of the common stochasticity observed here, the direct implication is an increased risk of the otherwise unlikely concurrent collapse of multiple species.
Our findings should be compared with the case studies of the American breeding bird populations [39] and Hinkley Point's fish community [40], for which growth rates have been analyzed.
For both systems the aggregated non-normalised growth rates were found to follow power-law distributions, but independent analysis of some of these data has shown that after rescaling (by the species standard deviation) to a variable equivalent to the normalised growth rate fluctuation t , the distribution becomes normal [41], which corroborates the conclusion that the growth rate distributions of individual species are short-tailed. This is consistent with the previous analysis of 544 long-term time-series from the global population dynamics database [42], including both aquatic and terrestrial populations, which demonstrated that the abundances of most species are either lognormally distributed or shorter-tailed than lognormally. If one neglects the (significant) year-to-year correlations in abundance, lognormal distributions for individual population abundances imply normal distributions for individual species' growth rates.
The results presented here, on the other hand, show that normalised growth rate fluctuations follow Laplace distributions and this is confirmed to remain true for individual species in the ecosystems we consider (supplementary material, figure S4). Had we not normalised the growth rates to eliminate heterogeneity across species, the resulting aggregated non-normalised growth rate fluctuations would be more fat-tailed than what we observe, without necessarily revealing a simple scaling behaviour (supplementary material, figure S5). Importantly, we have verified that the observed Laplace distributions for the normalised growth rates are not artefacts of the landing data used here as a proxy for population abundance, remaining valid for stock assessment data. This is demonstrated in figure S6 (supplementary material) for the RAM Legacy Stock Assessment Database [3], which is the most complete marine stock assessment catalogue available. The estimates of population abundance therein were acquired under controlled set-tings, using a variety of methodologies designed to avoid systematic biases, and in a variety of ecosystems, which substantiates the conclusion that Laplace statistics underlie the growth of marine populations in general. Note that different types of growth rate distributions are consistent with the lognormally-distributed abundances observed in previous studies; by Eq. 1, the (log) population abundance is (up to a constant) the sum of the growth rates in all previous years. Thus, by the Central Limit Theorem, any growth rate distribution (including Laplace) will eventually lead to normal distributions for the log abundances, provided the growth rates are independent and identically distributed with finite variance. However, direct analysis of the marine population abundances in this study reveals that their distributions are no closer to lognormal than to power laws (supplementary material, figure S7). Altogether, our results are significantly different from those suggested by previous studies, and they do not follow from existing ecological models.
Our demonstration that the growth rates of marine species are governed by the Laplace distribution, which has a heavier tail than a normal distribution with the same standard devi- , I 00 = {t|x t = 1 and x t+1 = 1}, and | · | denotes the number of elements in the set.
Model predictions and validation. For each time series {x t } of length T , we use its first T 0 data points as the training set for parameter estimation and the remaining T − T 0 data points for validation. We calibrate our model both without and with floor effect, and simulate both variants for 1, 000 independent runs. Each run starts with the same initial condition x T 0 and is computed for a total of T − T 0 time steps, with x t reset to 1 whenever its estimated value is below 1. For both the emprical and simulated data, we calculate the pseudo-extinction risk as the probability that x t ≤ θX , where θ ∈ 1 X , 1 is the pseudo-extinction threshold and X = max 1≤t≤T 0 x t . The special case of pseudo-extinction risk at the floor level is defined as the fraction of years for which the population abundance is at the floor (i.e., x t = 1). In this study, T = 57 years and we choose T 0 = 45 years.
Estimation of common stochasticity. We investigate the existence of common stochasticity among all species considered in each ecosystem by extending our model as follows. In a given year t, for those species that do not remain at the floor abundance (according to Eq. 3), we draw a vector ξ t = ξ (i) t of growth fluctuations from a multivariate Laplace distribution [45] characterised by (vector) mean 0, (vector) variance 1, and normalised covariance matrix Γ. We set the diagonal elements of Γ identically to 1 and the off-diagonal elements identically to α 2 , where α is a parameter ranging from 0 to 1. The resulting growth fluctuations for the individual species are distributed according to Eq. 2 but the correlation coefficient is α 2 between ξ (i) t and ξ (j) t for any pair i = j. As such, the parameter α can be interpreted as the proportion of common stochasticity. The value of α that best represents the empirical data is determined in figure 4.
Competing interests. The authors declare that they have no competing interests.  The inset shows the probabilities P (r t = 0| ln x t = 0) and P (r t > 0| ln x t = 0). When the abundance of a species reaches the floor level ln x t = 0, it has a high probability (89%) of staying at that level the following year, which reveals a strong relation between low population and nearly zero growth rate. We note that, while our addition of 1's to the reported populations will tend to deflate the apparent magnitude of the growth rate at small population abundances when calculated using Eq. 1, the stickiness to the floor reported here is a property of the raw data and is not influenced by this transformation. (D) Same as in panel C at the ecosystem level, where the 'population abundance' of an ecosystem is taken to be the weighted average over all of its 12 recorded species, with the weights defined by the inverse of the average abundance of the species. The exponent √ 2 is within the 95% confidence interval of the linear regressions.

Supplementary Material 1 Effects of Non-Constant Catchability
For a given species, the potential relationship between the reported catch (landing) in a given year, x t , and the underlying population abundance, n t , can be expressed in the general form where f is a function that then describes the catchability, or catch per unit population. For example, in the classical Type I fishery model, f (n t ) = qE, where q is catchability per unit effort and E is total effort [1]. Here we focus on more general functions f and study their effects on growth rates and population dynamics. We will make the assumption that catchability approaches a constant as population size grows large, i.e., f (n t ) → a as n t → ∞. Under this assumption, we explore the following three mostly commonly-considered scenarios.
Baseline. This is the simplest scenario, under which the amount caught is directly proportional to the actual population abundance (i.e., f (n t ) = a). Results reported in the main text were obtained under this assumption. In this case, the growth rates calculated from catch data are equal to the actual population growth rates since x t+1 /x t = n t+1 /n t , and hence the particular value of a does not affect the growth rate distribution nor the predictions of pseudo-extinction risk made by our stochastic model.
Hyperdepletion. In this scenario, catchability decreases at low population abundance, which is typically assumed to mean that f (0) = 0 and df /dn t > 0 for all n t ≥ 0 [2]. As a concrete example, we consider the following functional form: where b > 0 is a parameter. Note that we can use Eq. S4 to solve for population in terms of catch in this case as Hyperstability. In this scenario, catchability increases at low population abundance, which is typically taken to mean that f (0) > 0 and df /dn t < 0 for all n t ≥ 0 [2]. As a concrete example, we consider the following functional form: where c > 0 is an additional parameter. From this we obtain that Figure S1 illustrates the relationship between catch and underlying population abundance for the three scenarios above. In our analysis, we use the parameter values a = 0.1, b = 10 4 , and c = 10 2 . These values were chosen based on the typical scales of the catch data in the LME dataset. Note that for the purposes of this section, we will calculate growth rates and train/validate our model based on the population abundances, which are obtained by transforming the catch data in the LME dataset according to one of the three catchability relations defined above. Figure S2 shows the distribution of growth rate fluctuations under the hyperstability and hyperdepletion scenarios. In both scenarios, we see that the distribution is indistinguishable from the distribution of growth rate fluctuations in the baseline scenario, which we have shown to be Laplace ( Fig. 1 and Fig. S4C-D). As such, Laplace statistics should still form the basis of our predictive stochastic model under other conceivable catchability scenarios. But how do the model's predictions of pseudo-extinction risk fare when based on population rather than reported catch? Figure S3 shows the model-predicted pseudo-extinction risk as a function of pseudo-extinction threshold compared to the empirically-observed risk for the hyperstability and hyperdepletion scenarios. As for the baseline scenario of constant catchability ( Fig. 2A), the inclusion of the floor effect significantly improves our model's prediction of pseudo-extinction risk at low pseudo-extinction thresholds. However, in the case of hyperstability, both models (with and without floor effect) systematically underestimate pseudo-extinction risk.

Analysis of Stock Assessments
We use catch (landing) data for the analysis in the main text because it lends itself to highquality statistics, being available for a large number of years and for a large number of species.
Indeed, catch data is the only available indicator of the population abundance of most marine species. Nonetheless, one must be cautious in using reported landings as a direct proxy for population abundance, since there are factors that can affect catches other than changes in underlying population abundance, such as extreme weather events and changes in market demand or fishing effort.
To address the possibility that these artefacts may have affected our results, we have repeated our statistical analysis on the RAM Legacy Stock Assessment Database [3], which is the largest and most up-to-date catalogue of marine stock assessments available. For each assessed stock (comprising a specific species and geographical location), the relevant data consist of yearly estimates of the total biomass. These data integrate multiple independent sources of information beyond catch data, such as species-specific biological information and the results of research surveys, and are consequently regarded as more accurate estimates of population abundance.
To obtain meaningful statistics and facilitate comparison to results in the main text, we focus on those assessments that have at least 30 consecutive years of data and that have total biomass estimates measured in units of tonnes. Out of the 331 assessments in the database, 199 satisfy these criteria.     Assuming, as usual, that the hypothesis is not rejected for p-values larger than 0.05, it follows that the normalised growth rate fluctuations of most species are consistent with an exponential distribution with exponent close to √ 2. This confirms that the scaling for the aggregate data remains valid for individual species. We have verified that this conclusion also holds using the Akaike Information Criterion (AIC), which is another independent test of goodness of fit. According to the AIC, over 70% of species' growth rate distributions more plausibly follow a Laplace distribution compared to the null hypothesis of a normal distribution. The distribution of non-normalised growth rates (panel A) has a heavier tail than the Laplace distribution representing normalised growth rates; this difference is due to heterogeneity across species, which is mainly reflected in the relatively broad distribution of the standard deviations (panel B). Thus, the normalisation introduced in our analysis is an important step to reveal the scaling behaviour identified in this study.  for the abundances, where the standard deviations σ = σ(ln x t ) are determined by maximum likelihood estimation. (E) Distributions of the average log-abundance per species, ln x t , and of the full time-series of the log-abundances, ln x t , aggregated over all species. Panel E indicates that the relative species abundances are mainly determined by differences between the average populations of different species rather than by their time variations. On the other hand, the distributions of the abundances of individual species over time are generally no closer to lognormal functions (inset of panel D) than they are to power-law functions (inset of panel C).