Significance The conventional parametric approach to modeling relies on hypothesized equations to approximate mechanistic processes. Although there are known limitations in using an assumed set of equations, parametric models remain widely used to test for interactions, make predictions, and guide management decisions. Here, we show that these objectives are better addressed using an alternative equation-free approach, empirical dynamic modeling (EDM). Applied to Fraser River sockeye salmon, EDM models (i) recover the mechanistic relationship between the environment and population biology that fisheries models dismiss as insignificant, (ii) produce significantly better forecasts compared with contemporary fisheries models, and (iii) explicitly link control parameters (spawning abundance) and ecosystem objectives (future recruitment), producing models that are suitable for current management frameworks.

Abstract It is well known that current equilibrium-based models fall short as predictive descriptions of natural ecosystems, and particularly of fisheries systems that exhibit nonlinear dynamics. For example, model parameters assumed to be fixed constants may actually vary in time, models may fit well to existing data but lack out-of-sample predictive skill, and key driving variables may be misidentified due to transient (mirage) correlations that are common in nonlinear systems. With these frailties, it is somewhat surprising that static equilibrium models continue to be widely used. Here, we examine empirical dynamic modeling (EDM) as an alternative to imposed model equations and that accommodates both nonequilibrium dynamics and nonlinearity. Using time series from nine stocks of sockeye salmon (Oncorhynchus nerka) from the Fraser River system in British Columbia, Canada, we perform, for the the first time to our knowledge, real-data comparison of contemporary fisheries models with equivalent EDM formulations that explicitly use spawning stock and environmental variables to forecast recruitment. We find that EDM models produce more accurate and precise forecasts, and unlike extensions of the classic Ricker spawner–recruit equation, they show significant improvements when environmental factors are included. Our analysis demonstrates the strategic utility of EDM for incorporating environmental influences into fisheries forecasts and, more generally, for providing insight into how environmental factors can operate in forecast models, thus paving the way for equation-free mechanistic forecasting to be applied in management contexts.

One of the fundamental challenges of environmental science is to understand and predict the behavior of complex natural ecosystems. This task can be especially difficult when multiple drivers (e.g., species interactions, environmental influences) interact in a nonlinear state-dependent way to produce dynamics that appear to be erratic and nonstationary (1). In the standard parametric approach, which implicitly assumes that the selected model and its equations are essentially correct, the equations (really just mechanistic hypotheses) can lack the flexibility to describe the nonlinear dynamics that occur in nature. Consequently, these parametric models tend to perform poorly as descriptions of reality, with little explanatory or predictive power (2, 3), and limited usefulness for prediction and management.

