var model research paper

• R(eal) Basics
• Econometric Methods
• Time Series Topics
• Network Analysis
• Reproduction Projects
• About and Disclaimer

An Introduction to Vector Autoregression (VAR)

Since the seminal paper of Sims (1980) vector autoregressive models have become a key instrument in macroeconomic research. This post presents the basic concept of VAR analysis and guides through the estimation procedure of a simple model. When I started my undergraduate program in economics I occasionally encountered the abbreviation VAR in some macro papers. I was fascinated by those waves in the boxes titled impulse responses and wondered how difficult it would be to do such reseach on my own. I was motivated, but my first attempts were admittedly embarrassing. It took me quite a long time to figure out which kind of data can be analysed, how to estimate a VAR model and how to obtain meaningful impulse responses. Today, I think that there is nothing fancy about VAR models at all once you keep in mind some points.

Univariate autoregression

VAR stands for vector autoregression . To understand what this means, let us first look at a simple univariate (i.e. only one dependent or endogenous variable) autoregressive (AR) model of the form \(y_{t} = a_1 y_{t-1} + e_t\) . In this model the current value of variable \(y\) depends on its own first lag, where \(a_1\) denotes its parameter coefficient and the subscript refers to its lag. Since the model contains only one lagged value the model is called autoregressive of order one, short AR(1), but you can easily increase the order to p by adding more lags, which results in an AR(p). The error term \(e_t\) is assumed to be normally distributed with mean zero and variance \(\sigma^2\) .

Stationarity

Before you estimate such a model you should always check if the time series you analyse are stationary, i.e. their means and variances are constant over time and do not show any trending behaviour. This is a very important issue and every good textbook on time series analysis treats it quite – maybe too – intensively. A central problem when you estimate models with non-stationary data is, that you will get improper test statistics, which might lead you to choose the wrong model.

There is a series of statistical tests like the Dickey-Fuller, KPSS, or the Phillips-Perron test to check whether a series is stationary. Another very common practise is to plot a series and check if it moves around a constant mean value, i.e. a horizontal line. If this is the case, it is likely to be stationary. Both statistical and visual tests have their drawbacks and you should always be careful with those approaches, but they are an important part of every time series analysis. Additionally, you might want to check what the economic literature has to say about the stationarity of particular time series like, e.g., GDP, interest rates or inflation. This approach is particularly useful if you want to determine whether a series trend or difference stationary , which must be treated a bit differently.

At this point it should be mentioned that even if two time series are not stationary, a special combination of them can still be stationary. This phenomenon is called cointegration and so-called (vector) error correction models (VECM) can be used to analyse it. For example, this approach can significantly improve the results of an analysis of variables with known equilibrium relationships.

Autoregressive distributed lag models

Regressing a macroeconomic variable solely on its own lags like in an AR(p) model might be a quite restrictive approach. Usually, it is more appropriate to assume that there are further factors that drive a process. This idea is captured by models which contain lagged values of the dependent variable as well as contemporaneous and lagged values of other, i.e. exogenous, variables. Again, these exogenous variables should be stationary. For an endogenous variable \(y_{t}\) and an exogenous variable \(x_{t}\) such an autoregressive distributed lag , or ADL, model can be written as

\[y_{t} = a_1 y_{t-1} + b_0 x_{t}+ b_{1} x_{t-1} + e_t.\]

In this ADL(1,1) model \(a_1\) and \(e_t\) are definded as above and \(b_0\) and \(b_1\) are the coefficients of the contemporaneous and lagged value of the exogenous variable, respectively.

The forecasting performance of such an ADL model is likely to be better than for a simple AR model. However, what if the exogenous variable depends on lagged values of the endogenous variable too? This would mean that \(x_{t}\) is endogenous too and there is further space to improve our forecasts.

Vector autoregressive models

At this point the VAR approach comes in. A simple VAR model can be written as

\[\begin{pmatrix} y_{1t} \\ y_{2t} \end{pmatrix} = \begin{bmatrix} a_{11} & a_{12} \\ a_{21} & a_{22} \end{bmatrix} \begin{pmatrix} y_{1t-1} \\ y_{2t-1} \end{pmatrix} + \begin{pmatrix} \epsilon_{1t} \\ \epsilon_{2t} \end{pmatrix}\]

or, more compactly,

\[ y_t = A_1 y_{t-1} + \epsilon_t,\] where \(y_t = \begin{pmatrix} y_{1t} \\ y_{2t} \end{pmatrix}\) , \(A_1= \begin{bmatrix} a_{11} & a_{12} \\ a_{21} & a_{22} \end{bmatrix}\) and \(\epsilon_t = \begin{pmatrix} \epsilon_{1t} \\ \epsilon_{2t} \end{pmatrix}.\)

Note: Yes, you should familiarise yourself with some (basic) matrix algebra (addition, subtraction, multiplication, transposition, inversion and the determinant), if you want to work with VARs.

Basically, such a model implies that everything depends on everything . But as can be seen from this formulation, each row can be written as a separate equation, so that \(y_{1t} = a_{11} y_{1t-1} + a_{12} y_{2t-1} + \epsilon_{1t}\) and \(y_{2t} = a_{21} y_{1t-1} + a_{22} y_{2t-1} + \epsilon_{2t}\) . Hence, the VAR model can be rewritten as a series of individual ADL models as described above. In fact, it is possible to estimate VAR models by estimating each equation separately.

Looking a bit closer at the single equations you will notice, that there appear no contemporaneous values on the right-hand side ( rhs ) of the VAR model. However, information about contemporaneous relations can be found in the so-called variance-covariance matrix \(\Sigma\) . It contains the variances of the endogenous variable on its diagonal elements and covariances of the errors on the off-diagonal elements. The covariances contain information about contemporaneous effects between the variables. Like the error variance \(\sigma^2\) in a single equation model \(\Sigma\) is essential for the calculation of test statistics and confidence intervals.

The covariance matrices of standard VAR models are symmetric , i.e. the elements to the top-right of the diagonal (the “upper triangular”) mirror the elements to the bottom-left of the diagonal (the “lower triangular”). This reflects the idea that the relations between the endogenous variables only reflect correlations and do not allow to make statements about causal relationships, since the effects are the same in each direction. This is the reason why this model is said to be not uniquely identified .

Contemporaneous causality or, more precisely, the structural relationships between the variables is analysed in the context of so-called structural VAR ( SVAR ) models, which impose special restrictions on the covariance matrix – and depending on the model on other matrices as well – so that the system is identified. This means that there is only one unique solution for the model and it is clear, how the causalities work. The drawback of this approach is that it depends on the more or less subjective assumptions made by the researcher. 1 For many researchers this is too much subjectiv information, even if sound economic theory is used to justify those assumptions.

In this article I consider a VAR(2) process of the form

\[\begin{pmatrix} y_{1,t}\\ y_{2,t} \end{pmatrix} = \begin{bmatrix} -0.3 & -0.4 \\ 0.6 & 0.5 \end{bmatrix} \begin{pmatrix} y_{1,t-1} \\ y_{2,t-1} \end{pmatrix} + \begin{bmatrix} -0.1 & 0.1 \\ -0.2 & 0.05 \end{bmatrix} \begin{pmatrix} y_{1,t-2} \\ y_{2,t-2} \end{pmatrix} + \begin{pmatrix} \epsilon_{1t} \\ \epsilon_{2t} \end{pmatrix}\]

with \(\epsilon_{1t} \sim N(0, 0.5)\) and \(\epsilon_{2t} \sim N(0, 0.5)\) . Note that for simplification the errors are not correlated. Models with correlated errors are described in a post on SVAR .

The artificial sample for this example is generated in R with

