## 2 Data and Methods
### Var  
VAR models are mainly used to predict and analyze the dynamic impact of random disturbances to the system, including the magnitude, positive and negative directions, and duration of the impact [19; 20]. We used the Granger causality test and impulse response analysis methods that are suitable for exploring the relationship of geographical variables to assess the temporal dynamics and cumulative effects of drought-ecosystem function interactions in the Songnen Plain from 2000 to 2020 based on SPEI, LAI, ET, and GPP variables. In a VAR model, each variable is a linear function of its own lagged value and the lagged values of the other variables. The general mathematical expression of the model is:  
\[y_{t}=\mu+A_{1}y_{t-1}+\cdots+A_{p}y_{t-p}+\mu_{t} \tag{3}\]  
where \(y_{t}\) denotes a vector of \(K\) variables (i.e., variables in the system that are explained by other variables, often called \"dependent variables\") assigned a spherical perturbation term of the same dimension at time \(t\) and \(u_{t}\) (i.e., a white noise process),\(\mu\) is a \(K\times 1\) column vector of constant terms, \(A_{1}\),..., \(A_{p}\) is a matrix of regression coefficients of dimension \(K\times K\), and p is the lag length [19; 42]. The \(p\) parameter determines how much prior information on all variables (i.e., lag time) is included in the model to predict future behavior [43]. Herein, Equation (3) is written in matrix form taking SPEI and LAI as an example:  
\[ \begin{bmatrix}SPEI_{t}\\ LAI_{t}\end{bmatrix} =\sum_{i=1}^{p} \begin{bmatrix}a^{\prime}_{11,i}&a^{\prime }_{12,i}\\ a^{\prime}_{21,i}&a^{\prime}_{22,i}\end{bmatrix}\begin{bmatrix}SPEI_{t-i}\\ LAI_{t-i}\end{bmatrix} +CD_{t}+ \begin{bmatrix}u_{1t}\\ u_{2t}\end{bmatrix}  \tag{4}\]  
Granger causality test is a statistical hypothesis testing method [44] that can be applied to time series data analysis. Assuming there exists variables \(X\) and \(Y\), if the prediction effect of variable \(Y\) is better than that of \(Y\) only predicted by the past information of \(Y\), namely, variable \(X\) helps to explain the future changes of variable \(Y\), then variable \(X\) is considered to be the Granger cause of variable \(Y\), and vice versa. Essentially, two regression models containing \(X\), \(Y\), and \(Y\) only are established, respectively. Then, the residual sum of squares of these two models is used to construct the F-statistic to determine the relationship with the original hypothesis, i.e., the existence of causality. Here, we use this method to assess whether there is a potential causal relationship between SPEI, LAI, ET, and GPP. Note that the Granger causality test requires the input time series to be stable, whereas physical geography data generally has seasonality. Thus, it is necessary to remove seasonality from the data. Data stability is judged according to the ADF (Augmented Dickey-Fuller) coefficients. The results are scientific only when all variables involved in the modeling are cointegrated with the same order. After making a judgment about the existence of causality between variables, impulse response analysis can help to understand the direction and strength of the relationships between the variables more deeply. The impulse response function is a causal analysis method to study the relationship between variables within a dynamic system, which has been popularized along with the development of the VAR model. The function is expressed based on the Wold moving average representation of the VAR(\(p\)) process, measuring the response process of current and future values of one variable over time after being subjected to a positive impulse of another variable [45], which is described as follows:  
\[y_{t}=\Phi_{0}u_{t}+\Phi_{1}u_{t-1}+\Phi_{2}u_{t-2}+\cdots \tag{5}\]  
Additionally, the cumulative impulse response function is used to examine the cumulative effects that a change in a variable in the system at a given moment will have on the system internally over a future period.