Parametric Models as Hypotheses A common problem when applying the parametric approach to nonlinear systems is that of ephemeral fitting. That is, although population models may assume that demographic parameters such as growth rate or carrying capacity are fixed constants, these quantities are often observed to vary in time or in relation to other variables (e.g., resource availability, changing climate regimes) when tested on actual data (4). This principle is illustrated in Fig. 1A, where the Ricker spawner–recruit relationship is fit to the early (1948–1976) and late (1977–2005) halves of the time series from the Seymour stock. Very different relationships emerge in these two time periods, conflicting with the assumption of a fixed equilibrium and constant parameter values. Indeed, Beamish et al. (5) found that the Ricker model fit better when constrained by climate regimes, suggesting that the spawner–recruit relationship does vary in time, a fact consistent with the general notion of nonlinear state dependence (6, 7). Fig. 1. Model output for the Ricker, extended Ricker, and multivariate EDM models. (A) Ricker curves for the Seymour stock of Fraser River sockeye salmon are quite different for the early (blue; 1948–1976 brood years) vs. later (red; 1977–2005) time segments (triangles are observed data). Even fit to the whole time series (gray line), large errors remain. (B and C) Model output (surfaces; points are observed data) from the extended Ricker model (B) and EDM (C) using spawner abundance and Pine Island April SST to forecast recruitment of Seymour sockeye salmon. Although the Ricker model varies smoothly, it can forecast recruitment to be many times higher than the historical maximum. In the EDM model, however, the relationship between temperature and spawners is defined empirically by the data and, thus, more realistically depicted. At its core, nonlinear dynamics [which are known to be ubiquitous in marine species (8, 9)] occur when variables have interdependent effects; this can be problematic when applying a reductionist approach to understand nonlinear systems. For example, in laboratory experiments, guppies (Poecilia reticulatus) preferentially eat Drosophila or tubificid worms depending on which prey is more abundant (10). Thus, the strength of predation on, say, Drosophila, will change depending on the abundance of tubificid worms. This prey-switching behavior typifies nonlinear state dependence, whereby different components cannot be treated independently, as would be the case in a linear system or even a nonlinear system approximated at equilibrium. Consequently, applying a model that assumes separability of effects [e.g., vector autoregression (11)] to a system that is actually nonlinear can give the appearance of nonstationarity or stochasticity even when the underlying mechanisms are unchanged and deterministic. Nonlinearity is also known to affect the correct identification of causal drivers—a key prerequisite for understanding and predicting system behavior. In nonlinear systems, because interacting variables can exhibit transient (mirage) correlations that change in magnitude or sign (6, 7), the use of correlation to identify causal environmental variables can be misleading, producing both false positives (i.e., correlation does not imply causation) and false negatives (i.e., lack of correlation does not imply a lack of causation). Given the prevalence of nonlinear interactions in ecology, mirage correlations can be misleading. Indeed, a metaanalysis examining the robustness of correlations between recruitment and the environment (12) found that only 28 out of 74 initially significant correlations were upheld when subsequent data were included. Even when causal variables are known, their inclusion into improperly formulated models can produce conflicting results. For example, with sockeye salmon in the Fraser River, although anomalous oceanic conditions experienced by juveniles are thought to be responsible for the low abundance of returning adults in 2009 (13⇓–15), extensions of the standard Ricker model that explicitly include environmental factors surprisingly show no significant improvements in the actual forecasts (16⇓–18). A simple explanation for this apparent contradiction is that the extended Ricker model does not accurately portray the relevant interaction between the oceanic environment and sockeye salmon. Indeed, the model naively assumes that the environment acts on recruitment dynamics independently with a constant multiplicative effect (e.g., a 1 °C decrease in temperature always doubles recruitment regardless of other factors important to the state of the system). Although temperature, in all likelihood, does affect recruitment, it probably does not follow this arbitrary form. We demonstrate this by fitting the model to Pine Island sea surface temperature (SST) and Seymour spawner–recruit data (Fig. 1B), finding that the model predicts unrealistically high recruitment (much higher than the historically observed maximum) for hypothetical (but plausible) conditions of high spawner abundance and low temperature. Thus, although the equation may appear reasonable as a hypothesis, it apparently does not incorporate the environment realistically.

Empirical Dynamic Modeling In contrast to fitting an assumed set of equations, empirical dynamic modeling (EDM) instead relies on time series data to reveal the dynamic relationships among variables as they occur (1, 6, 19⇓–21). By extracting these relationships empirically, EDM accommodates potentially complex and changing interactions that cannot be described in a simple set of equations. Thus, prediction accuracy with EDM is constrained by the quantity and quality of data rather than by the hypotheses represented in a set of equations [which may be subject to process error due to false or incomplete specification (22)]. Fundamental to EDM is the concept of a time series as an observation on a dynamic system. Broadly speaking, a dynamic system can be viewed as a set of “states” (d-dimensional vectors where each coordinate is a system variable) and deterministic rules (governing dynamics) for how the states evolve over time. Collectively, the set of states and their trajectories forms an “attractor manifold,” and projecting the motion on this manifold to a coordinate axis produces a time series of the corresponding variable (SI Appendix, Fig. S1A). For example, in a simple predator–prey system where the system evolves as a function of the two abundances, the system state could be represented as the ordered pair of predator and prey abundances. This system state can be projected onto the prey coordinate axis to produce a time series of prey abundance, although many other observation functions are also possible (e.g., predator abundance, average number of prey for each predator). In theory, with time series for all of the system variables, it would be possible to reconstruct the original attractor manifold by plotting each time series as a separate coordinate. In practice, however, we typically do not have these data or know the identity of all relevant variables. Fortunately, a fundamental mathematical result proves that information about the entire system is contained in any one variable (23, 24), meaning that a shadow version of the original attractor can be constructed from just a single time series. This is accomplished by substituting lags of that time series for the unknown or unobserved variables (SI Appendix, Fig. S1B). These essential mechanics of EDM are detailed in ref. 6 and crisply summarized in a short animation (Movie S1). Although a single time series is usually sufficient to reconstruct a system’s dynamics, there are exceptions (e.g., it is not a closed system). In the case of sockeye salmon, abundance alone may not skillfully predict future returns because they are influenced by external environmental factors. Here, the environment may be thought to act as stochastic external forcing, necessitating its inclusion as an additional coordinate in a multivariate reconstruction (1, 7, 24). We demonstrate this by using spawners and SST to predict recruitment (Fig. 1C). Unlike a parametric model in which a hypothesized interaction must be specified in advance (the extended Ricker model; Fig. 1B), the empirical surface in Fig. 1C makes no assumptions about the relationship between variables, but instead captures the interaction between density dependence and environmental conditions as revealed by the data: ocean temperatures have a stronger effect on recruitment when spawner abundance is low.