The estimation of the parameters and the covariance matrix of a simple VAR model is straightforward. For \(Y = (y_{1},..., y_{T})\) and \(Z = (z_{1},..., z_{T})\) with \(z\) as a vector of lagged valus of \(y\) and possible deterministic terms the least squares estimator of the parameters is \(\hat{A} = YZ(ZZ')^{-1}\) . The covariance matrix is then obtained from \(\frac{1}{T-Q}(Y-\hat{A}Z) (Y-\hat{A}Z)'\) , where \(Q\) is the number of estimated parameters. These formalas are usually already programmed in standard statistics packages for basic applications.

In order to estimate the VAR model I use the vars package by Pfaff (2008). The relevant function is VAR and its use is straightforward. You just have to load the package and specify the data ( y ), order ( p ) and the type of the model. The option type determines whether to include an intercept term, a trend or both in the model. Since the artificial sample does not contain any deterministic term, we neglect it in the estimation by setting type = "none" .

Model comparison

A central issue in VAR analysis is to find the number of lags, which yields the best results. Model comparison is usually based on information criteria like the AIC, BIC or HQ. Usually, the AIC is preferred over other criteria, due to its favourable small sample forecasting features. The BIC and HQ, however, work well in large samples and have the advantage of being a consistent estimator of the true order, i.e. they prefer the true order of the VAR model - in contrast to the order, which yields the best forecasts - as the sample size grows.

The VAR function of the vars package already allows to calculate standard information criteria to find the best model. In this example we use the AIC:

Note that instead of specifying the order p , we now set the maximum lag length of the model and the information criterion used to select the best model. The function then estimates all five models, compares them according to their AIC values and automatically selects the most favourable. Looking at summary(var.aic) we see that the AIC suggests to use an order of 2 which is the true order.

Looking at the results more closely we can compare the true values - cointained in object A - with the parameter estimates of the model:

All the estimates have the right sign and are relatively close to their true values. I leave it to you to look at the standard errors of summary(var.aic) to check whether the true values fall into the confidence bands of the estimates.

Impulse response

Once we have decided on a final VAR model its estimated parameter values have to be interpreted. Since all variables in a VAR model depend on each other, individual parameter values only provide limited information on the reaction of the system to a shock. In order to get a better intuition of the model’s dynamic behaviour, impulse responses (IR) are used. They give the reaction of a response variable to a one-time shock in an impulse variable. The trajectory of the response variable can be plotted, which results in those wavy curves that can be found in many macro papers.

In R the irf function of the vars package can be used to obtain an impulse response function. In the following example, we want to know how Series 2 behaves after a shock to Series 1. After specifying the model and the variables for which we want an impulse response we set the time horizon n.ahead to 20. The plot gives the response of series 2 for the periods 0 to 20 to a shock in series 1 in period 0. The function also automatically calculates so-called bootstrap confidence bands. (Bootstrapping is a common procedure in impulse response analysis. But you should keep in mind that it has its drawbacks when you work with structural VAR models though.)

Note that the ortho option is important, because it says something about the contemporaneous relationships between the variables. In our example we already know that such relationships do not exist, because the true variance-covariance matrix – or simply covariance matrix – is diagonal with zeros in the off-diagonal elements. However, since the limited time series data with 200 observations restricts the precision of the parameter estimates, the covariance matrix has positive values in its off-diagonal elements which implies non-zero contemporaneous effects of a shock. To rule this out in the IR, we set ortho = FALSE . The result of this is that the impulse response starts at zero in period 0. You could also try out the alternative and set ortho = TRUE , which results in a plot that start below zero. I do not want to go into more detail here, but suffice it so say that the issue of so-called orthogonal errors is one of the central problems in VAR analysis and you should definitely read more about it, if you plan to set up your own VAR models.

Sometimes it is interesting to see what the long-run effects of a shock are. To get an idea about that you can also calculate and plot the cumulative impulse response function to get an idea of the overall long-run effect of the shock:

We see that although the reaction of series 2 to a shock in series 1 is negative during some periods, the overall effect significantly positive.

Lütkepohl, H. (2007). New Introduction to Multiple Time Series Analyis . Berlin: Springer.

Bernhard Pfaff (2008). VAR, SVAR and SVEC Models: Implementation Within R Package vars. Journal of Statistical Software 27(4).

Sims, C. (1980). Macroeconomics and Reality. Econometrica , 48(1), 1–48.

If you want to impress your professor, you can also use the term Wold-ordering problem to refer to this issue. ↩

VAR Models in Macroeconomics – New Developments and Applications: Essays in Honor of Christopher A. Sims: Volume 32

Table of contents, var models in macroeconomics - new developments and applications: essays in honor of christopher a. sims, advances in econometrics, copyright page, list of contributors, introduction, the relationship between dsge and var models.

This article reviews the literature on the econometric relationship between DSGE and VAR models from the point of view of estimation and model validation. The mapping between DSGE and VAR models is broken down into three stages: (1) from DSGE to state-space model; (2) from state-space model to VAR( ∞ ); (3) from VAR( ∞ ) to finite-order VAR. The focus is on discussing what can go wrong at each step of this mapping and on critically highlighting the hidden assumptions. I also point out some open research questions and interesting new research directions in the literature on the econometrics of DSGE models. These include, in no particular order: understanding the effects of log-linearization on estimation and identification; dealing with multiplicity of equilibria; estimating nonlinear DSGE models; incorporating into DSGE models information from atheoretical models and from survey data; adopting flexible modeling approaches that combine the theoretical rigor of DSGE models and the econometric model’s ability to fit the data.

Do DSGE Models Forecast More Accurately Out-Of-Sample than VAR Models? ☆ The views expressed in this article are those of the authors.

Recently, it has been suggested that macroeconomic forecasts from estimated dynamic stochastic general equilibrium (DSGE) models tend to be more accurate out-of-sample than random walk forecasts or Bayesian vector autoregression (VAR) forecasts. Del Negro and Schorfheide (2013) in particular suggest that the DSGE model forecast should become the benchmark for forecasting horse-races. We compare the real-time forecasting accuracy of the Smets and Wouters (2007) DSGE model with that of several reduced-form time series models. We first demonstrate that none of the forecasting models is efficient. Our second finding is that there is no single best forecasting method. For example, typically simple AR models are most accurate at short horizons and DSGE models are most accurate at long horizons when forecasting output growth, while for inflation forecasts the results are reversed. Moreover, the relative accuracy of all models tends to evolve over time. Third, we show that there is no support to the common practice of using large-scale Bayesian VAR models as the forecast benchmark when evaluating DSGE models. Indeed, low-dimensional unrestricted AR and VAR forecasts may forecast more accurately.

Unit Roots, Cointegration, and Pretesting in Var Models ☆ The views expressed here are the authors and not necessarily those of the Federal Reserve Bank of Atlanta or the Federal Reserve System.

This article investigates the robustness of impulse response estimators to near unit roots and near cointegration in vector autoregressive (VAR) models. We compare estimators based on VAR specifications determined by pretests for unit roots and cointegration as well as unrestricted VAR specifications in levels. Our main finding is that the impulse response estimators obtained from the levels specification tend to be most robust when the magnitude of the roots is not known. The pretest specification works well only when the restrictions imposed by the model are satisfied. Its performance deteriorates even for small deviations from the exact unit root for one or more model variables. We illustrate the practical relevance of our results through simulation examples and an empirical application.

Evaluating the Accuracy of Forecasts from Vector Autoregressions ☆ The views expressed herein are solely those of the authors and do not necessarily reflect the views of the Federal Reserve Bank of Cleveland, Federal Reserve Bank of St. Louis, Federal Reserve System, or any of its staff.

This article surveys recent developments in the evaluation of point and density forecasts in the context of forecasts made by vector autoregressions. Specific emphasis is placed on highlighting those parts of the existing literature that are applicable to direct multistep forecasts and those parts that are applicable to iterated multistep forecasts. This literature includes advancements in the evaluation of forecasts in population (based on true, unknown model coefficients) and the evaluation of forecasts in the finite sample (based on estimated model coefficients). The article then examines in Monte Carlo experiments the finite-sample properties of some tests of equal forecast accuracy, focusing on the comparison of VAR forecasts to AR forecasts. These experiments show the tests to behave as should be expected given the theory. For example, using critical values obtained by bootstrap methods, tests of equal accuracy in population have empirical size about equal to nominal size.

Identifying Structural Vector Autoregressions Via Changes in Volatility ☆ This article was written while the author was a Bundesbank Professor at the Freie Universität Berlin. An earlier version of the paper was published as DIW Discussion Paper 1259 – http://www.diw.de/sixcms/detail.php?id=diw_0.1.c.412678.de

Identification of shocks of interest is a central problem in structural vector autoregressive (SVAR) modeling. Identification is often achieved by imposing restrictions on the impact or long-run effects of shocks or by considering sign restrictions for the impulse responses. In a number of articles changes in the volatility of the shocks have also been used for identification. The present study focuses on the latter device. Some possible setups for identification via heteroskedasticity are reviewed and their potential and limitations are discussed. Two detailed examples are considered to illustrate the approach.

Panel Vector Autoregressive Models: A Survey ☆ The views expressed in this article are those of the authors and do not necessarily reflect those of the ECB or the Eurosystem.

This article provides an overview of the panel vector autoregressive models (VAR) used in macroeconomics and finance to study the dynamic relationships between heterogeneous assets, households, firms, sectors, and countries. We discuss what their distinctive features are, what they are used for, and how they can be derived from economic theory. We also describe how they are estimated and how shock identification is performed. We compare panel VAR models to other approaches used in the literature to estimate dynamic models involving heterogeneous units. Finally, we show how structural time variation can be dealt with.

Mixed-Frequency Vector Autoregressive Models ☆ This views expressed herein are solely those of the authors and do not necessarily reflect the views of the Norges Bank. The usual disclaimers apply.

The development of models for variables sampled at different frequencies has attracted substantial interest in the recent literature. In this article, we discuss classical and Bayesian methods of estimating mixed-frequency VARs, and use them for forecasting and structural analysis. We also compare mixed-frequency VARs with other approaches to handling mixed-frequency data.

Thresholds and Smooth Transitions in Vector Autoregressive Models ☆ The views expressed in this article are those of the authors and should not be interpreted as reflecting the views of the European Central Bank.

This survey focuses on two families of nonlinear vector time series models, the family of vector threshold regression (VTR) models and that of vector smooth transition regression (VSTR) models. These two model classes contain incomplete models in the sense that strongly exogeneous variables are allowed in the equations. The emphasis is on stationary models, but the considerations also include nonstationary VTR and VSTR models with cointegrated variables. Model specification, estimation and evaluation is considered, and the use of the models illustrated by macroeconomic examples from the literature.

Nonparametric Vector Autoregressions: Specification, Estimation, and Inference

For over three decades, vector autoregressions have played a central role in empirical macroeconomics. These models are general, can capture sophisticated dynamic behavior, and can be extended to include features such as structural instability, time-varying parameters, dynamic factors, threshold-crossing behavior, and discrete outcomes. Building upon growing evidence that the assumption of linearity may be undesirable in modeling certain macroeconomic relationships, this article seeks to add to recent advances in VAR modeling by proposing a nonparametric dynamic model for multivariate time series. In this model, the problems of modeling and estimation are approached from a hierarchical Bayesian perspective. The article considers the issues of identification, estimation, and model comparison, enabling nonparametric VAR (or NPVAR) models to be fit efficiently by Markov chain Monte Carlo (MCMC) algorithms and compared to parametric and semiparametric alternatives by marginal likelihoods and Bayes factors. Among other benefits, the methodology allows for a more careful study of structural instability while guarding against the possibility of unaccounted nonlinearity in otherwise stable economic relationships. Extensions of the proposed nonparametric model to settings with heteroskedasticity and other important modeling features are also considered. The techniques are employed to study the postwar U.S. economy, confirming the presence of distinct volatility regimes and supporting the contention that certain nonlinear relationships in the data can remain undetected by standard models.

Testing for Common Cycles in Non-Stationary VARs with Varied Frequency Data

This article proposes a new approach to detecting the presence of common cyclical features when several time series are sampled at different frequencies. We generalize the common-frequency approach introduced by Engle and Kozicki (1993) and Vahid and Engle (1993) . We start with the mixed-frequency VAR representation investigated in Ghysels (2012) for stationary time series. For non-stationary time series in levels, we show that one has to account for the presence of two sets of long-run relationships. The first set is implied by identities stemming from the fact that the differences of the high-frequency I (1) regressors are stationary. The second set comes from possible additional long-run relationships between one of the high-frequency series and the low-frequency variables. Our transformed vector error-correction model (VECM) representations extend the results of Ghysels (2012) and are important for determining the correct set of variables to be used in a subsequent common cycle investigation. This fact has implications for the distribution of test statistics and for forecasting. Empirical analyses with quarterly real gross national product (GNP) and monthly industrial production indices for, respectively, the United States and Germany illustrate our new approach. We also conduct a Monte Carlo study which compares our proposed mixed-frequency models with models stemming from classical temporal aggregation methods.

Multivariate Dynamic Probit Models: An Application to Financial Crises Mutation

In this article we propose a multivariate dynamic probit model. Our model can be viewed as a nonlinear VAR model for the latent variables associated with correlated binary time-series data. To estimate it, we implement an exact maximum likelihood approach, hence providing a solution to the problem generally encountered in the formulation of multivariate probit models. Our framework allows us to study the predictive relationships among the binary processes under analysis. Finally, an empirical study of three financial crises is conducted.

All feedback is valuable

Please share your general feedback

Report an issue or find answers to frequently asked questions

Contact Customer Support

Implications of Dynamic Factor Models for VAR Analysis

This paper considers VAR models incorporating many time series that interact through a few dynamic factors. Several econometric issues are addressed including estimation of the number of dynamic factors and tests for the factor restrictions imposed on the VAR. Structural VAR identification based on timing restrictions, long run restrictions, and restrictions on factor loadings are discussed and practical computational methods suggested. Empirical analysis using U.S. data suggest several (7) dynamic factors, rejection of the exact dynamic factor model but support for an approximate factor model, and sensible results for a SVAR that identifies money policy shocks using timing restrictions.

• Acknowledgements and Disclosures

MARC RIS BibTeΧ

Download Citation Data

More from NBER

In addition to working papers , the NBER disseminates affiliates’ latest findings through a range of free periodicals — the NBER Reporter , the NBER Digest , the Bulletin on Retirement and Disability , the Bulletin on Health , and the Bulletin on Entrepreneurship  — as well as online conference reports , video lectures , and interviews .

An official website of the United States government

The .gov means it’s official. Federal government websites often end in .gov or .mil. Before sharing sensitive information, make sure you’re on a federal government site.

The site is secure. The https:// ensures that you are connecting to the official website and that any information you provide is encrypted and transmitted securely.

• Publications
• Account settings

Preview improvements coming to the PMC website in October 2024. Learn More or Try it out now .

• Advanced Search
• Journal List
• HHS Author Manuscripts

Vector Autoregressive (VAR) Models and Granger Causality in Time Series Analysis in Nursing Research: Dynamic Changes Among Vital Signs Prior to Cardiorespiratory Instability Events as an Example

Eliezer bose.

Assistant Professor, School of Nursing, University of Pittsburgh, Pittsburgh, PA

Marilyn Hravnak

Professor, School of Nursing, University of Pittsburgh, Pittsburgh, PA

Susan M. Sereika

Professor School of Nursing, University of Pittsburgh Pittsburgh, PA

Associated Data

Patients undergoing continuous vital sign monitoring (heart rate [HR], respiratory rate [RR], pulse oximetry [SpO 2 ]) in real time display inter-related vital sign changes during situations of physiologic stress. Patterns in this physiological cross-talk could portend impending cardiorespiratory instability (CRI). Vector autoregressive (VAR) modeling with Granger causality tests is one of the most flexible ways to elucidate underlying causal mechanisms in time series data.

The purpose of this article is to illustrate development of patient-specific VAR models using vital sign time series (VSTS) data in a sample of acutely ill, monitored, step-down unit (SDU) patients, and determine their Granger causal dynamics prior to onset of an incident CRI.

CRI was defined as vital signs beyond stipulated normality thresholds (HR = 40–140/minute, RR = 8–36/minute, SpO 2 < 85%) and persisting for 3 minutes within a 5-minute moving window (60% of the duration of the window). A 6-hour time segment prior to onset of first CRI was chosen for time series modeling in 20 patients using a six-step procedure: (a) the uniform time series for each vital sign was assessed for stationarity; (b) appropriate lag was determined using a lag-length selection criteria; (c) the VAR model was constructed; (d) residual autocorrelation was assessed with the Lagrange Multiplier test; (e) stability of the VAR system was checked; and (f) Granger causality was evaluated in the final stable model.

The primary cause of incident CRI was low SpO 2 (60% of cases), followed by out-of-range RR (30%) and HR (10%). Granger causality testing revealed that change in RR caused change in HR (21%) (i.e., RR changed before HR changed) more often than change in HR causing change in RR (15%). Similarly, changes in RR caused changes in SpO 2 (15%) more often than changes in SpO 2 caused changes in RR (9%). For HR and SpO 2 , changes in HR causing changes in SpO 2 and changes in SpO 2 causing changes in HR occurred with equal frequency (18%).

Within this sample of acutely ill patients who experienced a CRI event, VAR modeling indicated that RR changes tend to occur before changes in HR and SpO 2. These findings suggest that contextual assessment of RR changes as the earliest sign of CRI is warranted. Use of VAR modeling may be helpful in other nursing research applications based on time series data.

Many nursing phenomena are time-dependent. Beginning early in their training years and beyond, nurses are taught to measure changes in behavior, patient outcomes, and responses occurring along a delineated time interval. Ever since Nurse Linda Richard’s invention of bedside records ( Holder, 2004 ), flow sheet usage and time series data have been ubiquitous and important in nursing. More recently, real-time data capture of vital signs using telemetry and portable monitors has created the opportunity to time stamp longitudinal data. Such time series data have potential to reveal patterns in manifestations of the variables at specific time intervals ( Chatterjee & Price, 2009 ). However, accurate and timely human-only interpretation of temporal patterns or trends in voluminous time series data is impossible. Finding patterns in time series data for nursing research purposes, particularly in identifying temporal pattern emergence prior to critical events, involves proper understanding and use of appropriate mathematical models to study changes across time at the level of the individual.

Models for time series analysis are idiographic in nature. This means that they are tools for analyzing unique and patient-specific fluctuations within a time series. As such, time series approaches provide a framework for analyzing future changes for an individual based on past trends and patterns. Predicting how a particular patient’s time series data could behave based on past trends has numerous implications for bed-side nursing practice. It provides a framework allowing the uniqueness of each patient to exist as a basis for assessing change over time—not as a deviation from a predetermined pattern generalizable to all patients, but as an alteration in a personal pattern.

Time series data occur when sets of observations from an individual case are arranged in temporal order. The primary interest is the relationship of the values from one point in time to the next for individual cases. In order to develop reasonable models for illuminating interrelationships among the points in a time series, large datasets with observations distributed densely over time are needed. With real-time data acquisition systems now in place and “big data” ubiquitous, researchers can now define “sampling windows.” Sampling windows are shorter researcher-defined epochs within an overall time period—defined by a duration of time (in minutes, hours, or days, as appropriate for the problem under study) and number of observations within the time period.

Time series models account for the fact that successive observations are correlated. The correlation among elements in a time series is called autocorrelation. Higher values signify greater association between an observation on a variable X (like heart rate) at some time t and at later times t + k , where k indicates the next 1, 2, up to N observations within the window, where N ∈ Z (the set of integers) is the maximum number of observations ( Chatfield, 2003 ).

Multivariate time series analysis involves more than one variable (like heart rate, respiratory rate, blood pressure, and SpO2). Strength of associations among different variables across time lags (shifts in time) are indexed with cross-correlations ( Box, Jenkins, & Reinsel, 2008 ). Several possible lags can be used to examine cross-correlations. For instance, if X t is the value of the variable in period t , lag( X t , 1) represents the value of the variable in period t-1 (lagged one period back in time).

The vector autoregression (VAR) model is one of the most commonly employed multivariate regression time series analytic techniques. The VAR model is advantageous because it can explain past and causal relationships among multiple variables over time, as well as predict future observations. Explanation and prediction of future observations in a time series is dependent upon correctly postulating a VAR model and estimating its parameters ( Lütkepohl, 2005 ).

Physiologic nursing research often depends on the collection of multivariate time series data to be used in various statistical models. For instance, a multivariate time series model was used to forecast daily attendances at the hospital emergency department of an acute care hospital; autoregressive integrated moving average models were applied to three different acuity categories (most acute, acute and least acute) of patients seen in the emergency department, and used as tools for predicting emergency department workload both for staff roster as well as for resource planning ( Sun, Heng, Seow, & Seow, 2009 ). In another study, a vector autoregressive time series model was used in a set of 450 intensive care unit patients to systematically learn and identify a collection of time series dynamics that were recurrent within each patient and that were shared across the entire cohort ( Lehman et al., 2015 ). The time series dynamics in combination with baseline acuity measures were used to better characterize the physiological state of a patient and predict in-hospital mortality ( Lehman et al., 2015 ).

Vital sign time series (VSTS) involves monitored physiologic vital signs (VS), such as heart rate (HR), respiratory rate (RR) or pulse oximetry (SpO 2 ). (See Figure 1 .) Physiologically, dynamic interactions among the VS are indicative of a patient’s efforts to compensate and attain homeostasis in response to values that exceed thresholds marking normal limits ( Scharf, Pinsky, & Magder, 2001 ), where failure to compensate leads to cardiorespiratory instability (CRI). Patient-specific multivariate time series modeling approaches ( Schulz et al., 2013 ) can be used to study the causal characteristics of the changes in physiological variables as they progress in time towards the development of CRI, which is a well-defined event that occurs or does not occur within some period of time. The vector autoregressive (VAR) model is one of the most flexible models for the analysis of causality in multivariate time series.

Time series plots for about 1 hour showing heart rate (HR), respiratory rate (RR) and pulse oximetry (SpO 2 ). The shaded portion represents incident cardiorespiratory instability (CRI) which this patient developed due to crossing of HR threshold limits, with changes in RR and SpO 2 occurring prior to CRI. Vector autoregressive (VAR) models are able to look at these evolving patterns over time for this specific patient to assess VS changes prior to incident CRI.

Granger causality is a concept of causality derived from the notion that causes may not occur after effects and that if one variable is the cause of another, knowing the status on the cause at an earlier point in time can enhance prediction of the effect at a later point in time ( Granger, 1969 ; Lütkepohl, 2005 , p. 41). The VAR model has been widely employed in econometric analyses ( Granger & Newbold, 1986 ) and in neurobiology ( Tang, Bressler, Sylvester, Shulman, & Corbetta, 2012 ) to elucidate underlying mechanisms using Granger causality. To the best of our knowledge, VAR has never been implemented in studying causality in physiologic time series variables of HR, RR and SpO 2 as multivariate inputs to the development of CRI.

We therefore provide an example of using the VAR approach. The purpose of this study was to develop patient-specific multivariate time series VAR models using HR, RR and SpO 2 in a sample of acutely ill monitored step-down unit (SDU) patients. VAR model provided the framework in order to study the Granger casual dynamics among the VS leading up to CRI. Since CRI can occur multiple times during the hospitalization of unstable patients, we decided to consider the sampling window as only their first instance of CRI, which we called CRI 1 .

The examples involve prospectively-collected physiologic data in patients who underwent continuous noninvasive VS monitoring. Heart rate (HR), respiratory rate (RR; bioimpedance), peripheral oximetry (SpO 2 ; plethysmography) were recorded once every 20 seconds (Δ t = 20 seconds) for the entire monitoring period of their hospitalization on the SDU. Instability threshold limits for VS deviations outside of normal were set (HR: ≤ 40 or ≥140 beats per minute; RR: ≤ 8 or ≥36 breaths per minute; and SpO 2 : ≤ 85%). Instability epochs were denoted by persistence of any deviation in VS for 3 minutes (60% duty cycle, which is the percentage of period the signal is active within a sampling window) of a 5 minute moving window. The instability epochs were then visually annotated by two expert critical care clinician reviewers. The instability epochs were coded as actual CRI or monitoring artifacts. The incident (first) CRI was denoted CRI 1 . After rejecting artifact, we created a subset of data for CRI 1 instances and VSTS up to a maximum of 6 hours prior to the incident CRI and saved it in a flat file using comma-separated value (.csv) format. (Flat files using .csv formats enable compact storage of large amounts of data making it easy to store, retrieve and import across applications.)

Sample and Setting

Patients were recruited from a 24-bed step-down unit (SDU) located in a tertiary academic center. Those who had at least two hours of continuous noninvasive monitoring of HR, RR, and SpO 2 were eligible. Data were extracted from input obtained from noninvasive beside monitors (model M1204, Philips Medical, Bothell, WA) that included continuous HR (3-lead electrocardiogram), RR (bioimpedance signal), and plethysmographic peripheral blood oxygen saturation (SpO 2 ) signal (model M1191B, Philips Medical Systems). Additional inclusion criteria were admission to a monitored bed on the SDU and age ≥21 years. All patients were admitted to the study unit following the usual standard of care for monitored bed admission and utilization. A convenience sample of 20 patients took part; ages ranged from 35 to 92 years; men and women were equally represented; and most were White (80%). The study protocol was approved by the University Institutional Review Board (IRB).

Data Set and Data Preprocessing

The start of the CRI 1 epoch was indicated by the first point in the patient’s streaming data at which any of a patient’s VS first crossed an instability threshold. The end of the epoch was indicated by the point at which all VS returned to normal values. Using the beginning of the CRI 1 epoch, a 6-hour window period prior to that time was chosen for VAR analysis. This was because Hravnak et al. (2008) , in the same patient population, determined that VS changes occur at least 6 hours prior to the onset of CRI.

In the current study, the first step was to preprocess the data to verify that a data value existed at every 20 seconds. Gaps were noticed within certain portions of the data stream (about 5–8% missingness). Missingness was likely due to poor signal capture. Gaps in the data stream were populated using linear interpolation between the points ( Lehman, Nemati, Adams, & Mark, 2012 ) to make the data stream continuous for each of the three vital sign variables. We used this strategy since linear interpolation assumes that the unknown value of any given physiologic variable lies on the line between the two known values, similar to that of normal physiological mechanisms. Once missing values had been replaced with interpolated values, the data set was ready for VAR modeling.

Vector Autoregressive Modeling

A VAR( p ) model for a multivariate time series is a regression model for outcomes at a specified time t and time lagged predictors, with p indicating the lag (e.g., p = 1 refers to the observation previous to t; p = 2 refers to two observations prior to t , and so on ). Key terms used in VAR modeling are defined in Table 1 . The VAR model for vital sign time series data was developed using the following steps:

Glossary of Terms

• Stationarity of the individual VSTS (HR, RR, SpO 2 ) is tested.
• Lag is determined using lag-length selection criteria.
• A VAR model with appropriate lags is built.
• Residual autocorrelation is assessed with the Lagrange Multiplier (LM) test.
• Stability of the VAR system is assessed with the autoregressive (AR) roots graph.
• The Granger causality test is performed.

Testing Stationarity

Stationary time series are detrended series, without periodic fluctuations. Stationarity is critical to development of a VAR model because in its absence, a model’s statistics such as means and correlations will not accurately describe the time series signal. The augmented Dickey-Fuller test is used to test stationarity. The null hypothesis is that the time series is nonstationary and the alternative is that the series is stationary. Rejection of the null hypothesis indicates that the series does not need transformation to achieve stationarity and modeling can proceed to the lag length selection step.

Selecting Lag Length

Lag length selection refers to the number of previous observations in a time series that will be used as predictors in the VAR model. Typically, a large number of lags will be used to generate a model and then a restriction applied to select a more parsimonious model. Lütkepohl (2005) indicated that using too few lags can result in autocorrelated errors whereas using too many lags results in over-fitting, causing an increase in mean-square-forecast errors of the VAR model. Selection of an appropriate lag is critical to inference in VARs. The lag length for the VAR ( p ) model can be determined using model selection criteria. The most common approach is to fit VAR ( p ) models with orders p = 0, 1, … p max and choose the value of p which minimizes some model selection criteria. The three most commonly used criteria used are the Akaike Information Criterion (AIC), Schwarz Bayesian Information Criterion (BIC), and the Hannan-Quinn criterion (HQ) (Lütkepohl, Section 8.1; Vrieze, 2012 ). The lag associated with the minimum value of a criterion is selected ( Lütkepohl, 2005 ).

Building the VAR Model

The general structure of the VAR model used in multivariate time series is that each variable, at a point in time, is a linear function of the most recent lag of both itself and the other variables. The general model is described in Supplemental Digital Content 1 . To illustrate, a trivariate VAR(1) model for HR, RR and SpO 2 for one case has the usual matrix form for a regression equation with multivariate outcomes, such that

or, in scalar notation,

As shown in Equations 1 – 4 , the regressors (predictors) of each outcome are the same, which are the lagged values of HR, RR and SpO 2 . Intercept terms are indicated by the c terms, regression coefficients are indicated by the subscripted A values, and error in prediction of each outcome at time t is indicated by the e terms. The equations can be solved using ordinary least squares (OLS) estimation. Since within the VAR ( p ), each equation has the same explanatory variables, each equation may be estimated separately. Recommended sources for information about estimation in VAR modeling are Lütkepohl (2005) and Hamilton (1994) .

Testing for Residual Autocorrelation

Once a VAR model has been developed, the next step is to determine if the selected model provides an adequate description of the data. In familiar regression models, this is performed by examining the residuals, which are differences between the actual observations and model-fitted values. In time series models, autocorrelation of the residual values is used to determine the goodness of fit of the model. Autocorrelation of the residuals indicates that there is information that has not been accounted for in the model. The Lagrange Multiplier (LM) test is a standard tool for checking residual autocorrelation in VAR models. The null hypothesis is that there is no residual autocorrelation; the alternative is that residual autocorrelation exists ( Lütkepohl, 2005 , Section 4.4.4). Technical detail is provided in Supplemental Digital Content 1 .

Assessing Stability of the VAR Model

Stability refers to checking whether the model is a good representation of how the time series evolved over the sampling window period. Technically, stability of a VAR system is evaluated using the roots of the characteristic polynomial of the coefficient matrix A , as in the example in Equation 1 . Stability in a VAR model is indicated by roots that are all less than 0, and are typically shown in a graph. See Supplemental Digital Content 1 for more detail.

Evaluating Granger Causality

VAR models describe the joint generation process of a number of variables over time, so they can be used for investigating relationships between the variables. Granger causality is one type of relationship between time series ( Granger, 1969 ). The basic idea of Granger causality (GC) can be stated as if the prediction of one time series is improved by incorporating the knowledge of a second time series, then the latter is said to have a causal influence on the first . Specifically, two autoregressive models are fitted to the first time series—with and without including the second time series—and the improvement of the prediction is measured by the ratio of the variance of the error terms. The null hypothesis for GC is that no explanatory power is added by jointly considering the lagged values of y and x as predictors. The null hypothesis that x does not cause y is rejected if coefficients for the lagged values of x are significant; i.e., Granger called a variable x causal for a variable y if the lagged values of x are helpful for improving forecasts of y ( y at future times). The VAR framework is flexible and provides an environment for implementing this type of analysis.

In the current study, all data analysis was performed using EViews 8 Student version, a Windows-based econometric program commonly employed in the financial industry for econometric analysis, forecasting and simulation. Time series data for each of the individual variables of HR, RR and SpO 2 were provided as input after importing .csv files of time series data.

Application in Example Data

The time series data for each of the individual VS variables were first visualized using a graphical line plot of each of the variables over time ( Figure 1 ). Graphs were assessed for trends and stationarity. From a visual perspective, the graphs verified that the series looked flat, without any upward or downward trend and with no periodic fluctuations over time.

Stationarity

The test for stationarity was performed individually for HR, RR, and SpO 2 for each patient. EViews 8 outputs an augmented Dickey Fuller (ADF) t -statistic value and its corresponding p -value for each variable. A p -value <.05 indicated that the time series was stationary for that particular variable. If the VS time series was not stationary, the series was differenced and the ADF test applied again on the differenced time series to check for stationarity. Among the 60 series (20 patients x 3 vital signs), 12 series were stationary when first tested and the rest were differenced once to achieve stationarity.

Once stationarity of a series was determined, we proceeded with lag-length selection using AIC and HQ. Table 2 shows an example case, were the 12 th lag was associated with the minimum AIC value. In selecting the 12 th lag, all the lags from 1 to 12 are subsequently included in the models. Lags ranged from 8 to 20 among the 20 cases for the HR, RR, and SpO 2 time series.

Lag Order Selection, Assessment of Residual Correlation, and Granger Causality for Heart Rate for One Patient

Note. Lag order begins with an initial lag order set high (in this case, lags). The “best” lag is the lag with smallest AIC value, associated with a lag of 12. The LGM evaluates serial correlation (autocorrelation) for a single variables; in this case, none of the LGM tests were significant, indicating that past values of HR did not affect future values for this patient, up to the specified lag order. Granger causality treats each variable as dependent and other variables with all of their lags as independent; p -value <.05 is suggestive of causality between the dependent variable and independent variable. Using this criterion in this time series, only the effect of RR on HR was significant (p <.04). AIC = Akaike Information Criterion; df = degrees of freedom; HQ = Hannan-Quinn Information Criterion; HR = heart rate; LGM = Lagrange Multiplier; RR = respiratory rate; SpO 2 = peripheral capillary oxygen saturation.

VAR model with appropriate lags

VAR models were estimated to include the number of lags from 1 until the specified lag as described above. Table 2 shows the example of a VAR model with lags from 1 to 12 for HR for one patient. Models were also estimated for RR and SpO 2. Significance of the coefficients was assessed to decide whether another VAR model should be constructed so as to minimize the nonsignificant coefficients. In all, only 5 models were re-estimated because coefficients were nonsignificant.

Residual autocorrelation

We next applied the LM test to check for residual serial correlation for each patient. In general, after fitting a VAR model, the residuals should be white noise and should have no autocorrelation. If autocorrelation among the residuals is noticed, then it implies that there was some information which was not accounted for by the model, such as insufficient lags. The null hypothesis of the LM test is that there is no serial correlation up to the specified lag order (past values of a variable do not affect the future value of the variable until the specified lag value). If a proper VAR system is obtained from the previous step, then the coefficients of the LM statistic are not significant, thus accepting the null hypothesis up to the specified lag. This is also illustrated in Table 2 for one patient.

Stability of the VAR system was evaluated using the roots of the characteristic polynomial with the variables of HR, RR and SpO 2 and the lag specification as per above criteria. Figure 2 shows stability evaluation using a graph of the characteristic roots for one patient. This evaluation was performed on a patient-by-patient basis. Of the 20 cases, models for 18 cases were stable ( Table 3 ), based on inspection of the graphs.

The stability of the VAR system is visually inspected by evaluating roots of the above VAR system. The dots represent the roots of characteristic polynomial. If no roots lie outside the unit circle, then the VAR system satisfies stability condition, indicating model stability.

Granger Causality for Each Patient

Note. Unidirectional relationships are indicated with (→) and bidirectional relationships are shown with (←→). CRI 1 = incident cardiorespiratory instability; HR = heart rate; ID = patient identification code; RR = respiratory rate; SpO2 = peripheral capillary oxygen saturation; VS = vital sign.

Granger causality

Granger causality assessment for one patient is also shown in Table 2 . In this case, all of the RR lags caused HR. This demonstrates that there was unidirectional causality running in one direction, from RR to HR, such that changes in RR causing changes in HR in the patient’s time series prior to CRI 1 .

Table 3 summarizes Granger causality assessment for all 20 cases. No causal structure was identified for three cases, and VAR models were unstable for two cases. The first VS over threshold to begin the CRI epoch was SpO 2 in 60% cases, followed by RR (30%) and HR (10%). Out of the 20 cases, assessment of Granger causality revealed that RR caused HR (i.e., RR changed before HR changed) more often than HR caused RR (21% vs 15%, respectively). Similarly, changes in RR causing changes in SpO 2 was more common than changes in SpO 2 causing RR (15% vs 9%, respectively). For HR and SpO 2 , changes in HR causing changes in SpO 2 and changes in SpO 2 causing changes in HR occurred with equal frequency (18%). There were three cases with bidirectional causality (#8, #14 and #18).

To demonstrate the utility of VAR modeling in nursing research, we carried out an analysis with the purpose of developing a stable patient-specific multivariate time series VAR model— using HR, RR and SpO 2 in a sample of SDU patients—in order to study the Granger casual dynamics among the monitored vital signs leading up to a first CRI. Our results suggest that using this information may be helpful to determine causality in VS threshold deviations, and determine a physiologic cause. For example in our 20 cases, SpO 2 falling below threshold was the first indicator of CRI. SpO 2 itself was generally not the cause for changed RR or HR, but rather RR and HR tended to be the cause for the fall in SpO 2 . This suggests that SpO 2 is a later VS change, and is preceded by more subtle changes in other VS which in turn caused changes in SpO 2 .

From an analysis standpoint, two cases had data only for 2 and 4 hours (cases #6 and #11, respectively) in the 6-hour window preceding CRI 1 . Also, there were two cases with VSTS from only two variables prior to CRI 1 —HR and RR—(# 3) and HR and SpO 2 (# 8). We also encountered two cases (#9 and #16) with intermittent monitoring (the reasons for the patients being off the monitor were not known) for an overall missing window period of at least 3 hours, within the 6-hour period prior to CRI 1 . In these two cases, we used the Markov Chain Monte Carlo (MCMC) approach to impute missing values to the dataset, since it was not possible to fill the data stream with linear interpolation. Briefly, for a continuous variable with missing values, MCMC uses the non-missing values to find the sample mean and standard deviation for the variable and then fills in the missing values with random draws from a normal distribution with mean and standard deviation equal to the sample values, limited within the range of the observed minimum and maximum values ( Geyer, 2011 ; Young, Weckman, & Holland, 2011 ).

No smoothing or other filtering approaches were used to remove artifacts in our dataset—this was done by prior expert review before data preprocessing. Filtering approaches disturb the information content, leading to spurious and missed causalities ( Florin, Gross, Pfeifer, Fink, & Timmermann, 2010 ). Hence, no filtering approaches were used, but the series was differenced to maintain stationarity, which is a prerequisite for inferring Granger causality. Finally, there were two cases (#12 and #20) where the VAR model could not be stabilized for reasons that were unclear. Interestingly, HR crossed the normality threshold first in both cases. It is possible that HR is less likely to induce compensatory changes than other VS, but this would require further exploration.

Limitations

A limitation to our study was the absence of continuous streaming data from all SDU patients for the entire 6-hour period. This is common in clinical care for a variety of reasons such as hallway ambulation, bathroom privileges, off-unit for testing or therapy, however. Linear interpolation or MCMC could have been avoided with totally uninterrupted data streams for all cases. We implemented these imputation procedures in order to maintain completeness of the VSTS for VAR modeling. Further, even though the VAR model is quite rich in capturing the temporal dynamics of highly complex hemodynamic responses, Granger causal relationships should not be over-interpreted as truly causal. An advanced discussion about Granger causality and causal inference with multiple time series is available in Eichler (2013) .

Future Directions

Our study provides an example demonstrating that the VAR modeling approach is able to expose Granger causal dynamics involving vital signs in evolving CRI. Future studies with larger samples and longer periods of VSTS monitoring to study casual dynamics are warranted. Such information could enhance target surveillance monitoring, enabling nurses to better recognize impending CRI.

Supplementary Material

Supplemental data file _doc_ pdf_ etc._.

Supplemental Digital Content 1. Technical details are provided in this document. .pdf

Acknowledgments

The research was supported by NIH Grant 1R01NR013912, National Institute of Nursing Research.

The authors have no conflicts of interest to declare.

Contributor Information

Eliezer Bose, Assistant Professor, School of Nursing, University of Pittsburgh, Pittsburgh, PA.

Marilyn Hravnak, Professor, School of Nursing, University of Pittsburgh, Pittsburgh, PA.

Susan M. Sereika, Professor School of Nursing, University of Pittsburgh Pittsburgh, PA.

• Box GEP, Jenkins GM, Reinsel GC. Time series analysis: Forecasting and control. 4. Hoboken, NJ: Wiley; 2008. [ Google Scholar ]
• Chatfield C. The analysis of time series: An introduction. 6. Boca Raton, FL: Chapman & Hall/CRC; 2003. [ Google Scholar ]
• Chatterjee S, Price A. Healthy living with persuasive technologies: Framework, issues, and challenges. Journal of the American Medical Informatics Association. 2009; 16 :171–178. doi: 10.1197/jamia.M2859. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Eichler M. Causal inference with multiple time series: Principles and problems. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences. 2013; 371 (1997) doi: 10.1098/rsta.2011.0613. [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Florin E, Gross J, Pfeifer J, Fink GR, Timmermann L. The effect of filtering on Granger causality based multivariate causality measures. NeuroImage. 2010; 50 :577–588. doi: 10.1016/j.neuroimage.2009.12.050. [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Geyer CJ. Introduction to Markov Chain Monte Carlo. In: Brooks S, Gelman A, Jones GL, Meng X-L, editors. Handbook of Markov chain Monte Carlo. Boca Raton, FL: Chapman & Hall/CRC; 2011. pp. 3–48. [ Google Scholar ]
• Granger CWJ. Investigating causal relations by econometric models and cross-spectral methods. Econometrica. 1969; 37 :424–438. doi: 10.2307/1912791. [ CrossRef ] [ Google Scholar ]
• Granger CWJ, Newbold P. Forecasting economic time series. 2. New York, NY: Academic Press; 1986. [ Google Scholar ]
• Hamilton JD. Time series analysis. Princeton, NJ: Princeton University Press; 1994. [ Google Scholar ]
• Holder VL. From handmaiden to right hand—The infancy of nursing. AORN Journal. 2004; 79 :374, 376–380, 382, 385–388, 390. doi: 10.1016/S0001-2092(06)60614-5. [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Hravnak M, Edwards L, Clontz A, Valenta C, DeVita MA, Pinsky MR. Defining the incidence of cardiorespiratory instability in patients in step-down units using an integrated electronic monitoring system. Archives of Internal Medicine. 2008; 168 :1300–1308. doi: 10.1001/archinte.168.12.1300. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Lehman LW, Adams RP, Mayaud L, Moody GB, Malhotra A, Mark RG, Nemati S. A physiological time series dynamics-based approach to patient monitoring and outcome prediction. IEEE Journal of Biomedical and Health Informatics. 2015; 19 :1068–1076. doi: 10.1109/JBHI.2014.2330827. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Lehman LW, Nemati S, Adams RP, Mark RG. Discovering shared dynamics in physiological signals: Application to patient monitoring in ICU. Engineering in Medicine and Biology Society, 2012 Annual International Conference of the IEEE; 2012. pp. 5939–5942. [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Lütkepohl H. New introduction to multiple time series analysis. Berlin, Germany: Springer Berlin Heidelberg; 2005. [ Google Scholar ]
• Scharf SM, Pinsky MR, Magder S, editors. Respiratory-circulatory interactions in health and disease. New York, NY: Marcel Dekker; 2001. [ Google Scholar ]
• Schulz S, Adochiei FC, Edu IR, Schroeder R, Costin H, Bär KJ, Voss A. Cardiovascular and cardiorespiratory coupling analyses: A review. Philosophical Transactions of the Royal Society A. 2013; 371 (1997):1–25. doi: 10.1098/rsta.2012.0191. [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Sun Y, Heng B, Seow YT, Seow E. Forecasting daily attendances at an emergency department to aid resource planning. BMC Emergency Medicine. 2009; 9 (1):1–9. doi: 10.1186/1471-227X-9-1. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Tang W, Bressler SL, Sylvester CM, Shulman GL, Corbetta M. Measuring Granger causality between cortical regions from voxelwise fMRI BOLD signals with LASSO. PLOS Computational Biology. 2012; 8 (5):1–14. doi: 10.1371/journal.pcbi.1002513. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Vrieze SI. Model selection and psychological theory: A discussion of the differences between the Akaike information criterion (AIC) and the Bayesian information criterion (BIC) Psychological Methods. 2012; 17 :228–243. doi: 10.1037/a0027127. [ PMC free article ] [ PubMed ] [ CrossRef ] [ Google Scholar ]
• Young W, Weckman G, Holland W. A survey of methodologies for the treatment of missing values within datasets: Limitations and benefits. Theoretical Issues in Ergonomic Science. 2011; 12 :15–43. doi: 10.1080/14639220903470205. [ CrossRef ] [ Google Scholar ]

•   Cadmus Home
• EUI Working Papers
• ECO Working Papers

Vector autoregressive models

Retrieved from Cadmus, EUI Research Repository

Show full item record

Files associated with this item

Collections

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here .

Loading metrics

Open Access

Peer-reviewed

Research Article

Choosing between AR(1) and VAR(1) models in typical psychological applications

Contributed equally to this work with: Fabian Dablander, Oisín Ryan, Jonas M. B. Haslbeck

Roles Conceptualization, Formal analysis, Writing – original draft, Writing – review & editing

* E-mail: [email protected]

Affiliation Department of Psychological Methods, University of Amsterdam, Amsterdam, Netherlands

Affiliation Department of Methodology and Statistics, Utrecht University, Utrecht, Netherlands

• Fabian Dablander, 
• Oisín Ryan, 
• Jonas M. B. Haslbeck

• Published: October 29, 2020
• https://doi.org/10.1371/journal.pone.0240730
• Peer Review
• Reader Comments

Time series of individual subjects have become a common data type in psychological research. The Vector Autoregressive (VAR) model, which predicts each variable by all variables including itself at previous time points, has become a popular modeling choice for these data. However, the number of observations in typical psychological applications is often small, which puts the reliability of VAR coefficients into question. In such situations it is possible that the simpler AR model, which only predicts each variable by itself at previous time points, is more appropriate. Bulteel et al. (2018) used empirical data to investigate in which situations the AR or VAR models are more appropriate and suggest a rule to choose between the two models in practice. We provide an extended analysis of these issues using a simulation study. This allows us to (1) directly investigate the relative performance of AR and VAR models in typical psychological applications, (2) show how the relative performance depends both on n and characteristics of the true model, (3) quantify the uncertainty in selecting between the two models, and (4) assess the relative performance of different model selection strategies. We thereby provide a more complete picture for applied researchers about when the VAR model is appropriate in typical psychological applications, and how to select between AR and VAR models in practice.

Citation: Dablander F, Ryan O, Haslbeck JMB (2020) Choosing between AR(1) and VAR(1) models in typical psychological applications. PLoS ONE 15(10): e0240730. https://doi.org/10.1371/journal.pone.0240730

Editor: Miguel Angel Sánchez Granero, Universidad de Almeria, SPAIN

Received: February 7, 2020; Accepted: October 1, 2020; Published: October 29, 2020

Copyright: © 2020 Dablander et al. This is an open access article distributed under the terms of the Creative Commons Attribution License , which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: All data are available from the Github archive at https://github.com/jmbh/ARVAR .

Funding: Fabian Dablander was supported by the NWO Vici grant 016.Vici.170.083. Oisín Ryan was supported by a grant from the Netherlands Organisation for Scientific Research (NWO; Onderzoekstalent Grant 406-15-128). Jonas Haslbeck was supported by the Eropean Research Council Consolidator Grant no. 647209.

Competing interests: The authors have declared that no competing interests exist.

Introduction

Time series of individual subjects have become a common data type in psychological research since collecting them has become feasible due to the ubiquity of mobile devices. First-order Vector Autoregressive (VAR) models, which predict each variable by all variables including itself at the previous time point, are a natural starting point for the analysis of dependencies across time in such data and are already used extensively in applied research [ 1 – 5 ].

A key question that arises when using these models is: how reliable are the estimates of the single-subject VAR model, given the typically short time series in psychological research (i.e., n ∈ [30, 200])? To be more precise, we would like to know how large the estimation error is in this setting. Estimation error is defined as the distance between the estimated parameters and the parameters in the true model, assuming that the true model has the same parametric form as the estimated model. If estimation error is large, it might be possible to obtain a smaller estimation error by choosing a simpler model, even though it is less plausible than the more complex model [ 6 ]. A possible simpler model in this setting is the first-order Autoregressive (AR) model, in which each variable is predicted only by itself at the previous time point. While the AR model introduces a strong bias by setting all interactions between variables to zero, it can have a lower estimation error when the number of available observations is small. When analyzing time series in psychological research it is therefore important to know (a) in which settings the AR or the VAR model has a lower estimation error, and (b) how to choose between the two models in practice.

Bulteel et al. [ 7 ] identified these important and timely questions, and offered answers to both. They investigated question (a) regarding the relative performance of AR and VAR models by selecting three empirical time series data sets, each consisting of a number of individual time series with the same data structure. For each of these data sets, they approximate the out-of-sample prediction error with out-of-bag cross-validation error for both the AR and the VAR model and their mixed model versions. The authors make a valuable contribution by assessing which of the many cross-validation schemes available for time series approximates prediction error best in this context. Using the approximated prediction error obtained via cross-validation, they find that the prediction error for AR is smaller than for VAR, and that the prediction error of mixed AR and mixed VAR is similar. In a last step, they link prediction and estimation error by stating that “[…] the number of observations T [here n ] that is needed for the VAR to become better than the AR is the same for the prediction MSE [mean squared error] as well as for the parameter accuracy [estimation error]” [ 7 , p. 10]. Although the latter statement implies that the estimation error of mixed AR and mixed VAR models are similar, Bulteel et al. [ 7 ] conclude that “[…] it is not meaningful to analyze the presented typical applications with a VAR model” (p. 14) when discussing both mixed effects (i.e., multilevel models with random effects) and single-subject models.

Using their statement about the link between prediction error and estimation error together with a preference towards parsimony, Bulteel et al. [ 7 ] also offer an answer to question (b) on how to choose between the AR and VAR models in practice: they suggest using the “1 Standard Error Rule”, according to which one should select the AR model if its prediction error is not more than one standard error above the prediction error of the VAR model, and select the model with lowest prediction error otherwise [ 8 , p. 244].

In this paper, we provide an extended analysis of the problems studied by Bulteel et al. [ 7 ]. First, regarding question (a) on the relative performance of the AR and VAR models: when the goal is to determine the estimation error in a given setting, one can obtain it directly with a simulation study. A simulation study allows for a more extensive analysis of this problem for three reasons. First, we do not need to make any claim about the relation between prediction error and estimation error, which—as we will show—turns out to be non-trivial. Second, in a simulation study we can average over sampling variance which allows us to compute the expected value of estimation (and prediction) error. While the approach of Bulteel et al. [ 7 ] in using three empirical datasets has the benefit of ensuring the models considered mirror data from psychological applications, these empirical datasets are naturally subject to sampling variation. And third, a simulation study allows us to map out the space of plausible VAR models and base our conclusions on this large set of VAR models instead of the VAR models estimated from the three data sets used by Bulteel et al. [ 7 ]. We perform such a simulation study, which allows us to give a direct answer to the question of how large the estimation errors of AR and VAR models are in typical psychological applications.

Regarding question (b) on choosing between AR and VAR models in practice, Bulteel et al. [ 7 ] base their “1 Standard Error Rule” (1SER) on the idea that the n at which the estimation errors of the AR and VAR models are equal is (approximately) the same n at which the prediction errors of those models are equal, combined with a preference towards the more parsimonious model. While the 1SER is used as a heuristic in the statistical learning literature [ 8 ], it is not clear whether this heuristic would perform better in the present problem than simply selecting the model with the lowest prediction error. We show that when choosing between AR and VAR models, the n at which the prediction errors become equal is not necessarily the same as the n at which estimation errors become equal: in fact, there is a substantial degree of variation in how the prediction and estimation errors of both models cross. Using the relationship between estimation and prediction error we are able to show via simulation when the 1SER is expected to perform better than selecting the model with lowest prediction error. This extended analysis of the problem studied by Bulteel et al. [ 7 ] provides a more complete picture for applied researchers about when the VAR model is appropriate in typical psychological applications, and how to select between AR and VAR models in practice.

When does VAR outperform AR?

In this section we report a simulation study which directly answers the question of how large the estimation errors of AR and VAR models are in typical psychological applications. This allows the reader to get an idea of how many observations n e one needs, on average, for the VAR model to outperform the AR model. In addition, we will decompose the variance around those averages in sampling variation and variation due to differences in the VAR parameter matrix Φ . Finally, explaining the latter type of variation allows us to obtain n e conditioned on characteristics of Φ . The analysis code for the simulation study is available from https://github.com/jmbh/ARVAR .

Simulation setup

Since the AR model is nested under the more complex VAR model, we focus solely on the VAR as the true data-generating model. To obtain realistic VAR models, we use the following approach: first, we estimate a mixed VAR model to the “MindMaastricht” data [ 9 ], which consists of 52 individual time series with on average n = 41 measurements on p = 6 variables, and is the only publicly available data set used by Bulteel et al. [ 7 ]. In a second step, we sample stationary VAR models with a diagonal error covariance matrix from this mixed model.

We expect that the estimation (and prediction) errors of the AR and VAR model depend not only on the number of observations n , but also on the characteristics of the underlying p × p VAR model matrix Φ . We therefore stratify the sampling process from the mixed model by two characteristics of Φ . This procedure allows us to obtain a better picture of how the performance of AR and VAR may differ depending on the characteristics of the data generating model.

Ideally, we would stratify by sampling a fully crossed grid of D and O values. However, this is not possible since some combinations have an extremely small probability: For example, if a matrix has auto-regressive parameters close to one, it is unlikely to describe a stationary process if it also contains high positively-valued cross-lagged parameters. We therefore adopt the following approach: we divided the D - O -space in a grid by dividing each dimension into 15 equally spaced intervals (see S1 Fig ). We then include only those cells in the design in which at least one VAR model has been sampled. This procedure returned 74 non-empty cells. We then sample those 74 cells until each of them contains 100 VAR models. We keep the cell size constant to render the results comparable across cells (see Supporting Information for a detailed description of this procedure).

This procedure returns a set of 74 × 100 = 7400 VAR models that includes essentially any stationary VAR model with p = 6 variables, and allows us to describe each model in the dimensions O and D . For each of these VAR models, we generate 100 independent time series, each with n = 500 observations and with a burn-in period of n burn = 100. We then estimate both the AR and the VAR model on the first n = {8, 9, …, 499, 500} observations of those time series. This yields a simulation study with 7400 × 493 (parameters × sample size) conditions, and for each of those conditions we have 100 replications. For each model, and each n , we compute the expected estimation error for both the AR model (EE AR ) and the VAR model (EE VAR ) model by averaging over the 100 replications. This means that while EE AR and EE VAR have different values depending on n and the underlying model, we have averaged over the sampling variation.

Simulation results

The simulation described above allows us to investigate the relative performance of AR and VAR models across different samples, sample sizes, and data-generating models. We define the estimation error as the mean squared error of the estimated parameters to the true parameters, and quantify the relative performance with two measures: the difference between the estimation errors of the AR and VAR models at a particular sample size, EE Diff = EE AR − EE VAR ; and, n e , the sample size at which the VAR model outperforms the AR model (EE AR > EE VAR ). In the following we examine the mean and variance of EE Diff and subsequently study n e and its dependence on the characteristics of the true VAR model.

Fig 1(a) shows the mean and standard deviation of EE Diff as a function of n , across all 7400 VAR models and 100 replications. The dashed line at EE Diff = 0 indicates the point at which the estimation errors of the two models are equal. Below that line, the AR model performs better, that is, its parameter estimates are closer to the parameters of the true VAR model than the parameter estimates of the VAR model. We see that, across all models, we obtain a median n e = 89. Note that, out of all 740,000 simulated data sets, in only 23 cases the estimation error curves did not yet cross with an n of 500. Notably, the variance around the difference in estimation error is substantial for all n . In the following we decompose this variance in variance due to sampling error, and variance due to differences in VAR matrices.

• PPT PowerPoint slide
• PNG larger image
• TIFF original image

Panel (a) shows EE Diff averaged over replications and models, and the band shows the standard deviation over replications and models; panel (b) shows EE Diff for each model averaged across replications; and panel (c) shows the EE Diff averaged over replications for three specific models, and the bands show the standard deviation across 100 replications (sampling variation).

https://doi.org/10.1371/journal.pone.0240730.g001

Panel (b) of Fig 1 displays the mean EE Diff for each of the 7400 VAR models, averaged across 100 replications. We see that the lines differ considerably and that n e substantially depends on the characteristics of the true VAR model. This shows that one cannot expect reliable recommendations with respect to n e that ignore the characteristics of the generating model. To illustrate the extent of the sampling variation of the models, we have chosen three particular VAR models (see coloured lines). Fig 1(c) shows that they exhibit considerable sampling variation. Note that, as the variance in (b) is due to differences in mean performance across VAR models, it does not decrease with n . In contrast, the variance in (c) depends on n as it pertains to the sampling variance of a single VAR model, which decreases with the square root of the number of observations. While the mean EE Diff (shown in Fig 1(a) ) gives a clear answer to the question of which n is required for the VAR model to outperform the AR model on average , both types of variation (see Fig 1(b) and 1(c) ) show that for any particular VAR model it is difficult to determine which model performs better with the sample sizes typically available in psychological applications. However, we see that the sampling variation across replications is smaller than the variation across VAR models for most n . This means that if one has information about the parameters of the data-generating model, one can make much more precise statements about the sample size necessary for the VAR model to outperform the AR model.

The large degree of variation around EE Diff also highlights the potential pitfalls of generalizing the findings of Bulteel et al. [ 7 ] beyond the empirical data sets, which consist of 28, 52, and 95 individual time-series with an average number of 41, 70 and 70 time points, analyzed by the original authors. This is because (i) it is unlikely that their (in total) 175 time series appropriately represent the population of all plausible VAR matrices, (ii) their sample is subject to a substantial amount of sampling variation, and (iii) the absence of systematic variations of n does not allow a comprehensive answer to how relative performance relates to sample sizes.

Above we suggested that the relative performance of AR and VAR models (quantified by EE Diff ) depends on the characteristics D and O of the true VAR parameter matrix. In Fig 2(a) we show the median (across models in cells) n at which the estimation error of VAR becomes smaller than the estimation of AR (i.e., EE Diff > 0), depending on the characteristics D and O . We see that the larger the average off-diagonal elements O , the lower the n at which VAR outperforms AR. This is what one would expect: when O is small (as indicated by the lowest rows of cells in Fig 2(a) ), the true VAR model is actually very close to an AR model. In such a situation, the bias introduced by the AR model by setting the off-diagonal elements to zero leads to a relatively small estimation error. This trade-off between a simple model with high bias but low variance and a more complex model with low bias but high variance is well-known in the statistical literature as the bias-variance trade-off [ 8 ]. It therefore takes a considerable amount of observations until the variance of the VAR estimates becomes small enough for it to outperform the AR model. When O is large (indicated by the upper rows of cells), the bias of the AR model leads to a comparatively larger estimation error. Finally, we can also see that the size of the diagonal elements D is not as critical in determining n e as the size of the off-diagonal elements: Picking any row of cells in Fig 2(a) , we can see that there is only a very small variation across columns, with larger D values appearing to lead to very slight decreases in n e in general. Note that the O characteristic also largely explains the vertical variation of the estimation error curves shown in Fig 1(b) : the curves on top (small n e ) have high O , while the curves at the bottom (large n e ) have low O . Fig 2(b) collapses across these values and illustrates the sampling distribution of n e , taking into account the likelihood of any particular VAR matrix (as specified by the mixed model estimated from the “MindMaastricht” data).

Left: n e , the n at which estimation error becomes lower for the VAR than for the AR model, as a function of D and O . Right: Sampling distribution of n e , the n at which the expected estimation error of the VAR model becomes lower than the expected estimation error of the AR model. The dashed line indicates the median of 89.

https://doi.org/10.1371/journal.pone.0240730.g002

In summary, we used a simulation study to investigate the relative performance of AR and VAR models in a much larger space of plausible data-generating VAR models in psychological applications than considered by Bulteel et al. [ 7 ]. Next to investigating the average relative performance as a function of n , we also looked into the variation around averages. We showed that there is substantial variation both due to sampling error and differences in VAR matrices, which means that for a particular time series at hand it is difficult to select between AR and VAR with the n available in typical psychological applications. Finally, we found that the size of the off-diagonal elements influences the relative performance of the VAR model more strongly than the size of the diagonal elements.

Choosing between VAR and AR based on prediction error

In the previous section, we directly investigated the estimation errors of the AR and the VAR model in typical psychological applications and showed that the n at which VAR becomes better than AR depends substantially on the characteristics of the true model. In practice, the true model is unknown, so we can neither look up the n at which VAR outperforms AR in the above simulation study, nor can we compute the estimation error on the data at hand. Thus, to select between these models in practice, we may choose to use the prediction error which we can approximate using the data at hand, for instance by using a cross-validation scheme as suggested by Bulteel et al. [ 7 ]. However, since we are interested in estimation error, we require a link between prediction error and estimation error. In the remainder of this section we investigate this link and discuss the implications of this link for the model selection strategy suggested by Bulteel et al. [ 7 ], who use the “1 Standard Error Rule” (1SER) to select the model with lowest estimation error. Finally, we use our simulation study from above to directly compare the performance of the 1SER with model selection based only on the minimum prediction error.

The relation between prediction error and estimation error

Bulteel et al. [ 7 ] suggest that the link between prediction error and estimation error is relatively straightforward: “[…] the number of observations T [here n ] that is needed for the VAR to become better than the AR is the same for the prediction MSE [mean squared error] as well as for the parameter accuracy [estimation error]” [ 7 , p. 10]. More formally, this claim states that if n e is the number of observations at which the estimation errors of the AR and VAR model are equal, and if n p is the number of observation at which the prediction errors of the AR and VAR model are equal, and n gap = n e − n p , then n gap = 0. Bulteel et al. [ 7 ] do not specify the exact conditions under which this statement should hold, and elsewhere in the text suggest that this should be considered an approximate rather than an exact relationship. If this relationship were indeed approximate, it would still be interesting to study in which settings n gap > 0 or n gap < 0, as this bears on model selection, and so we will focus our investigation on quantifying n gap and investigating any potential systematic deviations from zero through simulation. Clearly, it would be unreasonable to expect that n gap = 0 for any data set, since the observations in a given data set are subject to sampling error. We therefore interpret the statement of Bulteel et al. [ 7 ] as a statement about the expectation over errors of any given VAR model.

Assessing n gap through simulation

We now use the results of the simulation study from the previous section to check whether indeed n gap = 0 on average for all VAR models. To compute prediction error, we generate a test-set time series consisting of n test = 2000 observations (using a burn-in period of n burn = 100) for each of the 7400 VAR models described in the previous section. For each of the 100 replications of model and sample size condition, we average over the prediction errors which are obtained when estimated model parameters are evaluated on the test set. This is the out-of-sample prediction error (i.e., the expected generalization error) that Bulteel et al. [ 7 ] approximate with out-of-bag cross-validation error. We define prediction error as the mean squared error (MSE) of the predicted values relative to the true values in the test data set.

Fig 3 shows the estimation (solid lines) and prediction (dashed lines) errors for both the AR (black lines) and VAR (red lines) models as a function of n , averaged across the replications, for model A with D = 0.068 and O = 0.092 (left panel) and model B with D = 0.337 and O = 0.051 (right panel). For model A, we see that n gap < 0, which shows that n gap = 0 for all VAR models is incorrect. What consequences does this gap have for model selection? The negative gap implies that if the prediction errors for the AR and VAR model are the same, the VAR model should be selected, because its estimation error is smaller. In contrast, for model B we observe n gap > 0. In this situation, if the prediction errors are equal, one should select the AR model because it incurs smaller estimation error. Clearly, n gap differs between the two models, and this difference matters for model selection.

Scaled Mean Squared Error (MSE) of estimation (solid lines) and prediction errors (dashed lines) for both the AR (black lines) and VAR (red lines) models as a function of n , separately for model A with D = 0.068 and O = 0.092 (left panel) and model B with D = 0.337 and O = 0.051 (right panel). The red and green shaded area indicates the median n gap , and the grey shaded area shows the 20% and 80% quantiles across the 100 replications per model.

https://doi.org/10.1371/journal.pone.0240730.g003

So far we only investigated n gap for two individual VAR models. Fig 4(a) shows the distribution of the expected n gap across all VAR models, computed by averaging over 100 replications. Note that for 31 out of 7400 models the curves of prediction errors and estimation errors did not cross within n ∈ {8, 9, …, 499, 500}. The results in Fig 4 are therefore computed on 7369 models.

Panel (a) displays the distribution of the expected n gap across all 7369 VAR models, computed by averaging over 100 replications, and weighted by the probability defined by the original mixed model. Panel (b) shows the distribution of non-zero EE comp across all n , 7369 VAR models, averaged across replications and weighted by the probability defined by the original mixed model.

https://doi.org/10.1371/journal.pone.0240730.g004

Each of the data points in the histogram in Fig 4(a) corresponds to the expected n gap of one of the 7369 models. We see that the expected n gap has a right skewed distribution with a mode at zero. This allows us to make a precise statement regarding the crossing of estimation and prediction errors described above: while the most common value of n gap is zero, most expected n gap are not zero. In fact, n gap shows substantial variation across different VAR models. Explaining the variance of n gap is interesting, because n gap has direct consequences for model selection. If we can relate the n gap to characteristics of the Φ matrix, it is possible to make more specific statements with respect to when to apply a bias towards the AR or VAR model, when the prediction errors are the same or very similar. Note that such a function from Φ to n gap must exist, because the only way the 7400 models differ is in their entries of the VAR parameter matrix Φ . However, this function may be very complicated. For example, the Pearson correlation of n gap with D and O are 0.21 and −0.02, respectively. Predicting n gap by D and O including the interaction term with linear regression achieves R 2 = 0.048. This shows that a simple linear model including D and O is not sufficient to describe the relationship between n gap and Φ . Future research could look into better approximations of this relationship. If successful, one could build new model selection strategies on reliable predictions of n gap from empirical data.

Performance of the “1 Standard Error Rule”

Bulteel et al. [ 7 ] propose, in the words of Hastie et al., to “[…] choose the most parsimonious model whose error is no more than one standard error above the error of the best model.” [ 8 ], p. 244]. This model selection criteria is known as the “1 Standard Error Rule” (1SER) and is suggested by Hastie and colleagues as a method of choosing a model with the minimal out-of-sample prediction error (which is typically unknown), on the basis of out-of-bag prediction error (acquired with cross-validation techniques).

Making inferences from prediction error to estimation error requires a link between the two. Bulteel et al. [ 7 ] provide this link by suggesting that n gap = 0 (or n gap ≈ 0). However, they do not provide justification for why the 1SER should outperform simply selecting the model with the lowest prediction error. Above we showed that n gap = 0 does not hold for all VAR models. In fact, it is this result that explains why the 1SER can perform better than selecting the model with the lowest prediction error. Specifically, this is the case when n gap > 0, which characterizes the situation that the prediction error for VAR is lower than for AR while at the same time the estimation error of VAR is higher than for AR. In such a situation, a bias towards the AR model can be favorable. In contrast, if n gap < 0 and the prediction error of AR is lower than for VAR, even though the estimation error of VAR is lower than for AR, such a bias would be unfavorable. In the following, we assess the relative performance of the 1SER and simply selecting the model with lowest prediction error, both on average and as a function of n .

Fig 4(b) shows the distribution of non-zero EE comp across all 7400 VAR models, averaged over replications, and weighted by the probability given by the original mixed model. The only interesting cases when comparing model selection procedures are the cases in which they disagree. Therefore, we analyze only those cases for which EE comp ≠ 0. Note that for all but 2 of the 7400 models there is some n at which the two decision rules in question choose a different model. We find that using the 1SER is better in 50.1% of cases (where each case is weighted by the probability of the corresponding model). This would suggest that it makes essentially no difference whether we use the 1SER or select the model with lowest prediction error. However, these proportions average over the number of observations n and therefore cannot reveal differences in relative performance for different sample sizes.

Fig 5(a) shows EE comp as a function of n , averaged across all 7400 models. Because the VAR prediction error is huge for very small n , both model selection strategies choose the same model, resulting in EE comp = 0 for those n . However, from around n = 10 on until around n = 60, EE comp is substantially positive, indicating that the 1SER outperforms simply selecting the model with the lowest prediction error by a large margin. However, for n > 60 we see that EE comp approaches zero and then becomes slightly negative. The latter is also illustrated in panel (b), which displays the weighted proportion of models in which the 1SER is better (i.e., EE comp > 0). The explanation of this curve has three parts. First, n gap tends to be larger if the gap is located at a small n (Pearson correlation r = −0.15). If n gap is large (and therefore positive), the AR model has lower estimation error than the VAR model, even though the prediction errors are the same (compare Fig 5(b) ). In such situations, biasing model selection towards selecting the AR model is advantageous. Since the 1SER constitutes a bias towards the AR model, it performs better for small n . Second, this also explains why the 1SER performs worse than simply selecting the model with lowest prediction error for large n : here the gap is small (negative), indicating that if the prediction errors are the same, the VAR model performs better. Clearly, in such a situation, providing a bias towards AR is disadvantageous. Therefore, the 1SER performs worse. Finally, why does the curve get closer and closer to zero? The reason is that the standard error converges to zero with (the square root of) the number of observations, and therefore the probability that both rules select the same model approaches 1 as n goes to infinity.

Panel (a) displays EE comp averaged across 7400 models as a function of n (black line) and the standard deviation around the average (blue line). Panel (b) displays, for each n , the proportion of times that EE comp > 0 across 7400 models (i.e., the proportion of 1SER performing better).

https://doi.org/10.1371/journal.pone.0240730.g005

To summarize, we found that the 1SER is better than simply selecting the model with the lowest prediction error only in 50.1% of the cases in which the two rules do not select the same model. However, when looking at the relative performance as a function of n , we found that the 1SER is better than selecting the model with lowest prediction error until around n = 60, and worse above. Finally, we were able to explain the dependence of the relative performance on n with the fact that n gap is larger when it occurs at a smaller n . For applied researchers these results suggest that, for VAR models with p = 6 variables, the 1SER should be applied for n < 60.

In this paper we provided an extended analysis of the problem studied by Bulteel et al. [ 7 ] by using a simulation study to (a) map out the relative performance of AR and VAR models in typical psychological applications as a function of the number of observations n , and (b) investigate how to choose between AR and VAR models in practice. We found that, averaged over all models considered in our simulation, the VAR model outperforms the AR model for n > 89 observations in terms of estimation error. In addition, we show that and explain why the 1SE rule proposed by Bulteel et al. [ 7 ] performs better than selecting the model with the lowest prediction error when n is small.

Next to the average estimation errors of AR and VAR models, we also investigated the variance around those averages. We decomposed this variance in variance due to different true VAR models, and variance due to sampling. The variance across different VAR models showed that the relative performance, that is, the n at which VAR becomes better than AR ( n e ) depends on the characteristics of the true VAR parameter matrix Φ . For example, if the true VAR model is very close to an AR model, it takes more observations until the VAR model outperforms the AR model. This shows that one cannot expect reliable recommendations with respect to n e that ignore the characteristics of the generating model: n e critically depends on the size of the off-diagonal elements present in the data-generating model. The size of the sampling variation also indicates that, for many of the considered sample sizes, whether the VAR or AR model will have lower estimation error largely depends on the specific sample at hand. This implies that it is difficult to select the model with lowest estimation error with the sample sizes available in typical psychological applications.

The second question we investigated was: how should one choose between the AR and VAR model for a given data set? Bulteel et al. [ 7 ] suggest that, for any VAR model, the n at which the prediction errors of both models are equal is, in expectation, (approximately) the same n at which their estimation errors are equal (i.e., n gap ≈ 0). Combining this claim with a preference towards the more parsimonious AR model, they proposed using the “1 Standard Error Rule”, according to which one should select the AR model if its prediction error is not more than one standard error above the prediction error of the VAR model, and choose the model with lowest prediction error otherwise. We showed that the expected n gap varies as a function of the parameter matrix of the true VAR model. Using the relationship between estimation and prediction error we were able to explain when the 1SER is expected to perform better than selecting the model with lowest prediction error. In addition, we showed via simulation that the 1SER performs better than selecting the model with the lowest prediction error for n < 60, in cases where those decision rules select conflicting models. Our simulations also showed that as n → ∞, both decision rules converge to selecting the same model. This means that there is a relatively small range of sample sizes in which these decision rules lead to contradictory model selections for a given data-generating system. We recommend that researcher wishing to use prediction error to choose between these models examine both the 1SER and lowest prediction error rules, and in cases of conflict between the two, use the 1SER for low ( n < 60) sample sizes.

The relative performance of the AR and VAR model shown in our simulations can be understood in terms of the bias-variance trade-off. Because the AR model sets all off-diagonal elements to zero, it has a bias that is constant and independent of n . In contrast, the VAR model has a bias of zero, since the true model is a VAR model. This is why a VAR model will always perform better than (or at least as good as, if the all off-diagonal elements of the true VAR model are zero) an AR model as n → ∞. However, for finite sample sizes the variance of the estimates of the two models are different: while both variances converge to zero as n → ∞, for finite samples the variance of VAR parameters is much larger than the variance of AR parameters, especially for small n . This allows for the situation that the biased simpler model is showing a smaller error, even though the true model is in the class of the more complex model. This trade-off between bias and variance also explains the relative performance of AR and VAR models: From Fig 3 we saw that for small n , the variance of the VAR estimates is so large that the error is larger than the error of the AR model, despite the bias of the AR model. However, with increasing n , the variance of the estimates of both models approaches zero. This means that the larger n , the more the bias of the AR model contributes to its error. Thus, at some n the error of the VAR model becomes smaller than the error of the AR model. We agree with Bulteel et al. [ 7 ] that the fact that a simple (and possibly implausible) model can outperform a complex (and more plausible) model, even though the true model is in the class of the more complex model, is underappreciated in the psychological literature.

An interesting question we did not discuss in our paper is: which model should we choose if the AR and VAR models have equal estimation error? Since we defined the quality of a model by its estimation error, we could simply pick one of the two models at random. However, their model parameters are likely to be very different. The estimation error of the AR model comes mostly from setting off-diagonal elements incorrectly to zero, while the estimation error of the VAR model comes mostly from incorrectly estimating off-diagonal elements. In terms of the types of errors produced by the two models, the AR model will almost exclusively produce false negatives, while the VAR model will produce almost exclusively false positives. A specification of the cost of false positives/negatives in a given analysis may allow to choose between models when the estimation errors are the same or very similar. For example, in an exploratory analysis one might accept more false positives in order to avoid false negatives.

Throughout the paper we compared the AR model to the VAR model. However, we believe that it is unnecessarily restrictive to choose only between those extremes (all off-diagonal elements zero vs. all off-diagonal elements nonzero). The AR model, by imposing independence between processes, presents a theoretically implausible model for many psychological processes. Applied researchers who estimate the VAR model may be primarily interested in the recovery of cross-lagged effects rather than auto-regressive parameters, for example to determine which processes are dependent on one another (as evidenced by frequent discussions of Granger causality [ 11 ] In such settings, one could estimate VAR models with a constraint that limits the number of nonzero parameters or penalizes their size [ 12 , 13 ]. This would allow the recovery of large off-diagonal elements without the high variance of estimates in the standard VAR model. Similarly, one could estimate a VAR model and, instead of comparing it to an AR model and thus testing the nullity of the off-diagonal elements jointly, test the nullity of the off-diagonal elements of the VAR matrix individually. Further investigation of these alternatives in future studies would provide a more complete picture to applied researchers.

It is important to keep the following limitations of our simulation study in mind. First, we claimed that the 7400 models we sampled from the mixed model obtained from the “MindMaastricht” data represent typical applications in psychology. One could argue that there are sets of VAR models that are plausible in psychological applications that are not included in our set of models. While this is a theoretical possibility, we consider this extremely unlikely, since we heavily sampled the mixed model stratified by O and D . Any VAR model that is not similar to a model in our set of considered VAR models is therefore most likely non-stationary. When presenting our results we weighted all models by their frequency given the estimated mixed model in order to avoid giving too much weight to unusual VAR models. This means that it could be that the weighting obtained from the mixed model does not represent the frequency of VAR models in psychological applications well. While we consider this unlikely, we also used a uniform weighting across VAR models as a robustness check which left all main conclusions unchanged. A second limitation is that we only considered VAR models with p = 6 variables. While this is not a shortcoming compared to Bulteel et al. [ 7 ] who use VAR models with 6, 6, and 8 variables, the results shown in the present paper would likely change when considering more or less than six variables. Specifically, we expect that the n at which VAR outperforms AR becomes larger when more variables are included in the model, and smaller when less variables are included. This change may be nonlinear in nature: As we add variables to the model, we would expect the variance of the VAR model to grow much quicker than the variance of the AR model, since in the former case we need to estimate p 2 parameters, and in the latter only p . However, the bias of the AR model also grows with each new variable added, with p 2 − p elements set to zero in each case, and so again, this will largely depend on the data-generating system at hand. Similarly, we would expect that for models with more variables the 1SER outperforms selecting the model with lowest prediction error for sample sizes larger than 60. While the exact values will change for larger p , we expect that the general relationships between n , O , and D extend to any number of variables p .

Although Bulteel et al. [ 7 ] also consider mixed VAR and AR models, in the simulation studies presented above we focus exclusively on single-subject time-series for simplicity. Mixed models can be seen as a form of regularization, in which individual parameter estimates are shrunk towards the group-level mean if the number of observations n is small. One would expect that for small n , the use of mixed models would improve the estimation and prediction errors of both models, which is also what Bulteel et al. [ 7 ] report in their results. Indeed, mixed models are expected to improve the performance of VAR methods relative to AR, and thus may be a solution to the relatively poor performance of the VAR model we observe in sample sizes realistic for psychological applications. The reason is that the differential performance of AR and VAR models can be understood in terms of a bias-variance trade-off, where AR models are biased but have lower variance than VAR methods. The use of mixed VAR models should decrease this variance through shrinkage in small n settings [ 14 , 15 ]. The precise effect of using mixed models depends on the variance of parameters across individuals; however, we do not expect the general pattern of results reported here to change when moving from single-subject to mixed settings.

Future research could extend the analysis shown here to VAR models with less than or greater than six variables, which would allow to generalize the simulation results to more situations encountered in psychological applications. Another interesting avenue for future research would be to investigate the link between n gap and the VAR parameter matrix Φ . Since n gap has direct implications for model selection, such a link could possibly be used to construct improved model selection procedures. It would be useful to extend the simulation study in this paper to constrained estimation such as the LASSO, especially since those methods are already applied in practice [ 16 ]. Finally, it would be useful to study the performance of mixed VAR models in a simulation setting, and perhaps compare this approach to alternative methods of using group-level information in individual time-series analysis, such as GIMME, an approach originally developed for the analysis of neuroimaging data [ 17 ]. Early simulation studies have assessed the performance of mixed AR models in recovering fixed effects using Bayesian estimation techniques [ 18 ], but these analyses have yet to be extended to mixed VAR models or the recovery of individual-specific random effects.

To sum up, we used simulations to study the relative performance of AR and VAR models in settings typical for psychological applications. We showed that, on average, we need sample sizes approaching n = 89 for single-subject VAR models to outperform AR models. While this may seem like a relatively large sample size requirement, such longer time series are becoming more common in psychological research [ 19 , 20 ] Decomposing this variance showed that (i) one cannot expect reliable statements with respect to the relative performance of the AR and VAR models that ignore the characteristics of the generating model, and (ii) that choosing reliably between AR and VAR models is difficult for most sample sizes typically available in psychological research. Finally, we provided a theoretical explanation for when the “1 Standard Error Rule” outperforms simply selecting the model with lowest prediction error, and showed that the 1SER performs better when n is small.

Supporting information

S1 fig. d and o values for the initially sampled 10000 var models..

https://doi.org/10.1371/journal.pone.0240730.s001

Acknowledgments

We would like to thank Don van den Bergh, Riet van Bork, Denny Borsboom, Max Hinne, Lourens Waldorp, and two anonymous reviewers for their helpful comments on earlier versions of this paper.

• View Article
• PubMed/NCBI
• Google Scholar
• 8. Hastie T, Tibshirani R, Friedman J. The Elements of Statistical Learning. Springer; 2009.
• 10. Hamilton JD. Time series analysis. vol. 2. Princeton, NJ: Princeton University Press; 1994.
• 13. Hastie T, Tibshirani R, Wainwright M. Statistical learning with sparsity: the lasso and generalizations. CRC Press; 2015.

arXiv's Accessibility Forum starts next month!

Help | Advanced Search

Quantitative Finance > Statistical Finance

Title: large investment model.

Abstract: Traditional quantitative investment research is encountering diminishing returns alongside rising labor and time costs. To overcome these challenges, we introduce the Large Investment Model (LIM), a novel research paradigm designed to enhance both performance and efficiency at scale. LIM employs end-to-end learning and universal modeling to create an upstream foundation model capable of autonomously learning comprehensive signal patterns from diverse financial data spanning multiple exchanges, instruments, and frequencies. These "global patterns" are subsequently transferred to downstream strategy modeling, optimizing performance for specific tasks. We detail the system architecture design of LIM, address the technical challenges inherent in this approach, and outline potential directions for future research. The advantages of LIM are demonstrated through a series of numerical experiments on cross-instrument prediction for commodity futures trading, leveraging insights from stock markets.

Submission history

Access paper:.

• HTML (experimental)
• Other Formats

References & Citations

• Google Scholar
• Semantic Scholar

BibTeX formatted citation

Bibliographic and Citation Tools

Code, data and media associated with this article, recommenders and search tools.

• Institution

arXivLabs: experimental projects with community collaborators

arXivLabs is a framework that allows collaborators to develop and share new arXiv features directly on our website.

Both individuals and organizations that work with arXivLabs have embraced and accepted our values of openness, community, excellence, and user data privacy. arXiv is committed to these values and only works with partners that adhere to them.

Have an idea for a project that will add value for arXiv's community? Learn more about arXivLabs .

Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

• View all journals
• Explore content
• About the journal
• Publish with us
• Sign up for alerts
• Open access
• Published: 24 August 2024

Research on food security issues considering changes in rainfall

• SiMan Jiang 1 ,
• Shuyue Chen 2 ,
• Qiqi Xiao 2 &
• Zhong Fang 2  

Scientific Reports volume  14 , Article number:  19698 ( 2024 ) Cite this article

Metrics details

• Environmental economics
• Environmental social sciences

Ensuring food security is not only vital to the adequate supply of food in the world, but also the key to the economic development and social stability of various countries. Based on the panel data of 29 provinces in China from 2016 to 2020, this paper selects the number of foodborne diseases patients and iodine deficiency disease patients as reference objects, uses stunting rate of children under 5 years old, malnutrition rate of children under 5 years old, obesity rate of children under 5 years old, and newborn visit rate to measure improving nutrition, proposes Meta Entropy Two-Stage Dynamic Direction Distance Function (DDF) Under an exogenous Data Envelopment Analysis (DEA) model to measure the efficiency of hunger eradication, food security, and improving nutrition under the influence of exogenous variable rainfall. The research results indicate that the sustainability of China’s agricultural economy is insufficient, and the focus of attention should be different in different stages. In addition, the average efficiency of the three regions generally shows a decreasing level in the eastern, western, and central regions. In order to improve China's ability to guarantee food security, we must continue to strengthen the construction of agricultural infrastructure, increase policy support for green agricultural production, promote the diversification of agricultural production, and enrich people’s agricultural product consumption varieties.

Similar content being viewed by others

Using household economic survey data to assess food expenditure patterns and trends in a high-income country with notable health inequities

Global food nutrients analysis reveals alarming gaps and daunting challenges

The COVID-19 crisis will exacerbate maternal and child undernutrition and child mortality in low- and middle-income countries

Introduction.

The issue of food security is not only related to the livelihoods of countries but also to global development. In 2015, the United Nations General Assembly adopted the 2030 Agenda for Sustainable Development, putting forward 17 Sustainable Development Goals (SDGs), of which the second goal (SDG2) focuses on food security and commits to eradicating hunger, achieving food security, improving nutrition, and promoting sustainable agriculture by 2030, also known as the “Zero Hunger” goal. Food security is an important cornerstone and key issue for global sustainable development. Currently, Food production has made significant progress globally in eradicating hunger, food insecurity, and malnutrition. However, many people are still facing hunger and malnutrition due to the impact of various factors such as extreme weather, global COVID-19 pandemic, and geopolitical conflicts in recent years. In addition, the loss of arable land and urban expansion have adversely affected agricultural land and put enormous pressure on preventing the degradation of ecosystem service functions and adapting to climate change, which brings new uncertainties to global food security and new challenges to food security in China. To cope with the uncertainty of global food security, food security in China has become even more important.

Having entered a new stage of development, China has made significant achievements in food security. In the face of the global food crisis, China’s food production has achieved a good harvest for 19 consecutive years, the total food output has remained above 650 million tons for 8 consecutive years, the self-sufficiency rate of food rations has exceeded 100 percent and that of cereal foods has exceeded 95 percent, with the per capita food possession at approximately 480 kg, which is higher than the internationally acknowledged food security line of 400 kg, and China has achieved basic self-sufficiency in cereals and absolute security in rations. Using 9 percent of the world's arable land and 6 percent of its freshwater resources, China has been able to feed nearly 20 percent of its population, making a historic transition from hunger to subsistence to well-being. However, after a long period of sustained improvement, China's food security situation was reversed in 2015 due to multiple challenges, including agricultural environmental pollution and intensifying climate change. Climate change will have a negative impact on food production, which will increase the price of agricultural products and subsequently increase China’s food imports, which in turn will affect China’s level of food self-sufficiency. Currently, for every 0.1 °C increase in temperature, China’s yield per unit area of the three major food crops will decrease by about 2.6 per cent, and just a 1 per cent increase in precipitation will increase the yield per unit area by 0.4 per cent. In recent years, climate change has led to significant changes in China’s agroclimatic resources: From 1951 to 2021, the annual average surface temperature in China increased at a rate of 0.26 °C per decade; annual rainfall in China increased by an average of 4.9 mm per decade, showing a trend of “northern expansion of the rainfall belt”. The “double increase in water and heat” of climate change has led to significant changes in China’s agroclimatic resources, with the crop growing season lengthening by 1.8 days per decade. The impact of climate change on agricultural production is both negative and positive, but the negative impact of uneven rainfall and extreme weather on agriculture is significant and requires increased attention. The problem of uneven rainfall is reflected in the redistribution of global rainfall, with increased rainfall in some areas causing flooding and damage to crop roots and soil structure, thus reducing food production; reduced rainfall in some other areas leads to drought, which affects crop growth and development, and likewise reduces food production. Droughts used to exist in the northern regions of China, but seasonal droughts are now occurring in many southern regions, especially at critical times of crop growth, leading to significant reductions in crop yields. At present, China's food security still faces many risks and challenges, with new problems in both production and consumption, such as the contradiction between the basic balance of food supply and demand and structural scarcity, the contradiction between food production methods and the upgrading of food demand, and the contradiction between the international food market linkage and the volatility of domestic food prices, which has resulted in a potentially further deterioration of the food security situation. These food insecurity trends will ultimately increase the risk of malnutrition and further affect the quality of diets, affecting people’s health in different ways. Currently, with less than a decade to go before the achievement of the 2030 SDGs, the global food security situation is still spiraling downwards. Therefore, food security should always be a matter of crisis awareness.

Food security is affected by several factors, and rainfall is one of the major influences on food production. The regional impact of rainfall on production is complex and can have an impact on the total food production in China. Although the national rainfall has not shown a significant trend in the last 50 years, there are significant regional differences. In the scientific study of global change, there will be a long way to go to study the impact of rainfall changes on food production and food security in different regions of China and to propose effective countermeasures.

Literature review

The concept of food security was first officially introduced by the Food and Agriculture Organization of the United Nations (FAO) in 1974. It is defined as ensuring sufficient global availability of basic food supplies at all times, particularly in the case of natural disasters or other emergencies to prevent the exacerbation of food shortages, while steadily increasing food consumption in countries with low per capita intake to reduce production and price fluctuations. Make food security one of the basic rights of human life. This concept reflects people’s concerns about the occurrence of global food crisis at that time, recognizing that the decline of food supply plays a major role in promoting the expansion of hunger, while the instability of food prices caused by supply–demand imbalances exacerbates the severity of the hunger situation 1 . Although the early definition of food security primarily emphasized the quantity of food supply, namely the accessibility of food, and measures to address hunger mainly focused on expanding food production, there has been a growing recognition of the importance of food stability as a crucial aspect of food security 2 . As the world's economic situation evolves, people have gained a better understanding of food security, leading to an expanded conceptual framework. In 1982, the FAO revised the definition of food security to ensure sufficient food supply, stable food flows, and stable food sources for individuals or households. This new interpretation incorporates some micro considerations into the existing macro perspective, emphasizing the significance of balancing food supply and demand 3 . During the World Food Summit in 1996, the FAO updated the definition of food security to ensure that all individuals have physical and economic access to sufficient, nutritious and safe food at all times, and the effective utilization of these food nutrients, and defined four pillars of food security: availability, accessibility, utilization, and stability. In 2001, the FAO added the term “social” to the original definition of food security, which has become the most widely cited definition in current international food policies, that is, to ensure that all individuals at all times have physical, social and economic access to sufficient, nutritious and safe food to meet people's needs and preferences regarding food and promote people to lead positive and healthy lives.

Food security is closely related to people's lives, and it has always been the focus of academic attention. The existing research mainly analyzes the impact of resource endowment, climate change and government policy on food security, and then explores practical paths for various countries and regions to ensure food security in the future.

The literature primarily focuses on water resources, land resources, and human resources and other aspects to study the impact of resource endowment on food security. From the perspective of water resources, Kang et al. 4 summarized the evolution of irrigation water productivity in China over the past 60 years, studied the differences in food productivity under different planting patterns, fertilization levels, and irrigation water consumption, analyzed the current situation of water resources’ impact on food security and explored comprehensive measures to improve agricultural water use efficiency in the future; Chloe et al. 5 combining interviews and surveys from British farmers with the resilience theory to analyze the influencing factors of water scarcity risk and management strategies, found that farmers need to establish resilience by maintaining the buffer of water resources or increasing the availability of backup resources to minimize the negative impacts of water scarcity on food production and farmer’s economic income. From the perspective of land resources, Charoenratana and Shinohara 6 pointed out that land and its legal rights are crucial factors for farmer income and agricultural production, and sustainable food security can only be achieved if land is kept safe. Li et al. 7 indicated that while there has been a strong transition of cultivated land from non-staple food production to food production in the suburbs of Changchun after rapid urbanization, overall, the utilization diversity of suburban cultivated land in the black soil region of Northeast China has decreased, leading to a reduction in local supply of non-staple food. From the perspective of human resources, Yang et al. 8 found that the relationship between non-agricultural employment and food production presents an inverted U-shaped pattern, which means that in the case of a small supply of non-agricultural labor force, increasing non-agricultural employment will have a positive impact on food output, while in the case of a large supply of non-agricultural labor force, increasing non-agricultural employment is not conducive to food output increase. Abebaw et al. 9 investigated the impact of rural outmigration on food security in Ethiopia, and the results showed that rural outmigration significantly increased the daily calorie intake per adult by approximately 22%, reducing the gap and severity of food poverty by 7% and 4%, respectively.

Climate change. There is no consensus on the impact of climate change on food security. The majority of scholars assert that climate change will have significant negative effects on the availability, accessibility, and stability of food. Bijay et al. 10 argued that the ongoing global climate change has caused a range of issues, including increased carbon dioxide, frequent droughts, and temperature fluctuations, which pose significant obstacles to pest management, consequently impeding increased food production. Muhammad et al. 11 concluded through empirical analysis that climate change has a substantial adverse impact on irrigation water, agriculture, and rural livelihoods, and the latter three have a significant positive correlation with food security, suggesting that climate change is detrimental to food security. Atuoye et al. 12 examined the influence of gender, migration, and climate change on food security, and their findings revealed that as global climate changes, the impact of controlling carbon emissions on non-migrant food insecurity in Tanzania is reduced, while it exacerbates the impact on migrant food insecurity. However, some scholars contend that climate change can improve agricultural production conditions in certain regions, thereby facilitating increased food production and positively impacting food security 13 , 14 .

Government policy. Bizikova et al. 15 evaluated 73 intervention policies in a sample of 66 publications, of which 49 intervention policies had a positive impact on food security, 7 intervention policies had a negative impact, and 17 intervention policies had no impact. Chengyou et al. 16 used data such as mutual aid funds of impoverished villages in China to evaluate the effect of agricultural subsidies, and the empirical conclusion pointed out that agricultural subsidies can improve farmers' willingness to plant food, promote farmers in impoverished areas to increase the planting area, and help farmers improve their own food production capacity and economic income. Na et al. 17 proposed that food subsidies can increase the working time of part-time farmers in agricultural work, especially in food planting, and promote farmers to better switch between non-agricultural work and agricultural work. This subsidy effect is conducive to maintaining sufficient supply and sustainable development of food production.

The existing literature studies food security from different perspectives and draws reasonable conclusions and policy recommendations, but it fails to analyze the issue of food security under the comprehensive effect of resource endowment, climate change and government policy. Based on this, this paper proposes the Entropy Window two-stage DDF to measure the efficiency of hunger eradication, food security and improving nutrition in 29 provinces of China under the influence of exogenous variable rainfall. From the perspective of food security, the impact of resource endowment, climate change and government policy on food security is comprehensively considered. In addition, in terms of climate change, different from the existing research focusing on the negative effects of high temperature, low temperature and drought on food production, this paper focuses on the impact of extreme changes in rainfall on food security, providing a certain complement to the existing literature on food security research.

Research methods

The evolution of DEA methods has seen many discussions of the dynamic DEA model. Färe and Grosskopf 18 first established the concept of dynamic DEA, devised a form of dynamic analysis, and then proposed a delayed lag (carryover) variable for the dynamic model. Tone and Tsutsui 19 then extended it to a dynamic DEA approach based on weighted relaxation, including four types of connected activities: (1) desired (good); (2) undesired (bad); (3) discretionary (free); and (4) non-discretionary (fixed). Battese and Rao 20 and Battese et al. 21 next demonstrated that it is possible to compare the technical efficiencies of different groups using the Meta-frontier model. Portela and Thanassoulis 22 proposed a convex Meta-frontier concept that can take into account the technology of all groups, the state-of-the-art level of technological production, as well as the communication between groups and can be further extended to improve business performance. O’Donnell et al. 23 proposed a Meta-frontier model for defining technical efficiency using an output distance function, which accurately calculates group and Meta-frontier technical efficiencies and finds that the level of technology of all groups is superior to the level of technology of any one group.

In this paper, the evaluation performance based on DDF is better, which can provide more accurate estimation results. Therefore, this paper modifies the traditional DDF model, combines Dynamic DEA with Network Structure 19 , 24 and Entropy method 25 , and considers exogenous issues to construct Meta Entropy Two-Stage Dynamic DDF Under an Exogenous DEA Model in order to measure the efficiency of hunger eradication, food security, and improving nutrition in 29 provinces of China under the influence of rainfall.

The entropy method

In this model, the stage 2 (Hunger eradication and improving nutrition of sustainable stage) output item “Improving nutrition” covers four detailed indicators: (1) stunting rate of children under 5 years old; (2) malnutrition rate of children under 5 years old; (3) obesity rate of children under 5 years old; and (4) newborn visit rate. If these detailed indicators are put into DEA, then there will be problems that cannot be solved. Therefore, this model first uses the Entropy method and then finds the weights and output values of four detailed indicators of improving nutrition in stage 2. The Entropy method mainly includes the following four steps.

Step 1: Standardize the data of the four detailed indicators of improving nutrition in stage 2 in 29 provinces of China.

Here, \(r_{mn}\) is the standardized value of the \(n\) th indicator of the \(m\) th province; \(\mathop {\min }\limits_{m} x_{mn}\) is the minimum value of the \(n\) th indicator of the \(m\) th province; and \(\mathop {\max }\limits_{m} x_{mn}\) is the maximum value of the \(n\) th indicator of the \(m\) th province.

Step 2: Add up the standardized values of the four detailed indicators of improving nutrition in stage 2.

Here, \(P_{mn}\) represents the ratio of the standardized value of the \(n\) th indicator to the sum of the standardized values for the \(m\) th province.

Step 3: Calculate the entropy value ( \({\text{e}}_{{\text{n}}}\) ) for the \({\text{n}}\) th indicator.

Step 4: Calculate the weight of the \({\text{n}}\) th indicator \(\left( {{\text{w}}_{{\text{n}}} } \right)\) .

Using the above steps, we are able to find the weights and output values of the four detailed indicators of improving nutrition in stage 2.

Meta entropy two-stage dynamic DDF under an exogenous DEA model

Suppose there are two stages in each \(t \left( {t = 1, \ldots ,T} \right)\) time periods. In each time period, there are two stages, including agricultural production stage (stage 1), hunger eradication and improving nutrition of sustainable stage (stage 2).

In stage 1, there are \(b \left( {b = 1, \ldots ,B} \right)\) inputs \(x1_{bj}^{t}\) , producing \(a \left( {a = 1, \ldots , A} \right)\) desirable outputs \(y1_{aj}^{t}\) and \(o \left( {o = 1, \ldots , O} \right)\) undesirable outputs \(U1_{oj}^{t}\) . Stage 2 takes \(d \left( {d = 1, \ldots , D} \right)\) inputs \(x2_{dj}^{t}\) , creating \(s \left( {s = 1, \ldots ., S} \right)\) desirable outputs \(y2_{sj}^{t}\) and \(c \left( {c = 1, \ldots ., C} \right)\) undesirable outputs \(U2_{cj}^{t}\) ; the intermediate outputs connecting stages 1 and 2 are \(z_{hj}^{t} \left( {h = 1, \ldots ,H} \right)\) ; the carry-over variable is \(c_{lj}^{t} \left( {l = 1, \ldots ,L} \right)\) ; the exogenous variable is \(E_{vj}^{t} \left( {v = 1, \ldots ,V} \right)\) .

Figure  1 illustrates the framework diagram of the model. In stage 1, the input variables are agricultural employment, effective irrigation area and total agricultural water use, and the output variables are total agricultural output value and agricultural wastewater discharge. In stage 2, the input variable is local financial medical and health expenditure, and the output variables are the number of foodborne disease patients, the number of iodine deficiency disease patients, and improving nutrition. The link between stage 1 and stage 2 is the intermediate output: total agricultural output value. And the exogenous variable is rainfall.

Model framework.

Under the frontier, the DMU can choose the final output that is most favorable to its maximum value, so the efficiency of the decision unit under the common boundary can be solved by the following linear programming procedure.

Objective function

Efficiency of \({\text{DMUi}}\) is:

Here, \({\text{w}}_{1}^{{\text{t}}}\) and \({\text{w}}_{2}^{{\text{t}}}\) are the weights for stages 1 and stage 2, and \({ }\theta_{1}^{{\text{t}}}\) and \(\theta_{2}^{{\text{t}}}\) are the efficiency values for stages 1 and stage 2.

Subject to:

Stage 1: Agricultural production stage

Here, \({\text{q}}_{{{\text{bi}}1}}^{{\text{t}}}\) , \({\text{q}}_{{{\text{ai}}1}}^{{\text{t}}}\) , and \({\text{q}}_{{{\text{oi}}1}}^{{\text{t}}}\) denote the direction vectors associated with stage 1 inputs, desirable outputs, and undesirable outputs.

Stage 2: Hunger eradication and improving nutrition of sustainable stage

Here, \({\text{q}}_{{{\text{di}}2}}^{{\text{t}}}\) , \({\text{q}}_{{{\text{ci}}2}}^{{\text{t}}}\) , and \({\text{q}}_{{{\text{hi}}\left( {1,2} \right)}}^{{\text{t}}}\) denote the direction vectors associated with stage 2 inputs, undesirable outputs, and the intermediate outputs connecting stages 1 and 2.

The link of two periods

The exogenous variables

From the above results, the overall efficiency, the efficiency in each period, the efficiency in each stage, the efficiency in each stage in each period are obtained.

Input, desirable output, and undesirable output efficiencies

The disparity between the actual input–output indicators and the ideal input–output indicators under optimal efficiency represents the potential for efficiency improvement in terms of input and output orientation. This paper chooses the ratio of actual input–output values to the computed optimal input–output values as the efficiency measure for the input–output indicators. The relationship between the optimal value, actual value, and indicator efficiency is as follows:

If the actual input and undesirable output equals the optimal input and undesirable output, then the efficiencies of that input and undesirable output are equal to 1 and known as efficient. However, if the actual input exceeds the optimal input, then the efficiency of that input indicator is less than 1, which denotes being inefficient.

If the actual desirable output equals the optimal desirable output, then the efficiency of that desirable output is equal to 1 and is referred to as efficient. However, if the actual desirable output is less than the optimal desirable output, then the efficiency of that desirable output indicator is less than 1 and is considered inefficient. ME (Mean Efficiency) reflects the average efficiency of a certain region throughout the study period, with higher values indicating higher efficiency in that region.

Empirical study

Comparative analysis of total efficiency values considering and not considering exogenous variables.

As shown in Fig.  2 , in terms of the average total efficiency value for each region, without considering the exogenous variable rainfall, from 2016 to 2020, the average total efficiency values of the eastern, central, and western regions show a pattern of “eastern > western > central” in descending order. With the exogenous variable rainfall taken into account, the average total efficiency values for each region for each year were greater than the corresponding average total efficiency values without taking into account the exogenous variable rainfall, which may be attributed to the fact that rainfall plays a key role in irrigating the farmland and replenishing the soil moisture, which is an important factor in the process of agricultural production, and that the addition of rainfall has a more pronounced marginal effect on the increase in the total efficiency values. With the exogenous variable rainfall taken into account, the average total efficiency values for each region in each year are larger than the corresponding average total efficiency values without taking into account the exogenous variable rainfall, indicating that there is more room for improvement in the average total efficiency values without taking rainfall into account than in the efficiency values with rainfall taken into account. Except for 2016, when the average total efficiency value of the western region was greater than that of the eastern region and the central region, the average total efficiency values of the eastern, central, and western regions from 2017 to 2020 also showed a pattern of “eastern > western > central” from largest to smallest. It can be concluded that whether or not the exogenous variable rainfall is taken into account, the eastern region has a better overall efficiency in agricultural production and achieving food security than the western and central regions due to its better agricultural infrastructure, good economic base, and better educated labor force.

Average efficiency by region from 2016 to 2020.

The three regions of the East, Central and West maintain a similar fluctuating upward trend. The average efficiency in the eastern and western regions is relatively high, and the five-year fluctuation interval is small, ranging from 0.75 to 0.85. After considering the exogenous variable rainfall, the average total efficiency value in the central region increased from 0.62 to 0.66. However, compared with the eastern and western regions, the total efficiency in the central region is still at a lower level and the five-year fluctuation interval is larger, between 0.55 and 0.70, with the largest fluctuation interval in the average efficiency in 2017–2018, at − 0.11. This may be due to the downsizing of grain sowing area under the structural reform of the agricultural supply side, leading to a small decline in the total national grain output in 2018, which in turn affects the level of efficiency in eradicating hunger, guaranteeing food security and improving nutrition. From this, it can be concluded that the eastern and western regions should give full play to their original advantages and promote the modernization and sustainable development of agricultural production in order to accelerate the achievement of the three major goals of eradicating hunger, guaranteeing food security and improving nutrition, while the central region still has more room for improvement and needs to further play the role of agricultural policies to alleviate the people’s worries about food.

Table 1 Average efficiency by province and city from 2016 to 2020 demonstrates the average efficiency values for each province and city from 2016 to 2020 when rainfall is considered and not considered. From the point of view of the annual average total efficiency by province, after considering the exogenous variable rainfall, the efficiency value of most provinces has been improved. The average efficiency has also been improved from 0.6134 to 0.6189. Among them, the efficiency value of Qinghai increases from 0.8167 to 1, and the ranking also rises from 11th to 1st place. Qinghai is deep inland, with less rainfall throughout the year, and its agricultural and animal husbandry production is more sensitive to the changes of rainfall, and the addition of exogenous variable rainfall makes the average total efficiency more accurately portrayed, and achieves the DEA validity. Shandong’s ranking drops from 9 to 11th after considering the exogenous variable rainfall. As a major agricultural province, Shandong’s food production will be seriously affected by persistent heavy precipitation and other extreme weather events, which indicates that Shandong needs to take measures to strengthen the ability of its agricultural production to cope with extreme precipitation.

Two-stage average efficiency analysis

The average efficiency values of the two stages in both cases of considering exogenous variable rainfall and not considering exogenous variable rainfall are very similar, indicating that exogenous variable rainfall does not have much effect on the efficiency of stage 1 and stage 2, and therefore only the specific case with exogenous variable rainfall is discussed. Figures  3 and 4 show the efficiency values for Stage 1 and Stage 2 for each province and city for the years 2016–2020 when rainfall is considered. As shown in Fig.  3 , the difference between the efficiency values for Stage 1 and Stage 2 is still relatively significant. The efficiency of agricultural production in Stage 1 is significantly higher than that of hunger elimination, food security and nutritional improvement in Stage 2, and the fluctuation is relatively smooth, which indicates that there is still much room for improvement in China’s food production in terms of hunger elimination, food security and nutritional improvement, and that how to develop high-quality and high-efficiency agriculture and increase the output of food units is an urgent problem to be solved by each province.

Comparison of the average efficiency of the two phases by province from 2016 to 2020.

Average efficiency values for the two phases in each province from 2016 to 2020.

Specifically, there are large gaps in the efficiency of agricultural production in China's provinces, which can be roughly categorized into three types: the first type has an efficiency value of 1, realizing the DEA is effective, and is filled in green in Fig.  4 ; the second type has an efficiency value between 1 and the average, and is filled in yellow in Fig.  4 ; and the third type has an efficiency value below the average, and is filled in red in Fig.  4 .

In the first stage, the first category is Shanghai, Shandong, Tianjin, Beijing and other 15 provinces, whose agricultural production efficiency values are all 1, at the meta-frontier, and these provinces rely on a solid economic foundation and sound agricultural infrastructure to realize the optimal efficiency of effective inputs and outputs; the second category is Guangxi, Hubei, Sichuan, and Liaoning, whose agricultural production efficiencies are higher than the national average and close to the meta-frontier; the third category consists of 10 provinces such as Gansu, Inner Mongolia, Jilin, Heilongjiang, etc., whose economic development is relatively slow, meteorological conditions are poor, agricultural production is susceptible to meteorological disasters, and the efficiency of agricultural production is below the average level, among which the value of Gansu’s agricultural production efficiency is the lowest, 0.496.

In the second stage, the first category includes seven provinces, including Yunnan, Tianjin, Beijing, and Ningxia, which either have higher economic levels or better climatic conditions, and have the highest efficiency in eradicating hunger, achieving food security, and improving nutrition, with an efficiency value of 1; the second category includes eight provinces, including Shanghai, Chongqing, Jilin, and Shaanxi, which have an efficiency in eradicating hunger, achieving food security, and improving nutrition higher than the national average, and are close to the meta-frontier; the third category includes 14 provinces, including Fujian, Shanxi, Inner Mongolia, and Guangxi, which are below the national average, among which Sichuan has the lowest efficiency value of 0.1, which is evident that Sichuan, as a “Heavenly Grain Silo,” is more likely to speed up the realization of mechanization and digital development to improve comprehensive grain production capacity.

In summary, provinces with high efficiency values in agricultural production and in eradicating hunger, achieving food security and improving nutrition can be categorized into two groups, one of which is the developed and coastal provinces with good economic and climatic conditions, such as Beijing, Shanghai, Tianjin, and Hainan, can enhance agricultural sustainable efficiency and actively promote the sustainable development of the agricultural economy; the other category is the provinces with relatively backward economic development, including Yunnan, Ningxia, Qinghai, Heilongjiang and other central and western regions, although their development is relatively late and low, they have unique climatic conditions, geographic conditions, ecological conditions, and other resource advantages, which bring opportunities for sustainable agricultural development in the central and western regions. As for the provinces with lower efficiency values for agricultural production and hunger eradication, reaching food security and improving nutrition, they are not only affected by the level of economic development and ecological conditions such as climate and environment, but also by the level of urbanization, such as Fujian, Jiangsu, Zhejiang, Guangdong and other eastern coastal provinces with a high level of urbanization will also face pressure on the supply of agricultural products as the sown area of crops continues to decrease due to a combination of factors such as the occupation of arable land by construction sites as well as abandonment of land.

Comparative analysis of output indicator efficiency in the regions

Taking rainfall as an exogenous variable into account, the efficiency of the number of foodborne diseases patients and improving nutrition showed a higher pattern in the eastern and western regions and a lower pattern in the central region. Table 2 shows the efficiency values of each output indicator for 2016–2020. From 2016 to 2020, the efficiency of these two output indicators in the eastern and western regions showed an upward trend, while that in the central region showed a downward trend. It shows that the contribution of agricultural production to food security in the eastern and western regions is small, and more perfect institutional measures should be formulated to ensure food security; the contribution of agricultural production to improving nutrition in the central region is relatively small, and corresponding health expenditures need to be increased to improve people's own nutritional supplements. In terms of the efficiency of the number of iodine deficiency disease patients, the efficiency in the eastern and central regions was low and showed a downward trend from 2016 to 2020, while the efficiency in the western region was high and the fluctuation was relatively small. As people in the eastern and central regions can easily buy kelp, laver and other iodine-rich foods, local residents eat iodine-rich food at high frequency and in large amounts, while in the western region, which is far from the sea, daily eating may not meet the human body's daily demand for iodine. Therefore, in order to reduce the incidence of iodine deficiency diseases caused by geographical location and dietary habits, governments in the western region need to speed up the opening of transportation channels and purchasing channels for iodized salt and iodine-rich foods.

Conclusions and policy recommendation

The key to sustainable agricultural development lies in the organic integration of ecological sustainability, economic sustainability, and social sustainability, emphasizing the coordination between agroecological production capacity and human development. The conclusions of this paper are as follows.

First, in the total factor efficiency analysis, the average total efficiency values of the eastern, central, and western regions in each year when the exogenous variable rainfall is taken into account are higher than the corresponding average total efficiency values without considering exogenous variable rainfall. This may be due to the fact that rainfall is an important factor in the agricultural production process and the inclusion of rainfall has a more pronounced marginal effect on the increase in the total efficiency value. In addition, there is a certain difference between the average total efficiency values of the eastern and western regions regardless of whether exogenous variable rainfall is considered. Still, the difference is not very large, and all three regions maintain a similar trend of fluctuating upward. However, the average total efficiency value of the central region is still at a lower level compared to the eastern and western regions, and the fluctuations of the eastern and western regions over the 5 years are small, fluctuating between 0.75 and 0.85, while the average efficiency of the central region over the 5 years is low and fluctuates greatly, fluctuating between 0.55 and 0.70, and the fluctuations of the average efficiency in 2017–2018 are the largest, at − 0.11. Besides the average efficiency of the eastern region was slightly lower than that of the western region in 2016, the average efficiency of the three regions generally showed a decreasing hierarchy of eastern, western, and central regions one by one. In terms of the annual average total efficiency of each province, after considering the exogenous variable rainfall, the efficiency values of most provinces have improved, with Qinghai's average total efficiency rising to 1, achieving optimal input–output efficiency.. In contrast, Shandong's average efficiency ranking has declined.

Second, under the condition of considering the exogenous variable rainfall, the efficiency value in stage 1 (agricultural production stage) is significantly higher than the efficiency value in stage 2 (eliminating hunger, achieving food security and improving nutrition), and the fluctuation is relatively smooth, which suggests that China's food production still has a large room for improvement, and that the focus of attention should be different in different stages. Specifically, in stage 1, the provinces with lower agricultural production efficiency values belong to the central and western inland provinces with slower economic development and poorer meteorological conditions, while in stage 2, the provinces with lower efficiency value also include the more economically developed eastern coastal provinces, such as Fujian, Jiangsu, Zhejiang, Guangdong, etc. The rapid population growth in the developed eastern coastal areas, coupled with the impact of the construction of arable land and the impact of a combination of factors such as the abandonment of land, crop sowing area has been decreasing, resulting in per capita arable land area is lower than the national average level. This shows that although the developed eastern coastal provinces have a better foundation for agricultural development, they are also facing enormous pressure on the supply of agricultural products and increasingly fierce competition in the future industrial development.

Third, a comparative analysis of the efficiency of output indicators by region, taking into account the exogenous variable of rainfall, reveals that the efficiency of the number of foodborne diseases patients and improving nutrition are both high in the eastern and western regions and low in the central region and that the efficiency of these two output indicators shows a rising trend in the eastern and western regions and a declining trend in the central region in the period from 2016 to 2020. In terms of the efficiency of the number of iodine deficiency disease patients, the efficiency of the eastern and central regions is low and shows a similar downward trend over the five-year period, while the efficiency of the western region is high and fluctuates relatively little, with no significant trend of change.

Through the above empirical analysis, it can be seen that rainfall, an exogenous variable, has a significant impact on the average efficiency in the eastern, central and western regions. Therefore, this paper puts forward corresponding policy recommendations on hunger eradication, food security and improving nutrition. The specific recommendations are as follows:

First, continue to strengthen the construction of agricultural infrastructure and increase the per capita arable land. All regions, especially the central and western regions, need to continue to increase investment in agriculture, build agricultural infrastructure such as water conservancy facilities, transportation facilities, and electric power facilities, promote the transformation and upgrading of old agricultural infrastructure, help the rapid development of agricultural mechanization in China, further enhance the ability to resist natural disasters, and improve agricultural output and production efficiency. In this way, the contradiction between food production and the growing rigid demand for food can be alleviated.

Second, increase policy support for green agricultural production to ensure China's food security. Due to the developed industry and serious pollution, the eastern region should pay more attention to green agricultural production. Each province shall formulate corresponding subsidy plans for green agricultural production according to the specific conditions of the province, strengthen green technology to lead the green development of agriculture, increase the enthusiasm of farmers to carry out green agricultural production, promote the promotion of green agricultural production, decrease the use of harmful fertilizers, pesticides, agricultural film, etc., to reduce agricultural pollution, so as to increase the supply of green agricultural products on the market to decrease the prevalence of foodborne diseases.

Third, promote the diversification of agricultural production and enrich people's agricultural product consumption varieties. On the one hand, each province extends local agricultural production varieties according to climate conditions and resources, rationally layout the supply structure of agricultural products, and increase policies to encourage farmers to carry out diversified agricultural production. On the other hand, some regions are limited by resource endowments and cannot expand the types of agricultural production, so it is necessary to speed up the construction of infrastructure such as logistics and preservation, improve the system of connecting production and marketing of agricultural products, enrich the "vegetable basket" of people in these regions with poor agricultural resources, and meet people's diversified consumption demand for agricultural products. In addition, nutrition guidance, publicity and education should be strengthened to raise people's awareness of rational diet and nutritious diet.

Data availability

All data generated or analysed during this study are included in this published article [and its supplementary information files].

Upton, B. J., Cissé, D. J. & Barrett, B. C. Food security as resilience: Reconciling definition and measurement. Agr. Econ. 47 , 135–147 (2016).

Article   Google Scholar  

Jennifer, C. G. et al. Viewpoint: The case for a six-dimensional food security framework. Food Policy 106 , 102164 (2022).

Maxwell, S. Food security: A post-modern perspective. Food Policy 21 , 155–170 (1996).

Kang, S. et al. Improving agricultural water productivity to ensure food security in China under changing environment: From research to practice. Agr. Water Manage. 179 , 5–17 (2017).

Chloe, S., Jerry, K. & Tim, H. Managing irrigation under pressure: How supply chain demands and environmental objectives drive imbalance in agricultural resilience to water shortages. Agr. Water Manage. 243 , 106484 (2021).

Charoenratana, S. & Shinohara, C. Rural farmers in an unequal world: Land rights and food security for sustainable well-being. Land Use Policy 78 , 185–194 (2018).

Li, W. et al. Measuring urbanization-occupation and internal conversion of peri-urban cultivated land to determine changes in the peri-urban agriculture of the black soil region. Ecol. Indic. 102 , 328–337 (2019).

Yang, J., Wan, Q. & Bi, W. Off-farm employment and food production change: New evidence from China. China Econ. Rev. 63 , 135 (2020).

Abebaw, D. et al. Can rural outmigration improve household food security? Empirical evidence from Ethiopia. World Dev. 129 , 104879 (2020).

Bijay, S., Anju, P. & Samikshya, A. The impact of climate change on insect pest biology and ecology: Implications for pest management strategies, crop production, and food security. J. Agr. Food Res. 14 , 100733 (2023).

Google Scholar  

Muhammad, U. et al. Pathway analysis of food security by employing climate change, water, and agriculture nexus in Pakistan: Partial least square structural equation modeling. Environ. Sci. Pollut. Res. 30 , 88577–88597 (2023).

Atuoye, K. N. et al. Who are the losers? Gendered-migration, climate change, and the impact of large scale land acquisitions on food security in coastal Tanzania. Land Use Policy 101 , 105154 (2021).

Shi, W. et al. Has climate change driven spatio-temporal changes of cropland in northern China since the 1970s. Clim. Change 124 , 163–177 (2014).

Article   ADS   Google Scholar  

Tao, F. et al. Responses of wheat growth and yield to climate change in different climate zones of China, 1981–2009. Agr. Forest Meteorol. 189 , 91–104 (2014).

Bizikova, L., Jungcurt, S., McDougal, K. & Tyler, S. How can agricultural interventions enhance contribution to food security and SDG 2.1. Glob. Food Secur. 26 , 100450 (2020).

Chengyou, L. et al. The effectiveness assessment of agricultural subsidy policies on food security: Evidence from China’s poverty-stricken villages. Int. J. Environ. Res. Public Health 19 , 13797 (2022).

Na, X., Liqin, Z. & Xiyuan, L. Sustainable food production from a labor supply perspective: Policies and implications. Sustainability 14 , 15935 (2022).

Färe, R. & Grosskopf, S. Productivity and intermediate products: A frontier approach. Econ. Lett. 50 , 65–70 (1996).

Tone, K. & Tsutsui, M. Dynamic DEA: A slacks-based measure approach. Omega 38 , 145–156 (2010).

Battese, G. E. & Rao, D. P. Technology gap, efficiency, and a stochastic metafrontier function. Int. J. Bus. Econ. 1 (2), 87–93 (2002).

Battese, G. E., Rao, D. P. & O’Donnell, C. J. A metafrontier production function for estimation of technical efficiencies and technology gaps for firms operating under different technologies. J. Prod. Anal. 21 (1), 91–103 (2004).

Portela, M. C., Thanassoulis, E. & Simpson, G. Negative data in DEA: A directional distance approach applied to bank branches. J. Oper. Res. Soc. 55 (10), 1111–1121 (2004).

O’Donnell, C. J., Rao, D. P. & Battese, G. E. Metafrontier frameworks for the study of firm-level efficiencies and technology ratios. Empirical Econ. 34 (2), 231–255 (2008).

Charnes, A., Cooper, W. W. & Rhodes, E. Measuring the efficiency of decision making units. Eur. J. Oper. Res. 2 , 429–444 (1978).

Article   MathSciNet   Google Scholar  

Shannon, C. E. A mathematical theory of communication. Bell Syst. Tech. J. 27 , 379–423 (1948).

Download references

Author information

Authors and affiliations.

School of Marxism of Fujian Police College, Fuzhou, Fujian, 350007, People’s Republic of China

SiMan Jiang

School of Economics, Fujian Normal University, Fuzhou, Fujian, 350007, People’s Republic of China

Shuyue Chen, Qiqi Xiao & Zhong Fang

You can also search for this author in PubMed   Google Scholar

Contributions

SM.J and QQ.X. wrote the main manuscript text, SY.C and Z.F prepared figures and methodology. All authors reviewed the manuscript.

Corresponding author

Correspondence to Zhong Fang .

Ethics declarations

Competing interests.

The authors declare no competing interests.

Additional information

Publisher's note.

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Supplementary information 1., supplementary information 2., rights and permissions.

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/ .

Reprints and permissions

About this article

Cite this article.

Jiang, S., Chen, S., Xiao, Q. et al. Research on food security issues considering changes in rainfall. Sci Rep 14 , 19698 (2024). https://doi.org/10.1038/s41598-024-70803-x

Download citation

Received : 16 April 2024

Accepted : 21 August 2024

Published : 24 August 2024

DOI : https://doi.org/10.1038/s41598-024-70803-x

Share this article

Anyone you share the following link with will be able to read this content:

Sorry, a shareable link is not currently available for this article.

Provided by the Springer Nature SharedIt content-sharing initiative

• Food security
• Improving nutrition
• Entropy two-stage DDF

By submitting a comment you agree to abide by our Terms and Community Guidelines . If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.

Quick links

• Explore articles by subject
• Guide to authors
• Editorial policies

Sign up for the Nature Briefing: Anthropocene newsletter — what matters in anthropocene research, free to your inbox weekly.

Corporate Debt Maturity Matters for Monetary Policy

Joachim Jungherr

Matthias Meier

Timo Reinelt

Immo Schott

Download PDF (4 MB)

2024-30 | August 22, 2024

We provide novel empirical evidence that firms’ investment is more responsive to monetary policy when a higher fraction of their debt matures. In a heterogeneous firm New Keynesian model with financial frictions and endogenous debt maturity, two channels explain this finding: (1.) Firms with more maturing debt have larger roll-over needs and are therefore more exposed to fluctuations in the real interest rate (roll-over risk). (2.) These firms also have higher default risk and therefore react more strongly to changes in the real burden of outstanding nominal debt (debt overhang). Unconventional monetary policy, which operates through long-term interest rates, has larger effects on debt maturity but smaller effects on output and inflation than conventional monetary policy.

Suggested citation:

Jungherr, Joachim, Matthias Meier, Timo Reinelt, and Immo Schott. 2024. “Corporate Debt Maturity Matters for Monetary Policy.” Federal Reserve Bank of San Francisco Working Paper 2024-30. https://doi.org/10.24148/wp2024-30

Business Cycle and Early Warning Indicators for the Economy of Hong Kong– Challenges of Forecasting Work amid the COVID-19 Pandemic

• Research Paper
• Open access
• Published: 22 August 2024

Cite this article

You have full access to this open access article

• Sharon Pun-wai Ng 1 ,
• Eddie Ming-lok Kwok   ORCID: orcid.org/0000-0002-7096-2534 1 ,
• Brian Chi-yan Cheng 1 &
• Alfred Yiu-po Yuen 1  

16 Accesses

Explore all metrics

This paper discusses the use of the OECD’s framework to identify early warning indicators for building a Composite Leading Indicator (CLI) and the use of the Vector Autoregressive Model (VAR) for constructing a short-term forecasting model of economic growth of Hong Kong. With the onset of the COVID-19 pandemic, this paper further evaluates the performance of the CLI and forecasting model which were built based on pre-COVID-19 parameters and identify further adjustments to enhance the model performance.

Similar content being viewed by others

Does the Real Business Cycle Help Forecast the Financial Cycle?

Estimating a time-varying financial conditions index for South Africa

Forecasting the South African Financial Cycle: A Linear and Non-Linear Approach

Avoid common mistakes on your manuscript.

1 Introduction

With a view to predicting the turning points of the economy after the outbreak of the global financial crisis in 2008, the Census and Statistics Department (C&SD) of the Hong Kong Special Administrative Region conducted the study on the trial compilation of the CLI for Hong Kong and short-term GDP forecasting model (Ng et al., 2010 ), which was discussed at the International Seminar on Early Warning and Business Cycle Indicators organised by the United Nations Statistics Division in 2010. Subsequently, internal reviews had been carried out over time to assess the performance of the CLI and the forecasting model.

Starting from early 2020, the outbreak of COVID-19 pandemic has caused serious disruptions to economic activities. Producing early warning indicators and economic forecasts has become more challenging, as the empirical relationships between the pre-COVID-19 parameters and the business cycle have changed drastically in the midst of the pandemic. Thus, the performance of the early warning indicators has to be reviewed critically and the formation of the CLI and forecasting model has to be adjusted, which is the aim of this present research.

2 Literature Review

Numerous studies have been performed on the cyclical movement and short-term forecasting of output. The empirical research of business cycle started with Burns and Mitchell (1946), where a comprehensive list of composite leading, coincident, and lagging indicators of business cycles was developed using various economic variables for the economy of the United States. Since 1970s, OECD has developed a system for CLI and is now compiling it for 32 out of 38 member countries and 6 non-member countries. The methodology of OECD is based on various landmark work on business cycle, such as Bry and Boschan ( 1971 ), Prescott ( 1986 ), Kydland and Prescott ( 1990 ), Hodrick and Prescott ( 1997 ), Christiano-Fitzgerald ( 1999 ). The Conference Board has also compiled composite indices of leading, lagging, and coincident indicators for several major economies, using the methodology developed by the Department of Commerce of the United States (Green & Beckman, 1993 ). For short-term forecasting of output, various international studies (Rathjens & Robins, 1993 ; Miller & Chin, 1996 ; Stark, 2000 ; Koenig et al., 2001 ; Angelini et al., 2008 ) have focused on using monthly data to produce quarterly model forecasts, done by the traditional bridge equation, Principal Component Analysis (PCA), Vector Autoregressive (VAR) model (Sims, 1980 ) or other time series models.

For the case of Hong Kong, Gerlach and Yiu ( 2004 ) attempted a dynamic factor model for current quarter estimates of output. Genberg and Chang ( 2007 ) tried a VAR model to perform short-term forecasts of GDP using its components and other macroeconomic or monetary variables. On the other hand, Li and Kwok ( 2009 ) decomposed the output volatility of Hong Kong and analysed its cyclical movement across regions, making use of the Hodrick–Prescott (HP) filter and various Autoregressive Conditional Heteroskedasticity models. C&SD (Ng et al., 2010 ) conducted trial compilation of the Hong Kong’s CLI in 2010, utilising the high frequency indicators available and developed a high frequency macroeconomic model to forecast the one-quarter ahead real GDP, while an update (Ng et al., 2022 ) was done in 2021. Cheng et al. ( 2012 ) made similar attempt to compile the CLI of Hong Kong to track the local economy. On the other hand, it is noted that in recent years, the study of early warning indicators has been given more and more attention, and has found wider applications in the financial sector (Ito et al., 2014 ; Aldasoro et al., 2018 ). The emergence of COVID-19 pandemic has also made the forecasting of time series more challenging and Bobeica and Hartwig ( 2021 ) discussed the impact of the COVID-19 shock on time series models. Some recent studies have focused on the effect of COVID-19 on seasonal patterns in treatment of time series data, for example Eurostat ( 2020 ). Specifically, Lahiri and Yin ( 2024 ) compared additive outliers specified by Lucca and Wright ( 2021 ) method with those selected by automatic procedure.

3 Analytical Framework

3.1 cli of hong kong.

The method employed mostly follows the OECD’s framework (OECD, 2012 ), which is based on various landmark work on business cycle, such as Bry and Boschan ( 1971 ), Hodrick and Prescott ( 1997 ) and Christiano-Fitzgerald ( 1999 ). The first step involves the selection of seasonally adjusted real GDP as the reference indicator. Since real GDP figures are available quarterly but not monthly, the Denton method is applied to generate the monthly series of real GDP for compiling the CLI. Denton method is a temporal disaggregation commonly used to adjust quarterly series to monthly series (Denton, 1971 ). The quadratic match-sum method was adopted in this paper. For component indicators, they are selected based on the economic relevance and empirical relationship with GDP, as well as cross correlation analysis of the cyclical component obtained from the Hodrick-Prescott (HP) filter as shown in Eq. ( 1 ). The original time series and the trend component are indicated, respectively, by \(\:\{{\text{y}}_{t},\:t=\text{1,2},\dots\:,T\}\) and \({\tau _t},\;t = 1,2, \ldots,T\) . Thus, \({c_t} = {y_t} - {\tau _t}\) is the cyclical component of the series \(\:\:{\text{y}}_{t}\) . The HP filter gives a two-sided linear filter that computes the smoothed series \(\:{\tau\:}\) of \(\:\text{y}\) by minimising the variance of \(\:\text{y}\) around \(\:{\tau\:}\) , subject to a penalty that constrains the second difference of \(\:\:{\tau\:}.\)

The deviation of the original series from the trend (the first term of the equation) as well as the curvature of the estimated trend (the second term) is minimised through the HP filter. The trade-off between the two goals is controlled by the smoothing parameter λ. The higher the value of λ, the smoother is the estimated trend.

Cross correlation analysis shown in Eq. ( 2 ) is performed to verify if the component series has a leading or coincident relationship with the reference series. The indicator is said to be leading, or at least coincident, if the lag of the component series or the component series itself exhibits a strong correlation with the reference series.

where \(\:l=0,\pm\:1,\pm\:2,\dots\:\) , while

where T indicates the time period. Since the cyclical components are stationary, the mean and the variance of the cyclical series do not change over time. The cyclical component is obtained after subjecting the original series into the HP filter Footnote 1 , which is used to run the cross correlation between real GDP ( x ) and potential indicators ( y ) in Eq. ( 2 ) for selection of component indictors to form the CLI.

After selection of indicators, the data series of selected indicators are seasonally adjusted by the X-12-ARIMA method to remove possible seasonal patterns with additive outliers selected by automatic procedure, if the series are originally not seasonally adjusted and with high effect. A double HP filter Footnote 2 is adopted to extract cyclical component of the series, which are then normalised in order to express various indicators (which are originally measured in different units) in a common scale. Finally, the normalised series are aggregated by simple average (i.e. weighted equally) to form the CLI. The normalisation procedure introduces an implicit weighting of component series, however, with the series being weighted by the inverse of their mean absolute deviation (OECD, 2012 ).

3.2 Short-Term Forecast of Real GDP of Hong Kong

While CLI is useful in providing an early qualitative reference on turning points of the economy, it does not contain explanatory power in predicting output growth. For predicting short-run output growth, we attempt to make use of the potential high frequency component indicators identified in the above study of the CLI to form a monthly VAR model (Sims, 1980 ) (Eq.  6 ) so as to forecast the real GDP growth.

VAR (p) model for \(\:{\text{y}}_{\text{t}}\) in level form:

where y t refers to vector of variables; \(\mu \) refers to vector of constants; \(\:{\beta\:}_{k}\) refers to matrix of estimated coefficients; p refers to lag; and \(\:{\varepsilon\:}\) refers to vector of error terms.

Each equation in the unrestricted VAR can be estimated efficiently using ordinary least squares. Assuming the relationship identified in the past continues to hold in the forecast period, the period ahead forecasts can be obtained by simply substituting the estimated coefficients and their past values in the equations. In other words, the one-quarter ahead forecast can be produced when the data of component indicators in the last quarter are available. The model can also produce forecasts of current quarter real GDP growth with partial data of that quarter. To achieve this, recursive estimation, in which parameters are estimated recursively by constantly updating the information set when the forecast moves forward in time, is performed.

To further validate the VAR model, a benchmark simple autoregressive model with exogenous input, ARX(1,q-1) model, which does not depend on the lag of the indicator series, is built to compare the performance between the two models according to the root mean squared error (RMSE).

ARX(1,q-1) model for \(\:{y}_{t}\) in level form:

where \(\:{\text{y}}_{\text{t}}\) refers to the real GDP and \(\:{x}_{i,t}\) refers to the COVID-19 adjusted indicators series.

The out-of-sample approach is used to split the whole time series into two: (i) an initial fit period in which a model is trained; and (ii) a testing period of the latest three time points held out for estimating the loss of that model. In other words, the last three data of the time series are excluded in the model development and used for out-of-sample validation to assess the accuracy of the model. This validation methodology can preserve the temporal order of observations to deal with the dependency among observations and account for the potential temporal correlation between the consecutive values of the time series (Cerqueira et al. 2020 ).

4 Empirical Results

4.1 cli of real gdp of hong kong.

The CLI of Hong Kong for the pre-COVID-19 period included nine component indicators that are selected based on the cross correlation analysis, and taking into account the economic interpretation of potential indicators. These component indicators are retail sales, Hong Kong’s Purchasing Managers’ Index (PMI), Hang Seng Index, price of private domestic premises, transactions of the property market, M1, total merchandise trade, visitor arrivals and manufacturing PMI of the mainland of China. Details of the cross correlation analysis are given in Table  1 . A list of other potential indicators which was tested but found to have weak or no leading correlation with GDP is given in Table  4 in Appendix.

The COVID-19 pandemic, however, had generated some behavioural changes in Hong Kong’s economy with more reliance on local economic activities. In particular, the contribution to GDP of tourism reduced from around 5% from the pre-COVID-19 period to only 0.4% in 2020 and as low as 0.1% in 2021. As the inbound tourism was virtually halted after the outbreak of the pandemic, the number of visitor arrivals was no longer relevant in predicting GDP movement. Thus, it is taken out of the component indicator set. On the other hand, having further reviewed the indicator set and conducted several trials using relevant indicators that can potentially better correlate with the ups and downs in local economic activities along with changing COVID-19 situation, the number of passenger journey by taxis is chosen to replace the visitor arrivals after the fourth quarter of 2019, where its peak and trough are slightly ahead of the movement of real GDP (Fig. 1 ).

Normalised selected component series vs. normalised reference series

Statistically, the cyclical properties of CLI and real GDP can also be confirmed by cross correlation analysis. As shown in Table  2 , while strong leading correlation between CLI and real GDP is identified in both the CLI of pre-COVID-19 indicators and CLI of COVID-19 adjusted indicators, the latter has a lower standard deviation. Furthermore, in Fig. 2 the CLI using the COVID-19 adjusted indicators exhibits a much better leading property of real GDP than the CLI using pre-COVID-19 indicators since the latter is heavily dragged by the halt in inbound tourism after Q4 2019.

Cycle of CLI and real GDP, Jan 2018– March 2022

4.2 VAR Model for Short-Term Real GDP Forecasting

Regarding the forecasting of quarter-to-quarter percentage change of seasonally adjusted real GDP growth using the unrestricted VAR, all indicator series used in constructing CLI have been applied in the model. The results of current quarter model forecasts with two months’ actual data are shown in Fig. 3 . Given the higher volatility of GDP growth since Q1 2020, the forecast results are generally not as accurate as before. To examine the COVID-19 shock through the VAR models, impulse response functions (IRFs) analysed by with and without the epidemic period are generated (Bobeica & Hartwig, 2021 ), as presented in Figs. 4 and 5 respectively. The variable of visitor arrivals is clearly subject to structural break as its IRFs with the epidemic period far exceed the 95% confidence interval of the VAR model without that period (Fig. 4 j), while the IRFs of the variable of passenger journey by taxis are relatively more stable (Fig. 5 j). Replacing the variable of visitor arrivals by the taxi passenger journey to build the COVID-19 adjusted indicators, better forecast results are observed compared to the forecast results using the pre-COVID-19 indicators, with the absolute difference between the actual and predicted growth lowered by 2.3% points on average during the 9 quarters from Q1 2020 to Q1 2022.

Quarterly growth (%) of real GDP, Q1 2018– Q1 2022

IRFs of VAR using pre-COVID-19 indicators

IRFs of using COVID-19 adjusted indicators

Regarding the standard forecasting accuracy metrics, the VAR model performs better than the benchmark AR model, given its much lower values of RMSE for both the training and validation sets (Table 3 ).

5 Concluding Remark

The present study reviews the set of relevant statistical indicators that can exhibit good leading properties for real GDP in Hong Kong amid COVID-19 pandemic to form the CLI and the forecasting model using VAR. Necessary adjustment is made to take into account the “new normal” brought by the pandemic.

Looking ahead, there is no guarantee that the model will continue to be valid along with the evolving situation of the COVID-19 pandemic. Indeed, it should be noted that the forecasting ability of CLI and VAR model is based on past relationship between the occurrence of crises and relevant data. The historical patterns may not hold in future if an economy is impacted significantly by external shocks which have no historic references. This is particularly true for economies like Hong Kong that are small and open and can be heavily influenced by external shocks. The composition of the component indicators may need to be reviewed periodically to ensure their relevance. On the other hand, with the rapidly changing economic condition in the foreseeable future, further development and enhancement in the study of early warning indicators can be expected.

See Table  4 .

Quarterly data is applied to run the cross correlation and the default smoothing parameter value is applied (i.e. λ = 1 600) for the HP filter.

Based on the OECD framework, the double HP filter is applied to extract the cycle after de-trending and smoothing. In the first step, the long term trend is removed by applying the HP filter with a high λ value (i.e. λ = 133 107.94). Then, the HP filter with a small value of λ (i.e. λ = 13.93) is applied again to smooth the cycle component.

Aldasoro, I., Borio, C., & Drehmann, M. (2018). Early warning indicators of banking crises: Expanding the family. Bank for International Settlements Quarterly Review (March) , 29–45.

Angelini, E., Camba-Mendez, G., Giannone, D., Reichlin, L., & Rustler, G. (2008). Short-term forecasts of euro area GDP growth. European Central Bank Working Paper 949 , 1–31.

Bobeica, E., & Hartwig, B. (2021). The COVID-19 shock and challenges for time series models. European Central Bank Working Paper 2558 , 1–41.

Bry, G., & Boschan, C. (1971). Cyclical analysis of time series: Selected procedures and compute programs. NBER .

Cerqueira, V., Torgo, L., & Mozetič, I. (2020). Evaluating time series forecasting models: An empirical study on performance estimation methods. Machine Learning , 109 , 1997–2028.

Article   Google Scholar  

Cheng, M., Chung, L., Tam, C. S., Yuen, R., Chan, S., & Yu, I. W. (2012). Tracking the Hong Kong economy. Hong Kong Monetary Authority Occasional Paper , 03/2012 , 1–29.

Google Scholar  

Christiano, L. J., & Fitzgerald, T. J. (1999). The band pass filter. NBER Working Paper 7257 , 1–73.

Denton, F. T. (1971). Adjustment of monthly or quarterly series to annual totals: An approach based on quadratic minimization. Journal of the American Statistical Association , 66(333) , 99–102.

Eurostat (2020). Methodological note - guidance on time series treatment in the context of the COVID-19 crisis.

Genberg, H., & Chang, J. (2007). A VAR framework for forecasting Hong Kong’s output and inflation. Hong Kong Institute for Monetary Research Working Paper 2 , 1–25.

Gerlach, S., & Yiu, M. (2004). A dynamic factor model for current-quarter estimates of economic activity in Hong Kong. Hong Kong Institute for Monetary Research Working Paper 16 , 1–24.

Green, G. R., & Beckman, B. A. (1993). Business cycle indicators: Upcoming revision of the composite indexes. The Department of Commerce of the United States, Survey of Current Business, Oct 1993 , 44–51.

Hodrick, R. J., & Prescott, E. C. (1997). Postwar U.S. business cycles: An empirical investigation. Journal of Money Credit and Banking , 29 , 1–18.

Ito, Y., Kitamura, T., Nakamura, K., & Nakazawa, T. (2014). New financial activity indexes: Early warning system for financial imbalances in Japan. Bank of Japan Working Paper , 14 , 1–56.

Koenig, E., Dolmas, S., & Piger, J. M. (2001). The use and abuse of real-time data in economic forecasting. Federal Reserve Bank of St. Louis Working Paper 15 , 1–49.

Kydland, F. E., & Prescott, E. C. (1990). Business cycles: Real facts and a monetary myth. Federal Reserve Bank of Minneapolis Quarterly Review , 3–18.

Lahiri, K., & Yin, Y. (2024). Seasonality in U.S. disability applications, labor market and the pandemic echoes. Labour Economics , 87 .

Li, K. W., & Kwok, M. L. (2009). Output volatility of five crisis-affected East Asia economies. Japan and the World Economy , 21 , 172–182.

Lucca, D., & Wright, J. H. (2021). Reasonable seasonals? Seasonal echoes in economic data after COVID-19. Federal Reserve Bank of New York Liberty Street Economics .

Miller, P. J., & Chin, D. (1996). Using monthly data to improve quarterly model forecasts. Federal Reserve Bank of Minneapolis Quarterly Review Spring , 16–33.

Ng, P. W., Kwok, M. L., & Tam, K. M. (2010). Research study of business cycle and early warning indicators for the economy of Hong Kong. Contributed paper for the Third International Seminar on Early Warning and Business Cycle Indicators, November 2010 , 1–19.

Ng, P. W., Kwok, M. L., & Cheng, C. Y. (2022). Research study of business cycle and early warning indicators for the economy of Hong Kong– Challenges of forecasting work amid the COVID-19 pandemic. Proceedings 63rd ISI World Statistics Congress, 11–16 July 2021 , 1–4.

OECD (2012). OECD System of Composite Leading Indicators, April 2012. 1–18.

Prescott, E. C. (1986). Theory ahead of business cycle measurement. Federal Reserve Bank of Minneapolis Quarterly Review Fall, 9–22.

Rathjens, P., & Robins, R. (1993). Forecasting quarterly data using monthly information. Journal of Forecasting , 12 , 321–330.

Sims, C. (1980). Macroeconomics and reality. Econometrica , 48 , 1–48.

Stark, T. (2000). Does current-quarter information improve quarterly forecasts for the U.S. economy? Federal Reserve Bank of Philadelphia Working Paper 2 , 1–60 .

Download references

Author information

Authors and affiliations.

Census and Statistics Department, The Government of the Hong Kong Special Administrative Region of the People’s Republic of China, Hong Kong, Hong Kong SAR, China

Sharon Pun-wai Ng, Eddie Ming-lok Kwok, Brian Chi-yan Cheng & Alfred Yiu-po Yuen

You can also search for this author in PubMed   Google Scholar

Corresponding author

Correspondence to Eddie Ming-lok Kwok .

Ethics declarations

The authors have no financial or proprietary interests in any material discussed in this article.

Additional information

Publisher’s note.

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/ .

Reprints and permissions

About this article

Ng, S.Pw., Kwok, E.Ml., Cheng, B.Cy. et al. Business Cycle and Early Warning Indicators for the Economy of Hong Kong– Challenges of Forecasting Work amid the COVID-19 Pandemic. J Bus Cycle Res (2024). https://doi.org/10.1007/s41549-024-00098-4

Download citation

Received : 07 September 2022

Accepted : 02 July 2024

Published : 22 August 2024

DOI : https://doi.org/10.1007/s41549-024-00098-4

Share this article

Anyone you share the following link with will be able to read this content:

Sorry, a shareable link is not currently available for this article.

Provided by the Springer Nature SharedIt content-sharing initiative

• Leading indicators
• Forecasting
• Short-term forecast

JEL Classification

• Find a journal
• Publish with us
• Track your research

Suggestions or feedback?

MIT News | Massachusetts Institute of Technology

• Machine learning
• Sustainability
• Black holes
• Classes and programs

Departments

• Aeronautics and Astronautics
• Brain and Cognitive Sciences
• Architecture
• Political Science
• Mechanical Engineering

Centers, Labs, & Programs

• Abdul Latif Jameel Poverty Action Lab (J-PAL)
• Picower Institute for Learning and Memory
• Lincoln Laboratory
• School of Architecture + Planning
• School of Engineering
• School of Humanities, Arts, and Social Sciences
• Sloan School of Management
• School of Science
• MIT Schwarzman College of Computing

MIT engineers’ new theory could improve the design and operation of wind farms

Press contact :, media download.

*Terms of Use:

Images for download on the MIT News office website are made available to non-commercial entities, press and the general public under a Creative Commons Attribution Non-Commercial No Derivatives license . You may not alter the images provided, other than to crop them to size. A credit line must be used when reproducing images; if one is not provided below, credit the images to "MIT."

Previous image Next image

The blades of propellers and wind turbines are designed based on aerodynamics principles that were first described mathematically more than a century ago. But engineers have long realized that these formulas don’t work in every situation. To compensate, they have added ad hoc “correction factors” based on empirical observations.

Now, for the first time, engineers at MIT have developed a comprehensive, physics-based model that accurately represents the airflow around rotors even under extreme conditions, such as when the blades are operating at high forces and speeds, or are angled in certain directions. The model could improve the way rotors themselves are designed, but also the way wind farms are laid out and operated. The new findings are described today in the journal Nature Communications , in an open-access paper by MIT postdoc Jaime Liew, doctoral student Kirby Heck, and Michael Howland, the Esther and Harold E. Edgerton Assistant Professor of Civil and Environmental Engineering.

“We’ve developed a new theory for the aerodynamics of rotors,” Howland says. This theory can be used to determine the forces, flow velocities, and power of a rotor, whether that rotor is extracting energy from the airflow, as in a wind turbine, or applying energy to the flow, as in a ship or airplane propeller. “The theory works in both directions,” he says.

Because the new understanding is a fundamental mathematical model, some of its implications could potentially be applied right away. For example, operators of wind farms must constantly adjust a variety of parameters, including the orientation of each turbine as well as its rotation speed and the angle of its blades, in order to maximize power output while maintaining safety margins. The new model can provide a simple, speedy way of optimizing those factors in real time.

“This is what we’re so excited about, is that it has immediate and direct potential for impact across the value chain of wind power,” Howland says.

Modeling the momentum

Known as momentum theory, the previous model of how rotors interact with their fluid environment — air, water, or otherwise — was initially developed late in the 19th century. With this theory, engineers can start with a given rotor design and configuration, and determine the maximum amount of power that can be derived from that rotor — or, conversely, if it’s a propeller, how much power is needed to generate a given amount of propulsive force.

Momentum theory equations “are the first thing you would read about in a wind energy textbook, and are the first thing that I talk about in my classes when I teach about wind power,” Howland says. From that theory, physicist Albert Betz calculated in 1920 the maximum amount of energy that could theoretically be extracted from wind. Known as the Betz limit, this amount is 59.3 percent of the kinetic energy of the incoming wind.

But just a few years later, others found that the momentum theory broke down “in a pretty dramatic way” at higher forces that correspond to faster blade rotation speeds or different blade angles, Howland says. It fails to predict not only the amount, but even the direction of changes in thrust force at higher rotation speeds or different blade angles: Whereas the theory said the force should start going down above a certain rotation speed or blade angle, experiments show the opposite — that the force continues to increase. “So, it’s not just quantitatively wrong, it’s qualitatively wrong,” Howland says.

The theory also breaks down when there is any misalignment between the rotor and the airflow, which Howland says is “ubiquitous” on wind farms, where turbines are constantly adjusting to changes in wind directions. In fact, in an  earlier paper in 2022, Howland and his team found that deliberately misaligning some turbines slightly relative to the incoming airflow within a wind farm significantly improves the overall power output of the wind farm by reducing wake disturbances to the downstream turbines.

In the past, when designing the profile of rotor blades, the layout of wind turbines in a farm, or the day-to-day operation of wind turbines, engineers have relied on ad hoc adjustments added to the original mathematical formulas, based on some wind tunnel tests and experience with operating wind farms, but with no theoretical underpinnings.

Instead, to arrive at the new model, the team analyzed the interaction of airflow and turbines using detailed computational modeling of the aerodynamics. They found that, for example, the original model had assumed that a drop in air pressure immediately behind the rotor would rapidly return to normal ambient pressure just a short way downstream. But it turns out, Howland says, that as the thrust force keeps increasing, “that assumption is increasingly inaccurate.”

And the inaccuracy occurs very close to the point of the Betz limit that theoretically predicts the maximum performance of a turbine — and therefore is just the desired operating regime for the turbines. “So, we have Betz’s prediction of where we should operate turbines, and within 10 percent of that operational set point that we think maximizes power, the theory completely deteriorates and doesn’t work,” Howland says.

Through their modeling, the researchers also found a way to compensate for the original formula’s reliance on a one-dimensional modeling that assumed the rotor was always precisely aligned with the airflow. To do so, they used fundamental equations that were developed to predict the lift of three-dimensional wings for aerospace applications.

The researchers derived their new model, which they call a unified momentum model, based on theoretical analysis, and then validated it using computational fluid dynamics modeling. In followup work not yet published, they are doing further validation using wind tunnel and field tests.

Fundamental understanding

One interesting outcome of the new formula is that it changes the calculation of the Betz limit, showing that it’s possible to extract a bit more power than the original formula predicted. Although it’s not a significant change — on the order of a few percent — “it’s interesting that now we have a new theory, and the Betz limit that’s been the rule of thumb for a hundred years is actually modified because of the new theory,” Howland says. “And that’s immediately useful.” The new model shows how to maximize power from turbines that are misaligned with the airflow, which the Betz limit cannot account for.

The aspects related to controlling both individual turbines and arrays of turbines can be implemented without requiring any modifications to existing hardware in place within wind farms. In fact, this has already happened, based on earlier work from Howland and his collaborators two years ago that dealt with the wake interactions between turbines in a wind farm, and was based on the existing, empirically based formulas.

“This breakthrough is a natural extension of our previous work on optimizing utility-scale wind farms,” he says, because in doing that analysis, they saw the shortcomings of the existing methods for analyzing the forces at work and predicting power produced by wind turbines. “Existing modeling using empiricism just wasn’t getting the job done,” he says.

In a wind farm, individual turbines will sap some of the energy available to neighboring turbines, because of wake effects. Accurate wake modeling is important both for designing the layout of turbines in a wind farm, and also for the operation of that farm, determining moment to moment how to set the angles and speeds of each turbine in the array.

Until now, Howland says, even the operators of wind farms, the manufacturers, and the designers of the turbine blades had no way to predict how much the power output of a turbine would be affected by a given change such as its angle to the wind without using empirical corrections. “That’s because there was no theory for it. So, that’s what we worked on here. Our theory can directly tell you, without any empirical corrections, for the first time, how you should actually operate a wind turbine to maximize its power,” he says.

Because the fluid flow regimes are similar, the model also applies to propellers, whether for aircraft or ships, and also for hydrokinetic turbines such as tidal or river turbines. Although they didn’t focus on that aspect in this research, “it’s in the theoretical modeling naturally,” he says.

The new theory exists in the form of a set of mathematical formulas that a user could incorporate in their own software, or as an open-source software package that can be freely downloaded from GitHub . “It’s an engineering model developed for fast-running tools for rapid prototyping and control and optimization,” Howland says. “The goal of our modeling is to position the field of wind energy research to move more aggressively in the development of the wind capacity and reliability necessary to respond to climate change.”

The work was supported by the National Science Foundation and Siemens Gamesa Renewable Energy.

Share this news article on:

Related links.

• Unified momentum model software
• Howland Lab
• Department of Civil and Environmental Engineering

Related Topics

• Civil and environmental engineering
• Renewable energy
• Climate change

Related Articles

Michael Howland gives wind energy a lift

A new method boosts wind farms’ energy output, without new equipment

The answer may be blowing in the wind

Previous item Next item

More MIT News

Pursuing the secrets of a stealthy parasite

Read full story →

Study of disordered rock salts leads to battery breakthrough

Toward a code-breaking quantum computer

Uphill battles: Across the country in 75 days

3 Questions: From the bench to the battlefield

Duane Boning named vice provost for international activities

• More news on MIT News homepage →

Massachusetts Institute of Technology 77 Massachusetts Avenue, Cambridge, MA, USA

• Map (opens in new window)
• Events (opens in new window)
• People (opens in new window)
• Careers (opens in new window)
• Accessibility
• Social Media Hub
• MIT on Facebook
• MIT on YouTube
• MIT on Instagram

Information

• Author Services

Initiatives

You are accessing a machine-readable page. In order to be human-readable, please install an RSS reader.

All articles published by MDPI are made immediately available worldwide under an open access license. No special permission is required to reuse all or part of the article published by MDPI, including figures and tables. For articles published under an open access Creative Common CC BY license, any part of the article may be reused without permission provided that the original article is clearly cited. For more information, please refer to https://www.mdpi.com/openaccess .

Feature papers represent the most advanced research with significant potential for high impact in the field. A Feature Paper should be a substantial original Article that involves several techniques or approaches, provides an outlook for future research directions and describes possible research applications.

Feature papers are submitted upon individual invitation or recommendation by the scientific editors and must receive positive feedback from the reviewers.

Editor’s Choice articles are based on recommendations by the scientific editors of MDPI journals from around the world. Editors select a small number of articles recently published in the journal that they believe will be particularly interesting to readers, or important in the respective research area. The aim is to provide a snapshot of some of the most exciting work published in the various research areas of the journal.

Original Submission Date Received: .

• Active Journals
• Find a Journal
• Proceedings Series
• For Authors
• For Reviewers
• For Editors
• For Librarians
• For Publishers
• For Societies
• For Conference Organizers
• Open Access Policy
• Institutional Open Access Program
• Special Issues Guidelines
• Editorial Process
• Research and Publication Ethics
• Article Processing Charges
• Testimonials
• Preprints.org
• SciProfiles
• Encyclopedia

Article Menu

• Subscribe SciFeed
• Recommended Articles
• Google Scholar
• on Google Scholar
• Table of Contents

Find support for a specific problem in the support section of our website.

Please let us know what you think of our products and services.

Visit our dedicated information section to learn more about MDPI.

JSmol Viewer

Parametric optimization study of novel winglets for transonic aircraft wings.

Share and Cite

Padmanathan, P.; Aswin, S.; Satheesh, A.; Kanna, P.R.; Palani, K.; Devi, N.R.; Sobota, T.; Taler, D.; Taler, J.; Węglowski, B. Parametric Optimization Study of Novel Winglets for Transonic Aircraft Wings. Appl. Sci. 2024 , 14 , 7483. https://doi.org/10.3390/app14177483

Padmanathan P, Aswin S, Satheesh A, Kanna PR, Palani K, Devi NR, Sobota T, Taler D, Taler J, Węglowski B. Parametric Optimization Study of Novel Winglets for Transonic Aircraft Wings. Applied Sciences . 2024; 14(17):7483. https://doi.org/10.3390/app14177483

Padmanathan, Panneerselvam, Seenu Aswin, Anbalagan Satheesh, Parthasarathy Rajesh Kanna, Kuppusamy Palani, Neelamegam Rajan Devi, Tomasz Sobota, Dawid Taler, Jan Taler, and Bohdan Węglowski. 2024. "Parametric Optimization Study of Novel Winglets for Transonic Aircraft Wings" Applied Sciences 14, no. 17: 7483. https://doi.org/10.3390/app14177483

Article Metrics

Article access statistics, further information, mdpi initiatives, follow mdpi.

Subscribe to receive issue release notifications and newsletters from MDPI journals