This paper has two interrelated foci: (i) obtaining stable and efficient data-driven closure models by using a multivariate time series of partial observations from a large-dimensional system; and (ii) comparing these closure models with the optimal closures predicted by the Mori–Zwanzig (MZ) formalism of statistical physics. Multilayer stochastic models (MSMs) are introduced as both a generalization and a time-continuous limit of existing multilevel, regression-based approaches to closure in a data-driven setting; these approaches include empirical model reduction (EMR), as well as more recent multi-layer modeling. It is shown that the multilayer structure of MSMs can provide a natural Markov approximation to the generalized Langevin equation (GLE) of the MZ formalism. A simple correlation-based stopping criterion for an EMR–MSM model is derived to assess how well it approximates the GLE solution. Sufficient conditions are derived on the structure of the nonlinear cross-interactions between the constitutive layers of a given MSM to guarantee the existence of a global random attractor. This existence ensures that no blow-up can occur for a broad class of MSM applications, a class that includes non-polynomial predictors and nonlinearities that do not necessarily preserve quadratic energy invariants. The EMR–MSM methodology is first applied to a conceptual, nonlinear, stochastic climate model of coupled slow and fast variables, in which only slow variables are observed. It is shown that the resulting closure model with energy-conserving nonlinearities efficiently captures the main statistical features of the slow variables, even when there is no formal scale separation and the fast variables are quite energetic. Second, an MSM is shown to successfully reproduce the statistics of a partially observed, generalized Lotka–Volterra model of population dynamics in its chaotic regime. The challenges here include the rarity of strange attractors in the model’s parameter space and the existence of multiple attractor basins with fractal boundaries. The positivity constraint on the solutions’ components replaces here the quadratic-energy–preserving constraint of fluid-flow problems and it successfully prevents blow-up.

# Statistical methods

Data-driven non-Markovian closure models. Physica D: Nonlinear Phenomena. 2015;297 :33–55.Abstract

.
Interannual and interdecadal oscillation patterns in sea level. Climate Dynamics. 1995;11 (5) :255–278.Abstract

.
Weather regime prediction using statistical learning. Journal of the Atmospheric Sciences. 2007;64 (5) :1619–1635.

.
Spectral methods: What they can and cannot do for climatic time series. In: Decadal Climate Variability: Dynamics and Predictability. Springer-Verlag, Berlin/Heidelberg ; 1996. pp. 446–482.

.
An empirical stochastic model of sea-surface temperatures and surface winds over the Southern Ocean. Ocean Science [Internet]. 2011;7 (6) :755–770. Publisher's VersionAbstract

.
North Atlantic climate variability in coupled models and data. Nonlinear Processes in Geophysics. 2008;15 :13–24.

.
Spatio-temporal filling of missing points in geophysical data sets. Nonlinear Processes in Geophysics. 2006;13 (2) :151–159.Abstract

.