Fraser River Sockeye Salmon In this work, we perform a real-world test comparing EDM and the standard parametric paradigm, by forecasting returns for the nine most historically abundant stocks of sockeye salmon from the Fraser River system in British Columbia, Canada (Fig. 2), of significance to Canada’s iconic fisheries. Total returns in this system are highly variable and can span over an order of magnitude: a record low of 1.6 million in 2009 was followed by a record high of 28.3 million in 2010 (Fig. 3). Although some of this variability occurs because of cyclic dominance (25, 26), large interannual fluctuations in mortality and productivity (recruits-per-spawner) are difficult to predict, leading to considerable uncertainties in current parametric forecast models (27). This is suggestive of nonlinear dynamics in this fishery, and indeed, a Canadian federal inquiry (13, 14) concluded that recent declines in productivity could not be attributed to any single mechanism but were likely caused by the interaction of multiple stressors (e.g., predators, food availability, environment). Applying a simple S-map test (P = 0.002) (SI Appendix, Fig. S2), we confirm the presence of nonlinear dynamics among returns of Fraser River sockeye salmon. Fig. 2. Combined returns of Fraser River sockeye salmon. Total returns (Dataset S1) for Fraser River sockeye salmon combined across stocks (1954 cycle line in black). Although not all stocks exhibit cyclic dominance, and those that do are not synchronized, cycles are still visible in the aggregated returns. Fig. 3. Early ocean environment for Fraser River sockeye salmon. Upon exiting the Fraser River, juvenile sockeye salmon migrate north through the Strait of Georgia, spending up to a month moving through this ecosystem (31), before continuing through Queen Charlotte Strait and into Queen Charlotte Sound. Red labels for the nine stocks studied in this work are located at the approximate spawning sites. Blue triangles denote the locations of the two lighthouses where SST is recorded. Image courtesy of DFO. Thus, we apply EDM methods to unravel the mechanisms by which the environment may affect sockeye salmon recruitment. First, we compare the classical Ricker spawner–recruit model with equivalent EDM spawner–recruit models. With nearly all adults returning as age 4 or age 5 fish, we can consider the total returns in a single calendar year to be composed of age 4 and age 5 recruits from different spawning broods. Following ref. 16, we predict annual returns by first estimating total recruitment for each spawning brood year. This recruitment is then partitioned by age, and the age 4 and age 5 estimates from separate brood years combined appropriately to forecast returns (Materials and Methods). Note that the time series of spawning abundance and recruitment already account for the effects of the fishery (this information is contained within the time series; Materials and Methods), which enables us to focus on just the natural population dynamics. Second, to investigate the causal influence of the oceanic environment, we consider forecasts produced by the extended Ricker model and equivalent multivariate EDM formulations. In both cases, if the inclusion of environmental variables significantly improves forecasts (Materials and Methods), those variables are taken to have a causal influence on salmon recruitment. Last, to avoid arbitrary fitting and to obtain a robust measure of forecast skill, we apply a fourfold cross-validation scheme for each model: the model is fit to three-fourths of the data to predict the remaining one-fourth out-of-sample, and the procedure is repeated for each one-fourth segment of the time series.

Results Comparison of Spawner–Recruit Forecast Models. As a fair comparison with the standard Ricker model where spawner abundance is used to predict recruitment, we examine an equivalent EDM spawner–recruit model, but which actually has fewer fitted parameters (Materials and Methods). Fig. 4 shows that this simple EDM model has significantly higher accuracy (ρ, correlation between observations and predictions) than the Ricker model, with more accurate forecasts in eight of nine cases and significantly lower error overall [mean absolute error (MAE); SI Appendix, Fig. S3]. Nonetheless, predictions for several stocks (Birkenhead, Chilko, Stellako, and Weaver) are not very skillful (ρ < 0.3), suggesting that in these cases, there is no simple spawner–recruit relationship (parametric or otherwise). Instead, environmental factors (e.g., SST, food availability) may dominate, and better performance can be obtained by accounting for these external drivers. Fig. 4. Comparison of forecast accuracy. Comparisons between equivalent EDM and Ricker models show better forecast accuracy for the EDM models [simple EDM vs. Ricker, t (492) = 1.77, P = 0.039; multivariate EDM vs. extended Ricker, t (492) = 2.20, P = 0.014]. Additionally, including environmental data significantly improves accuracy for EDM [t (492) = 2.83, P = 0.0024], but not for the Ricker models [t (492) = 1.26, P = 0.10]. Incorporating Environmental Influences. As in the actual forecast models (16), we further consider three environmental variables [the Pacific Decadal Oscillation (PDO), SST, and Fraser River discharge] observed at different times and locations (12 time series in total). Each of these factors is believed to have a potential effect on recruitment, although significance has yet to be demonstrated in practice. For each stock, we compare the relative performance of the extended Ricker and corresponding multivariate EDM models that incorporate these environmental variables (Table 1; Materials and Methods). Fig. 4 shows that multivariate EDM is consistently and significantly better at forecasting than the extended Ricker model for all nine stocks, and is true for both accuracy and precision metrics (SI Appendix, Fig. S3). Here, the relevant causal influence of these environmental variables is verified by the fact that multivariate EDM models that include them perform significantly better than their simple EDM spawner–recruit counterparts. Table 1. Forecast skill of models incorporating the environment By contrast, the extended Ricker models show no significant improvement over the simple Ricker models in any of the stocks. The difference between EDM and Ricker is especially visible for Late Stuart, Quesnel, Stellako, and Weaver, indicating that these particular environmental factors (currently considered in assessments) can explain much of the variability in these stocks, provided they are incorporated reasonably (i.e., with the minimal assumptions of EDM). For Birkenhead and Chilko, however, multivariate EDM models performed no better than the simplified stock–recruit versions, hinting that variables other than these are required to understand the dynamics of those stocks.

Materials and Methods Data. We analyze yearly time series data for the nine historically most abundant stocks (Birkenhead, Chilko, Early Stuart, Late Shuswap, Late Stuart, Quesnel, Seymour, Stellako, and Weaver) of sockeye salmon from the Fraser River system (Dataset S2). Data span brood years 1948–2005, except for Late Stuart and Weaver, where data begin in 1949 and 1966, respectively. We consider only single-stock models, so notation and equations are given as for a single stock. S t is the number of effective female spawners in brood year t, and R t is the corresponding recruitment (returning adults). Recruitment is partitioned by age: R a,t is the number spawned in year t and returning at age a in year t + a. Following ref. 16, total recruitment is the sum of age 4 and age 5 recruits: R t = R 4,t + R 5,t . In contrast, total returns, N y , are the adults that return to spawn in calendar year y, and computed as N y = R 4,y-4 + R 5,y-5 . As explained below, recruitment is forecast from spawner abundance, and age 4 and age 5 recruits (from different brood years) are summed to estimate total returns in a given calendar year. Note that both recruitment and returns are computed as catch plus escapement plus en route loss, whereas spawner abundance is based on observations of escapement and egg production (27). Thus, both spawner abundance and recruitment account for the effects of catch, and the models we consider here focus just on the population dynamics of this system. We investigate three environmental variables: the PDO, SST, and Fraser River discharge (Dataset S3). For the PDO, one annual time series is constructed as the average of monthly values from November to March (32). SST measures are monthly averages from two lighthouse stations (Entrance Island: April to June; Pine Island: April to July). River discharge is measured at Hope; we include peak daily flow and monthly averages (April to June). Fraser River sockeye salmon enter the ocean at age 2, so the environmental data are lagged 2 years to line up with ocean entry time. Attractor Reconstruction. The goal of attractor reconstruction is to approximate the originating dynamic system using time series data. The simplest construction uses successive lags of a single time series (23, 44): given time series {x t }, E-dimensional vectors x t are composed of E lags of x, each separated by a time step τ : x t = 〈 x t , x t − τ , … , x t − ( E − 1 ) τ 〉 . Generalizations of Takens’ theorem (24, 45) permit attractor reconstructions using multiple time series. For example, with {x t } and {y t } observed from the same system, one possible reconstruction forms vectors as 〈 x t , y t , y t − τ 〉 . To account for different scaling between variables, each time series is first linearly transformed to have mean = 0 and variance = 1. Simplex Projection and S-Map. Simplex projection estimates the trajectory (i.e., forecasts) of a novel system state by computing a weighted average of the trajectories of that state’s nearest neighbors (19). Given an attractor reconstruction, and a novel state x s , we first find the b nearest neighbors (typically setting b = E + 1) that are closest to x s : these neighbors are the vectors x n(s,i) , where n(s,i) designates the time index of the ith closest neighbor to x s . So, x n(s,1) is the closest neighbor to x s , x n(s,2) is the second closest neighbor, etc. We then evolve the neighbors forward, and compute a weighted average of the forward evolutions (h time steps into the future) to estimate x s+h : x ^ s + h = ( ∑ i = 1 b w i ( s ) x n ( s , i ) + h ) / ∑ i = 1 b w i ( s ) . [1] The weights, w i (s), are based on the distance between x s and its ith neighbor, x n(s,i) , scaled to the distance to the nearest neighbor: w i ( s ) = exp ( − d ( x s , x n ( s , i ) ) / d ( x s , x n ( s , 1 ) ) ) , and d(x s , x t ) is the Euclidean distance between the vectors x s and x t . In most cases, we desire forecasts of a scalar value rather than of the full system state. This is possible when the variable to be forecast, y, is an observation on the same dynamic system. As such, there will be a correspondence between x t and the scalar value of y t , and we can adjust Eq. 1 to compute a weighted average of the corresponding values of y: y ^ s + h = ( ∑ i = 1 b w i ( s ) y n ( s , i ) + h ) / ∑ i = 1 b w i ( s ) . [2] The S-map procedure computes a local linear map between lagged-coordinate vectors and a target variable and is often used to test for nonlinear state dependence (22). It includes a tuning parameter, θ, that controls the weights associated with individual vectors: θ = 0 reduces the S-map to a linear autoregressive model of order E, whereas θ > 0 gives more weight to nearby states when computing the local linear map, thus allowing for nonlinear behavior. Following refs. 9 and 46, we test for nonlinearity by computing the decrease in forecast error (MAE) as θ is tuned to be greater than 0 (see SI Appendix for details). Model Descriptions. We formulate EDM models to forecast recruitment from spawner abundance, combining age 4 and age 5 recruits (from different brood years) to estimate total returns in a given calendar year. Acknowledging the persistent 4-year quasicycle, the time series of recruits and spawners are scaled so that each cycle line has mean 0 and variance 1: S t ′ = ( S t − μ k ( S ) ) / σ k ( S ) and R a , t ′ = ( R a , t − μ k ( R a ) ) / σ k ( R a ) , where k = 1, 2, 3, or 4, depending on cycle line and can be computed as k = 1 + ((t − 1) mod 4). σ k and μ k are the mean and SD, respectively, for the kth cycle line. The simple EDM model approximates the system state with 1 lag of the transformed spawner abundance: x t = 〈 S t ′ 〉 . [3] Forecasts of the age 4 and age 5 recruits, R 4 , t ′ and R 5 , t ′ , are made using simplex projection. Here, the two nearest neighbors of S t ′ are identified, and the corresponding values of R 4 , t ′ (or R 5 , t ′ ) are combined in a weighted average to produce a forecast. These forecasts are transformed back into raw values, R ^ a . t = R ^ a , t ′ ⋅ σ k ( R a ) + μ k ( R a ) , and age 4 and age 5 recruits are combined to produce a forecast of total returns, N ^ y = R ^ 4 , y − 4 + R ^ 5 , y − 5 . The multivariate models combine spawner data with up to two environmental indicators: x t = 〈 S t ′ , U t + 2 ′ 〉 x t = 〈 S t ′ , U t + 2 ′ , U t + 2 ″ 〉 , [4] where U t + 2 ′ (or U t + 2 ″ ) is one of the environmental time series described previously, normalized to have mean = 0 and variance = 1. Just as in the simple EDM model, forecasts of the age 4 and age 5 recruits are made using simplex projection and combined to produce a forecast of total returns. However, because including environmental variables increases the embedding dimension, three nearest neighbors are used for models that include one environmental coordinate, and four nearest neighbors for models that include two environmental coordinates. Following ref. 16, we use the standard Ricker model to estimate total recruitment and then partition it into age 4 and age 5 fish. Age 4 and age 5 fish from separate brood years are combined to forecast the number of returns: R ^ t = S t ⁡ exp ( α − β S t ) , N ^ y = R ^ y − 4 ⋅ p 4 + R ^ y − 5 ⋅ ( 1 − p 4 ) , [5] where p 4 is the average fraction of recruits that return as age 4 fish. The extended Ricker model is similar but includes an additional term in the exponent for an environmental covariate: R ^ t = S t exp ( α − β S t + γ U t + 2 ) . [6] Fitting Procedure and Performance Measures. To avoid testing all combinations of (and overfitting) the environmental variables in the EDM model, we sequentially add the environmental variable that most improves forecast accuracy (ρ, the correlation between observed and predicted values). If none of the variables improves forecasts when added, then no further environmental variables are included. Thus, the best EDM model for some stocks may have only 0 or 1 environmental variable (Table 1). Similarly, for the extended Ricker model, we choose the environmental variable that gives the highest ρ. The Ricker models were fit using R 3.0.2 (www.r-project.org/), the Rjags package (cran.r-project.org/web/packages/rjags/index.html), and JAGS 3.2.0 (Just Another Gibbs Sampler; mcmc-jags.sourceforge.net/) following the procedure outlined in ref. 16. Medians of the posterior distribution are used to obtain point estimates suitable for comparison. The EDM models were constructed using R 3.0.2 and the rEDM package (https://github.com/ha0ye/rEDM). The package can be installed with the following lines of R code: > library(devtools)

> install_github(“ha0ye/rEDM”) R scripts for the models can be found in Dataset S4. Data files can be found in Datasets S1–S3. All forecasts are made using a fourfold cross-validation procedure. To quantify model performance, we use Pearson’s correlation coefficient (ρ) between observed and predicted returns as a measure of accuracy and MAE as a measure of error. Comparisons of ρ between models uses a one-sided t test with SE calculated using the HC4 estimator from (47) and with adjusted degrees of freedom as suggested by ref. 48. Improvement in MAE is computed using a one-sided paired t test for the difference, treating each forecast as an independent sample. To compute an aggregate statistic combining all nine stocks, we first scale the observations and predictions for each stock so that the observed returns have mean = 0, variance = 1, and then combine the normalized values across stocks. The comparisons of ρ and MAE are done using this combined set of observations and predictions (Fig. 4 and SI Appendix, Fig. S3).

Acknowledgments We thank Michael Fogarty, Terry Beacham, Brad Werner, Ethan Deyle, Charles Perretti, and two anonymous reviewers for their feedback on this work. This research is supported by a National Science Foundation Grant DEB-1020372 (to G.S. and H.Y.), Foundation for the Advancement of Outstanding Scholarship and Ministry of Science and Technology of Taiwan (C.-h.H.), National Science Foundation–National Oceanic and Atmospheric Administration Comparative Analysis of Marine Ecosystem Organization (CAMEO) Program Grant NA08OAR4320894/CAMEO (to G.S.), a National Science Foundation Graduate Research Fellowship (to H.Y.), the Sugihara Family Trust (G.S.), the Deutsche Bank–Jameson Complexity Studies Fund (G.S.), and the McQuown Chair in Natural Science (G.S.).

Footnotes Author contributions: H.Y., S.M.G., C.-h.H., and G.S. designed research; H.Y. performed research; H.Y. analyzed data; and H.Y., R.J.B., S.M.G., S.C.H.G., C.-h.H., L.J.R., J.T.S., and G.S. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission.

See Commentary on page 3856.

This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1417063112/-/DCSupplemental.