1 Introduction

Vegetation-climate interaction, defined as the significant mutual influence between vegetation and climate from local to global scales (Foley et al. 1998; Strengers et al. 2010), is important to current climate change and future climate projections(Quillet et al. 2010; Wu et al. 2016; Cui et al. 2020). As a crucial component of the earth system, terrestrial vegetation plays an active role in regional and global climate through its direct and indirect influences on biophysical and biochemical feedback mechanisms, such as land surface albedo, evapotranspiration, roughness length, absorption and emissions of gases, and so on (Betts et al. 1997; Still et al. 2003; Brovkin et al. 2006; Bonan 2008; Zeng et al. 2017; Duveiller et al. 2018; Forzieri et al. 2020; Li et al. 2020; Chen et al. 2021; Wei et al. 2023). In turn, climate change has a substantial impact on the phenology and activity of terrestrial vegetation through variations in local meteorological variables (Wu et al. 2015; Afuye et al. 2021). In general, temperature is identified as the critical factor of vegetation growth in Arctic and high latitude of Eurasia (Tanja et al. 2003; Tape et al. 2006; Kong et al. 2017), while vegetation variation is strongly affected by precipitation anomalies in midlatitude arid and semiarid regions of Northern Hemisphere with time-lag effects via soil moisture (Djebou et al. 2015; Wu et al. 2015; Chen et al. 2020), and incoming solar radiation dominates vegetation growth in low latitude humid region such as Amazon rain forests and eastern and southern Asia (Nemani et al. 2003; Zhao et al. 2018). Hence, for better predicting vegetation responses and feedbacks to climate change, to comprehend the linkage between climate variability and vegetation activities is of great importance.

Featured by a northwest-southeast elevation gradient, the low-latitude highlands of China (21–29° N, 97–108° E, Cao et al. 2014, the rectangular black box in Fig. 1, hereafter abbreviated as CLLH) is situated in the subtropical southwest China with a mean altitude of 1500 m above sea level. Being adjacent to the tropical Pacific and tropical Indian Ocean, winter CLLH climate conditions are substantially affected by El Niño-Southern Oscillation (ENSO), the strongest interannual air-sea coupled mode in the tropics. For winter precipitation (Precip), generally speaking, reinforced Precip is usually observed over the CLLH when winter El Niño occurs via enhanced water vapor transport and/or upward vertical motion (Wang et al. 2000; Feng et al. 2010, 2014; Yuan et al. 2014), Gao et al. (2020) also suggested that winter El Niño (La Niña) events intensify (weaken) the Precip extremes over the CLLH. However, the decreased Precip anomaly over the CLLH is weak during the La Niña mature winters (Lin et al. 2023), indicating an asymmetric Precip response to the ENSO. Tan et al. (2017) and Gan et al. (2019) pointed out that the distinct Precip deficit over the entire CLLH during the winter 2009/2010 is caused by the strong central Pacific El Niño event based on the results of numerical experiments; on the contrary, Li et al. (2019b) showed that during the strong winter 1998/1999 central Pacific La Niña event the Precip anomalies over the CLLH are relatively small and exhibit a west-to-east dipole pattern. Besides, the decadal modulation of ENSO on winter Precip over the western part of the CLLH is also evident, which is enhanced after the mid-1990s (Zhang et al. 2021). For winter surface air temperature (SAT), the East Asian winter monsoon tends to be weak (strong) in the winter of El Niño (La Niña) (Ha et al. 2012; Ma and Chen 2021), which is crucial to the variation of CLLH SAT. Chen and Song (2018) elucidated that the La Niña-like abnormal sea surface temperature (SST) pattern over the tropical central eastern Pacific may contribute to formation of anomalous northerly wind anomalies over the eastern part of the CLLH, leading to the cold meridional temperature advection and hence the negative winter SAT anomalies there. Li et al. (2023) also demonstrated that the increasing of the number of winter extreme low temperature days in the eastern part of the CLLH is closely related to the cold sea surface temperature anomalies (SSTA) over the central equatorial Pacific. However, both the winter El Niño and La Niña event may lead to cold extremes in the western part of the CLLH, depending on the different configurations between the ENSO and large-scale atmospheric circulation systems in the mid-high latitudes, such as the Arctic Oscillation (AO) and the Siberian High (Wan et al. 2023). Similar to the Precip, the relationship between winter SAT in the western part of the CLLH and the ENSO is non-stationary, which is weakened after the mid-1990s (Yuan et al. 2022). Besides the winter mean SAT, the ENSO has a profound impact on the winter synoptic temperature variability in the eastern part of the CLLH, and an interdecadal shift of this relationship occurred in the late 1980s, changing from significant to insignificant (Jian et al. 2021).

Fig. 1
figure 1

Location and Terrain height (unit: km) of the CLLH, the black box indicates the CLLH (21–29° N, 97–108° E)

The climatological early spring (March–April) Normalized Difference Vegetation Index (NDVI, see Sect. 2) over the CLLH is depicted in Fig. 2a. It is shown that the NDVI is over 0.5 in most area of the CLLH and featured by a north-to-south gradient, indicating relatively healthy and dense vegetation of this area. The annual mean value of area-averaged NDVI over the CLLH is 0.64, however the monthly area-averaged NDVI reaches its minimum value (0.58) in early spring throughout the year (black solid-dotted line in Fig. 2b), implying the lowest vegetation vigor and photosynthetic activity over the CLLH within the year. Meanwhile, we also calculated the autocorrelations of the CLLH monthly area-average NDVI between early spring and other months of the year (red dash-dotted line in Fig. 2b). The autocorrelation coefficients between early spring and March, April are 0.86 and 0.81, respectively, which are significant at the 99% confidence level, whereas those with other months are less than 0.26 and are not significant at the 90% confidence level. Besides, the autocorrelations of the CLLH monthly area-average NDVI between March, April and other months of the year are also calculated, which exhibit similar features to those between early spring and other months of the year (green dash-dotted line and blue dash-dotted line in Fig. 2b). The autocorrelation coefficient between March and April is 0.62 which is significant at the 99% confidence level, on the contrary the autocorrelation coefficients between March, April and other months of the year are weak and insignificant. These results demonstrate that the vegetation growth over CLLH in March and April can be studied as a whole and the interannual variability of early spring vegetation over the CLLH is different from the rest of the year. Hence early spring might be a unique season of vegetation growth over the CLLH.

Fig. 2
figure 2

A Climatological early spring CLLH NDVI, b Monthly area-average NDVI values (black solid-dotted line) over the CLLH (21–29° N, 97–108° E), autocorrelation coefficients of the CLLH monthly area-average NDVI between early spring (March–April, red dash-dotted line), March (green dash-dotted line), April (blue dash-dotted line) and other months of the year during 1982–2015. The big dot denotes the autocorrelation coefficient significant at the 95% confidence level

Efforts have been devoted to explore the potential impact of vegetation changes over the CLLH. Firstly, based on numerical simulations, He et al. (2007) reported that the CLLH deforestation will lead to a sharp decreasing of rain fall area and precipitation while temperature and drought index will rapidly increase, which would also intensify the risk of drought in the future (Zhao et al. 2013). Secondly, by utilizing the eddy covariance technique and the biometric-based method, Fei et al. (2017) found that the CLLH vegetation acts as an effective carbon sink, and the CLLH vegetation carbon sequestration potential exhibits an upward trend in the future (Lü et al. 2023). Besides, in early spring, the CLLH is one of the earliest atmospheric heating sources in the Asian summer monsoon region (Wen and Cao 2023; Yang et al. 2023a, b) further pointed out that this early spring diabatic heating source is mainly dominated by near-surface sensible heating, thus the CLLH vegetation variation has the potential to regulate weather and climate over the Asian summer monsoon region through modulation of surface energy partition (Lan et al. 2021). In a word, interannual variability of the CLLH vegetation can influence the terrestrial ecosystem and the atmosphere both locally and remotely, hence more in-depth researches are required.

Researches have been made to investigate the influence of local meteorological conditions on vegetation growth over the CLLH. For instance, in the western part of the CLLH, Shi and Chen (2018) showed that NDVI exhibits a closer relationship with the annual temperature than with precipitation. However, Yan et al. (2021) demonstrated that extreme precipitation has remarkable impact on net primary productivity of vegetation and its spatial distribution. Han et al. (2023) pointed that relative humidity is the main driving factor of NDVI, compared to other meteorological factors. In the eastern part of the CLLH, vegetation growth is thought to be governed by temperature (Hou et al. 2015; Xue et al. 2023).

Although located in the humid climate zone, in recent decades, the CLLH has suffered from severe and persistent droughts in early spring (Wang et al. 2015). Ding and Gao (2020) demonstrated that the extreme spring drought of 2019 in the western part of the CLLH is closely related with the abnormally weakened India-Burma Trough, which results in anomalous low-level moisture divergence. Moreover, the extreme spring drought of 2021 is caused by decaying phase of a La Niña event with negative geopotential height anomalies over the Philippine Sea, which exbibits distinctive characteristic from the historical perspective (Liu et al. 2022). Recently, some studies noticed that spring droughts have profound impact on vegetation over the CLLH. For example, Zhang et al. (2012) indicated that the severe and sustained spring drought of 2010 substantially reduced regional annual gross primary productivity and net primary productivity, and vegetation over the CLLH did not fully recover from drought stress until August. Lin et al. (2023) stated that extreme spring drought events not only have a delayed effect on the vegetation over the CLLH, but also affect its stability and resilience to future drought events. These studies implies that remote large-scale climate factors, such as tropical SST, can significantly modulate vegetation over the CLLH. However, few researches have focused on the effects of remote large-scale climate factors on vegetation interannual variation over the CLLH.

These motivate us to investigate the impact of tropical SSTAs on the interannual variability of early spring vegetation over the CLLH. We attempt to answer the following questions: Is tropical SSTAs related to the interannual variability of early spring vegetation over the CLLH? If so, what are the corresponding physical processes?

The remainder of the paper is organized as follows. Section 2 addresses the data and methods employed in this research. The major results, including the main mode of the early spring vegetation over the CLLH, the relationship between the late winter (January–February) SST seesaw of equatorial central Pacific–Philippine Sea and early spring vegetation interannual variation over the CLLH as well as the candidate physical mechanisms, are presented in Sect. 3. The physical mechanism of the later winter SST seesaw between the equatorial central Pacific and the Philippine Sea, the impacts of the ENSO and AO on early spring CLLH vegetation interannual variation, the contrast of impacts between the late winter SST seesaw and the ENSO on vegetation, are discussed in Sect. 4. Finally, conclusions are given in Sect. 5.

2 Data and methods

2.1 Data

The NDVI, a spectral vegetation index derived from the visible (red) spectral band (0.58–0.68 μm) and the near-infrared spectral band (0.725–1.10 μm) of the surface reflectance (Tucker 1979; Myneni et al. 1995), is calculated from the following formula:

$$\text{NDVI} = \frac{(\text{NIR-RED)}}{(\text{NIR+RED)}}$$
(1)

where RED and NIR are the reflectance values of the red and near-infrared bands, respectively. Providing an efficient way to describe vegetation phenology, NDVI is widely used as a remote sensing indicator for monitoring terrestrial vegetation greenness, vegetation coverage and vegetation productivity (Ichii et al. 2002; Pettorelli et al. 2005; Piao et al. 2011; Jiang et al. 2019). The NDVI data vary from – 1.0 to + 1.0, with positive values indicating areas covered by vegetation canopy, and negative values indicating ice and snow fields. High NDVI values (approximately 0.6–0.9) imply dense vegetation such as tropical forests or crops at their peak growth stage. In this study, the third generation NDVI biweekly data supplied by the Global Inventory Modeling and Mapping Studies (GIMMS NDVI3g) is adopted, with an original spatial resolution of 1/12° and a time span of 1982–2015. First of all, to minimize the interference from other factors unrelated to vegetation change (i.e., solar zenith angles, clouds, aerosols, etc.), the maximum value composite technique (Holben 1986) is utilized to obtain the monthly NDVI data. Then, the monthly NDVI data is further transformed into a coarser 0.5° × 0.5° spatial resolution to suppress the local noise and enhance the robustness of the large-scale climate signals in vegetation activity (Jiang et al. 2019; Ji and Fan 2020; Yao et al. 2021). Finally, pixels with a mean NDVI value less than 0.1 are excluded as non-vegetated or snow‐covered regions (Zhou et al. 2003; Wang et al. 2019; Ji and Fan 2020).

For the SST data, monthly SST used in this study is obtained from the Hadley Centre Sea Ice and Sea Surface Temperature dataset (HadISST), with a horizontal resolution of 1° × 1° (Rayner et al. 2003). Monthly SAT, Precip data over land are provided by the Climatic Research Unit of the University of East Anglia with a spatial resolution of 0.5° (CRU-TS v4.05, Harris et al. 2020). The monthly Global Precipitation Climatology Project (GPCP) Precip data at a 2.5° × 2.5° horizontal resolution is also adopted (Adler et al. 2003).

For land surface variables, we use the monthly data from the land component of the fifth generation of European Center for Medium-Range.

Weather Forecasts Re-Analysis (ERA5-Land, Hersbach et al. 2020), with a spatial resolution of 0.5°, including volumetric soil moisture (SM) and soil temperature (ST) of three layers (0–7, 7–28, 28–100 cm depth) as well as surface heat fluxes (including latent heat flux, downward shortwave radiation flux). Due to local terrain complexity of the CLLH, the monthly SM and ST provided by the Global Land Data Assimilation System Noah model (GLDAS-Noah, Rodell et al. 2004) are also employed for mutual verification, with a spatial resolution of 0.25° and three layers of 0–10, 10–40, 40–100 cm depth as well as surface heat fluxes. Besides, the total SM content and ST of top 100 cm are also utilized, which are summed up and weighted by the thickness of each layer. The monthly total cloud cover data is derived from CRU-TS v4.05.

For atmospheric circulation data, monthly ERA5 reanalysis (Hersbach et al. 2020) is used to reveal the physical processes induced by tropical SSTAs, with a horizontal resolution of 1° × 1°, including zonal and meridional winds, geopotential height, vertical velocity, specific humidity and temperature on pressure levels as well as sea level pressure and surface pressure on single level. Top shortwave and longwave radiation flux, sea level pressure, surface pressure, surface zonal and meridional wind stress on single level with the same resolution are also used.

Studies have well documented the substantial influences of ENSO on vegetation growth (e.g. Erasmi et al. 2009; Hao et al. 2020; Li et al. 2017; Wang et al. 2022a, b). Hence, the potential impact of late winter ENSO on the early spring CLLH NDVI is investigated in this study. The ENSO variability is measured by different Niño indices. The Niño 3 index (INiño3), Niño 3.4 index (INiño3.4) and Niño 4 index (INiño4) are defined as the area-averaged SSTAs in the Niño 3 region (5°S–5° N, 150°–90° W), Niño 3.4 region (5° N–5°S, 170° W–120° W) and Niño 4 region (5°–5° N, 160°–150° W), respectively. Furthermore, the index of Central-Pacific (CP) type of ENSO (ICP_ENSO) and the index of Eastern-Pacific (EP) type of ENSO (IEP_ENSO) are obtained by a combination of the Niño3 and Niño4 indices proposed by Ren and Jin (2011). Besides the ENSO, the potential impact of AO, which is the main mode of the variability over the northern Hemisphere extratropical region, is also considered. The monthly mean AO index (IAO) is adopted from the website of Climate Prediction Center (https://www.cpc.ncep.noaa.gov/products/precip/CWlink/daily_ao_index/monthly.ao.index.b50.current.ascii).

In this study, early spring is defined as March–April and late winter as January–February. The target region for the NDVI is limited to the CLLH (21–29° N, 97–108° E). All the monthly data sets employed in this study share a common time period of 1982–2015 (the same time span as the NDVI data), with linear trends removed prior to the analysis, which is expected to eliminate the potentially deceptive relationships caused by long-term changes in large-scale climate factors and vegetation indices.

2.2 Methods

Statistical methods, including empirical orthogonal function (EOF) analysis, correlation and partial correlation analysis, linear regression analysis and the Liang-Kleeman information flow (IF, Liang 2014) analysis, are employed in this study.

The EOF analysis is carried out to acquire the leading mode of early spring NDVI over the CLLH. Regression and correlation analyses are utilized to investigate the relationship between early spring CLLH NDVI and meteorological factors of interest, such as late winter tropical SST.

To examine the independence of this relationship from other different meteorological factors, the partial correlation coefficients are calculated as follows:

$${r}_{ab,c}=\frac{{r}_{ab}-{r}_{ac}\times {r}_{bc}}{\sqrt{(1-{{r}_{ac}}^{2})\times (1-{{r}_{bc}}^{2})}}$$
(2)

where rab,c is the partial correlation coefficient of a and b independent of c, rab, rac and rbc indicates the correlation coefficient of a and b, a and c, and b and c, respectively. Partial regression coefficients could be obtained in a similar way to isolate impacts of different factors.

The high correlation does not suggest a causal relationship. Strong correlation often exists between two variables which are regulated by the same external factor while causal relationship does not essentially exist. Two events are causally connected if and only if they are connected by underlying physical mechanisms. Therefore, the causalities between the NDVI and local meteorological factors are analyzed using the Liang-Kleeman IF method quantitatively where causality is measured through the time rate of information flowing from one series to the other, to verify whether changes in local meteorological factors indeed cause variation in the NDVI. From Liang (2014), under the assumption of linear model, the maximum likelihood estimation form of the information flow from time series X2 (i.e., local meteorological factors) to time series X1 (i.e., the NDVI) can be defined as:

$${T}_{2\to 1}=\frac{{C}_{11}{C}_{12}{C}_{2,d1}-{{C}_{12}}^{2}{C}_{1,d1}}{{{C}_{11}}^{2}{C}_{22}-{C}_{11}{{C}_{12}}^{2}}$$
(3)

where Cij is the sample covariance between Xi and Xj and Ci,dj is the sample covariance between Xi and (Xj,n+1Xj,n)/ dt (Euler forward scheme), with dt indicating the time interval. According to the Liang-Kleeman IF theory, \({T}_{2\to 1}\) can be zero or nonzero, when \(|{T}_{2\to 1}|>0\), X2 is causal to X1, and vice versa. Apparently, two uncorrelated variables (C12 = 0) must be noncausal (\({T}_{2\to 1}=0\)), while the converse does not hold (when \({T}_{2\to 1}=0\), C12 may not be zero). To sum up, correlation does not imply causation, however, causation implies correlation. More detailed derivation and formula verification can be found in the research of Liang (2014).

Relying on the diagnosing of quasi-geostrophic omega equation, Feng et al. (2014) pointed out that the abnormal decent motion over the CLLH during winter drought is primarily maintained by the anomalous cold temperature advection processes in La Niña years. However, quasi-geostrophic omega equation is not applicable in the low latitude region considering the fact that the Coriolis force is relatively small in this area. In the tropical and monsoon regions, vertical motion is constrained by the moist static energy (MSE) budgets (Neelin and Held 1987), hence the MSE budget equation is diagnosed in this study to explore the mechanisms responsible for the anomalous vertical motion over the CLLH.

The MSE is the sum of geopotential energy, internal energy and latent heat energy:

$$MSE \, = \, gz + C_{P} T + L_{v} q$$
(4)

where g is the gravitational acceleration, z is height, T is the temperature, q is the specific humidity, CP and Lv are the specific heat at constant pressure and the latent heat of vaporization, respectively. h = CPT + Lvq is the moist enthalpy. The anomalous column-integrated MSE budget is given by:

$$\left\langle {\frac{\partial }{{\partial t}}MSE} \right\rangle ^{\prime } + \left\langle {\vec{V} \cdot \nabla MSE} \right\rangle ^{\prime } + \left\langle {\omega \frac{\partial }{{\partial p}}MSE} \right\rangle ^{{\prime }} = F_{{net}}^{'}$$
(5)

\(\overrightarrow{V}\) is the horizontal wind velocity, \(\omega\) is the vertical wind velocity, \(\nabla\) is the horizontal gradient, Fnet denotes the total of surface and top-of-atmosphere radiative fluxes and surface turbulent fluxes entering into the atmospheric column. The primes denote the deviation from monthly climatology, the angle brackets denote the mass integral through an entire atmospheric column (from surface to 100 hPa in this study).

On the interannual time scale, the tendency term \(\left\langle {\frac{\partial }{{\partial t}}MSE} \right\rangle ^{\prime }\) is relatively smaller than other terms and can be neglected and the advection terms can be decomposed into linear and nonlinear terms (Wu et al. 2017; Hu et al. 2021; Wang et al. 2022b; Lu et al. 2023), hence Eq. (5) can be simplified to:

$$\underbrace {{\left\langle {\omega ^{\prime}\frac{\partial }{{\partial p}}\overline{{MSE}} } \right\rangle }}_{A} = \underbrace {{ - \left\langle {\vec{V}^{\prime} \cdot \nabla \bar{h}} \right\rangle }}_{B}\underbrace {{\; - \left\langle {\overline{{\vec{V}}} \cdot \nabla h^{\prime}} \right\rangle }}_{C}\;\;\underbrace {{ - \left\langle {\bar{\omega }\frac{\partial }{{\partial p}}MSE^{\prime}} \right\rangle }}_{D} + \underbrace {{F_{{net}}^{'} }}_{E} + \underbrace {{NL}}_{F}$$
(6)

in which the NL denotes the sum of all nonlinear and transient terms, overbars denote the monthly climatology. Furthermore, Term C can be decomposed as follows:

$$\underbrace {{ - \left\langle {\overline{{\vec{V}}} \cdot \nabla h^{\prime}} \right\rangle }}_{C} = \underbrace {{ - \left\langle {\overline{{\vec{V}}} \cdot C_{p} \nabla T^{\prime}} \right\rangle }}_{{C_{1} }}\;\underbrace {{ - \left\langle {\overline{{\vec{V}}} \cdot L_{v} \nabla q^{\prime}} \right\rangle }}_{{C_{2} }}$$
(7)

Surface wind stress acting on the sea surface is regarded as a fundamental forcing for upper ocean circulation. Surface wind stress curl is calculated as follows:

$$\text{Curl(}\tau )=\frac{\partial {\tau }_{y}}{\partial x}-\frac{\partial {\tau }_{x}}{\partial y}$$
(8)

where \(\text{Curl(}\tau )\), \({\tau }_{x}\), \({\tau }_{y}\) is the surface wind stress curl, surface zonal and meridional wind stress, respectively.

The significance level of the regression and correlation coefficients is estimated by Student's t-test, and that of the Liang-Kleeman IF is based on Z‐test.

3 Results

3.1 Statistical linkage between late winter SST seesaw of equatorial central Pacific–Philippine Sea and early spring CLLH NDVI

  1. a.

    The dominant mode of early spring NDVI anomalies

We first extract the dominant mode of early spring NDVI interannual variation over the CLLH for the period 1982–2015 and the corresponding principal component (PC) time series via the EOF method. Figure 3a displays the spatial distribution of NDVI anomalies related to the first EOF mode (EOF1). To assess the significance of the result, the spatial pattern of EOF1 is obtained by regressing NDVI anomalies onto the corresponding normalized PC (PC1, green solid-dotted line in Fig. 3d). The EOF1 explains 33.6% of the total interannual variance of early spring CLLH NDVI anomalies and is well separated from the remainder of the EOF modes according to the criterion proposed by North et al. (1982). As the percentage variances of other EOF modes are relatively small, this study focuses on analyzing the EOF1.

Fig. 3
figure 3

Early spring NDVI anomalies over the CLLH obtained by regression on a the normalized principal components (PC1) corresponding to EOF1, b Iseesaw in late winter during 1982–2015. c as (b), but for the Liang-Kleeman IFs (unit: nat·month−1, divided by 10 for better presentation) from Iseesaw. e The time series of PC1 (green solid-dotted line), Iseesaw (blue solid-dotted line) in late winter. Stippling in (ac) indicates where the NDVI anomalies are significant at the 95% confidence level. The blue and red boxes in b and c are used to define the western and eastern part of the CLLH respectively as described in Sect. 3.2b

The EOF1 is obviously characterized by same-sign NDVI anomalies over the CLLH (Fig. 3a). In the positive EOF1 phase, significant positive NDVI anomalies greater than 0.02 are located in most of the CLLH, and vice versa. The PC1 exhibits significant interannual variations (green solid-dotted line in Fig. 3d), and the correlation coefficient between PC1 and early spring CLLH NDVI index (defined as the area-average NDVI during 1982–2015) is 0.98 (exceeding the 99% confidence level), indicating that PC1 can well represent early spring NDVI interannual variability in most of the CLLH.

  1. b.

    Related later winter tropical SSTAs

To reveal the statistical linkage between early spring CLLH NDVI and late winter tropical SST, Fig. 4 presents the later winter SSTAs regressed onto the PC1. In the positive EOF1 phase, significant warm SSTAs are observed over the equatorial central Pacific, accompanied with significant cold SSTAs over the Philippine Sea. These features imply that increased NDVI appears in the CLLH during early spring while late winter warm (cold) SSTAs appear in the equatorial central Pacific (Philippine Sea), and vice versa.

Fig. 4
figure 4

Regression pattern of late winter SSTA (unit: °C) with respect to PC1 during 1982–2015. The stippled regions indicate SSTA significant at the 95% confidence level. The black rectangles outline the regions (3°S–4° N, 168° E–170° W) and (12–19° N, 126–155° E) to define the CPI and PSI indices, respectively

Based on the above analyses, the equatorial central Pacific index (CPI) is defined as the normalized area-averaged SSTAs in the equatorial central Pacific (3°S–4° N, 168° E–170° W) and the Philippine Sea index (PSI) is defined as the normalized area-averaged SSTAs in the Philippine Sea (12°–19° N, 126°–155° E) during 1982–2015 in late winter, where the SSTAs of these regions exceed the 95% confidence level as shown in Fig. 4. The correlation coefficient between PC1 and CPI (PSI) is 0.36 (– 0.46), exceeding the 95% (99%) confidence level. These result indicate an in-phase (out-of-phase) relationship between early spring CLLH NDVI and late winter SST over equatorial central Pacific (Philippine Sea).

There exist, indeed, high negative correlations between late winter SST over equatorial central Pacific and Philippine Sea. Here we calculate the correlation coefficient between CPI and PSI and find it to be – 0.77 during 1982–2015 (exceeding the 99% confidence level), which suggests a possible close out-of-phase variation of late winter SST between equatorial central Pacific and Philippine Sea. Besides, the partial correlation coefficient between PC1 and CPI (PSI) independent of the PSI (CPI) is 0.27 (– 0.01), which is weak and not significant at the 90% confidence level. These results imply that it is the synergetic effect of late winter SST over equatorial central Pacific and Philippine Sea, rather than the individual effect, is tightly linked to the early spring CLLH NDVI interannual variability. We also calculated the Liang-Kleeman IFs between the CPI and the PSI, with the IF from the CPI to the PSI being – 0.65 nat/month (exceeding the 99% confidence level) while the IF from the PSI to the CPI being 0.08 nat/month (not significant at the 90% confidence level), suggesting that the CPI (PSI) exerts substantial (weak) influence on the PSI (CPI). Hence there is a causal relationship between the PSI and the CPI rather than the spurious correlation relationship. The physical mechanism that the PSI influences the CPI will be discussed in Sect. 4.

Therefore, we define the index of late winter SST seesaw between equatorial central Pacific and Philippine Sea (Iseesaw, blue solid-dotted line in Fig. 3d) as the normalized times series of CPI minus PSI. The Iseesaw is highly correlated with CPI (PSI) with a correlation coefficient of 0.94 (– 0.94). This indicates that the Iseesaw well captures the late winter SST covariation over equatorial central Pacific and Philippine Sea.

The correlation coefficient between PC1 and Iseesaw is 0.44 during 1982–2015 (exceeding the 99% confidence level), suggesting a notably coherent variation between late winter SST seesaw of equatorial central Pacific–Philippine Sea and early spring CLLH NDVI interannual variability. Meanwhile, regression maps of early spring CLLH NDVI anomalies onto the Iseesaw are shown in Fig. 3b. In the positive phase of Iseesaw, positive NDVI anomalies greater than 0.01 are located in most of the CLLH, which bear a close resemblance to the counterpart related to the PC1 although the Iseesaw related anomalies exhibit less magnitude and smaller significant area (Fig. 3a). Meanwhile, the Liang-Kleeman IF analysis is adopted to clarify the causal effects of late winter Iseesaw on early spring CLLH NDVI. It is depicted in Fig. 3c that obvious IFs locate in the western and eastern part of the CLLH associated with strong and significant NDVI anomalies (Fig. 3b). We also examined the Iseesaw regressed NDVI anomalies over the CLLH in March and April respectively (Supplementary Fig. S1). In the positive phase of Iseesaw, evident positive NDVI anomalies exist in the CLLH in both March and April. It is also noted that the positive NDVI anomalies are strengthened from March to April.

Finally, the possible impact of preceding SSTAs on early spring CLLH NDVI is examined as well. Correlation coefficients between the PC1 and the preceding monthly Iseesaw imply that the Iseesaw in January and February significantly influence the PC1 (0.39 and 0.47 respectively, exceeding the 99% confidence level) whereas the Iseesaw in other months has little impact on the PC1 (Supplementary Fig. S2). The preceding globally seasonal SSTAs regressed on the PC1 are presented in Fig. S3. The SSTAs over the Atlantic and the Indian Ocean are weak and insignificant in all seasons, while evident and significant SSTAs appear over the pacific in winter, indicating that the late winter seesaw SST pattern is the most critical physical factor affecting the early spring CLLH NDVI. Moreover, significant positive SSTAs over the equatorial central Pacific have already emerged in autumn and are intensified in winter. This result suggests the potential impact of ENSO on the early spring CLLH NDVI, which will be discussed in Sect. 4.

To sum up, the above results clearly confirm that late winter SST seesaw of equatorial central Pacific–Philippine Sea is significantly associated with the same-sign pattern of early spring CLLH NDVI. When the later winter Iseesaw is higher than normal, the NDVI in the following early spring tends to increase in most of the CLLH, and vice versa.

3.2 The physical mechanisms responsible for the linkage between late winter SST seesaw of equatorial central Pacific–Philippine Sea and early spring CLLH NDVI: a local perspective

In the previous section, we have demonstrated the statistical linkage between the late winter SST seesaw and early spring CLLH NDVI. Next, we will identify the key climate factors related to early spring CLLH NDVI and establish linkages with the late winter SST seesaw from the local perspective.

a. The key climate factors related to early spring CLLH NDVI

Before identifying the key climate factors related to early spring CLLH NDVI, the land surface condition climatology over the CLLH in early spring is examined at first (Fig. 5). The SAT increases from north to south (Fig. 5a), so does the ST (Fig. 5b, c). For the Precip and SM (Fig. 5d–f), the CLLH is characterized by an obvious west-to-east gradient, with the eastern (western) part being relatively wet (dry). The cloud cover significantly increases from west to east (Fig. 5g), correspondingly land surface downward shortwave radiation flux significantly decreases from west to east (Fig. 5h, i). This climatological land surface condition contrast between east and west of the CLLH is mainly induced by the Yunnan-Guizhou quasi-stationary front, which frequently occurs during the winter half year over the CLLH and is characterized by severe cold, cloudy weather to its east but mild sunny weather to its west (Cai et al. 2022). Hence, the land surface condition climatology over the CLLH in early spring is characterized by the western part being relatively warm-and-dry while the eastern part being cool-and-wet.

Fig. 5
figure 5

Early spring Climatological (1982–2015) a SAT (unit: °C), b top 100 cm ST derived from ERA5-Land (unit: K), c top 100 cm ST derived from GLDAS-Noah (unit: K), d Precip (unit: mm·month−1), e top 100 cm volumetric SM derived from ERA5-Land (unit: m3·m−3), f top 100 cm SM content derived from GLDAS-Noah (unit: kg·m−2), g cloud cover (unit: %), h surface downward shortwave radiation flux derived from ERA5-Land (unit: W·m−2), i surface downward shortwave radiation flux derived from GLDAS-Noah (unit: W·m−2). The Precip, SAT and cloud cover data use the CRU-TS v4.05 dataset

Then, Correlation and IF analyses are employed to explore the key climate factors related to early spring CLLH NDVI. First of all, several relevant local climate factors are identified through the investigation of correlation with PC1, including Precip, SAT, total SM and ST of top 100 cm (Fig. 6). Notably, the correlation between PC1 and simultaneous Pre as well as SAT are very weak and insignificant (Fig. 6a, b). In sharp contrast, PC1 is significantly and positively correlated with late winter Pre in most of the CLLH and SAT in the eastern part of the CLLH (Fig. 6c, d). The correlation maps of early spring SM (ST) with PC1 (Fig. 6e–h) are very similar to those of late winter Pre (SAT). Yang and Zhang (2016) pointed out that the memory length of winter ST anomaly in the eastern part of the CLLH is longer than 2 months, Gao et al.(2017) reported that the CLLH autumn SM anomaly exhibits strong memory that can persist to the following spring. The above results indicate that more vigorous vegetation over the CLLH in early spring is closely related with a wetter land surface condition over most of the CLLH and a warmer land surface condition over the eastern part of the CLLH simultaneously, which may result from stronger-than-normal Precip in most of the CLLH and higher-than-normal SAT in the eastern part of the CLLH in late winter. Additionally, the PC1 shows weak and insignificant correlation with early spring total cloud cover as well as surface downward shortwave radiation flux (Supplementary Fig. S4a,b,c), which suggests that early spring incoming solar radiation is not tightly linked with CLLH NDVI interannual variation.

Fig. 6
figure 6

Correlations of PC1 with early spring a Precip, b SAT, e top 100 cm SM derived from ERA5-Land, f top 100 cm ST derived from ERA5-Land, g top 100 cm SM derived from GLDAS-Noah, h top 100 cm ST derived from GLDAS-Noah and with late winter c Precip, d SAT during 1982–2015. The Precip and SAT data use the CRU-TS v4.05 dataset. The stippled regions indicate correlation coefficients significant at the 95% confidence level

Is interannual variation of early spring CLLH NDVI intrinsically induced by the above anomalies? To answer this question, the Liang-Kleeman IF analysis is applied to further clarify the causal effects of local climate factors on early spring CLLH NDVI in Fig. 7. For early spring Precip and SAT, the IFs are almost zero (Fig. 7a, b) associated with weak correlations (Fig. 7a, b), which confirms that changes in early spring Precip and SAT cannot lead to changes in early spring CLLH NDVI. However, significant nonzero IFs of late winter Precip and SAT can be detected (Fig. 7c, d), suggesting potential impacts of late winter Precip and SAT on early spring CLLH NDVI. Specifically, it is noted that late winter Precip shows weak and insignificant IFs (Fig. 7c) although exhibits strong and significant correlations (Fig. 6c) in the relatively cool-and-wet eastern part of the CLLH, where the correlations and IFs of early spring SM display similar results (Fig. 6e, g and 7e, g). However, late winter SAT shows significant correlations and IFs in the relatively cool-and-wet eastern part (Fig. 6d and Fig. 7d), so does early spring ST (Fig. 6f, h and 7g, h). In the relatively warm-and-dry western part of the CLLH, late winter Precip and early spring SM demonstrate strong and significant IFs associated with significant positive correlations (Fig. 6c, e, g and 7c, e, g). Therefore, the above results indicate that there are causal influences from early spring SM (ST) and late winter Precip (SAT) to early spring vegetation growth in the relatively warm-and-dry (cool-and-wet) western (eastern) part of the CLLH.

Fig. 7
figure 7

Same as in Fig. 6, but for the Liang-Kleeman IFs (unit: nat·month−1) to PC1

The IFs of early spring total cloud cover and surface downward shortwave radiation flux to PC1 are also investigated, the results of which are negligible (Supplementary Fig. S4d–f) corresponding to the weak and insignificant correlations between PC1 and total cloud cover as well as surface downward shortwave radiation flux. Climatologically, the monthly area-averaged cloud cover (57.8%, Fig. 5g) over the CLLH approximates to its minimum in early spring (Supplementary Fig. S5), meanwhile the monthly area-averaged surface downward shortwave radiation flux (171.6 and 182.5 W·m-2 for ERA5-Land and GLDAS-Noah respectively, Fig. 5h, i) reaches its maximum in early spring (Supplementary Fig. S5). These results suggest that there is sufficient incoming solar radiation over the CLLH in early spring for vegetation growth. As a result, variations of in situ incoming solar radiation may have relatively little impact on CLLH NDVI in early spring compared with other climate factors.

To brief, the above correlation and Liang-Kleeman IF analyses indicate that early spring interannual variabilities of NDVI in the relatively warm-and-dry (cool-and-wet) western (eastern) part of the CLLH is mainly modulated by local SM (ST) anomalies, which may be provoked by Precip (SAT) variation in late winter. The abnormal signals of late winter Precip (SAT) may be recorded by land surface and may persist to early spring owing to the SM (ST) memory. Accompanied with sufficient incoming solar radiation, wetter (warmer)-than-normal land surface condition in the relatively warm-and-dry (cool-and-wet) western (eastern) part of the CLLH is favorable for in situ vegetation growth in early spring. Hence the local SM and ST might act as an important bridge that links late winter atmospheric circulations to the following early spring NDVI over the CLLH.

b. The key climate factors affected by later winter SST seesaw of equatorial central Pacific–Philippine Sea

The next question then, is what is the relationship between the late winter SST seesaw and the key climate factors that influence the early spring CLLH NDVI. To answer this question, we regress later winter and early spring key climate factors on late winter Iseesaw and calculate the IFs from the Iseesaw over the CLLH during 1982–2015.

For the Precip and SAT anomalies in late winter, positive and significant Pre anomalies around 0.3 mm/day are observed over most of the CLLH (Fig. 8a) accompanied by pronounced IFs (Fig. 9a). Meanwhile, significant increased SAT anomalies around 0.8 °C exist in the eastern part of the CLLH (Fig. 8b), alongside with apparent IFs there (Fig. 9b). These results demonstrate that the late winter SST seesaw of equatorial central Pacific–Philippine Sea indeed has a profound impact on CLLH, and the related physical mechanisms will be analyzed in Sect. 3.3.

Fig. 8
figure 8

Regression patterns of a Precip anomalies (unit: mm·day−1) and b SAT anomalies (unit: °C) with respect to Iseesaw in late winter during 1982–2015. The Precip and SAT data are derived from GPCP and CRU-TS v4.05, respectively. The black box indicates the CLLH (21–29° N, 97–108° E). The blue box indicates the southern Bay of Bengal (3–15° N, 75–100° E) as described in Sect. 3.3. The stippled regions indicate regression coefficients significant at the 95% confidence level

Fig. 9
figure 9

The Liang-Kleeman IFs (unit: nat·month−1) from Iseesaw to a Precip, b SAT in late winter and c top 100 cm SM, d ST derived from ERA5-Land, e top 100 cm SM, f ST derived from GLDAS-Noah in early spring during 1982–2015. The stippled regions indicate the IFs significant at the 95% confidence level

For the SM anomalies, significant increased SM appears in the eastern part of the CLLH in late winter (Fig. 10a, b). The wet soil anomalies are in fact associated with the increased rainfall corresponding to the positive phase of previous early winter Iseesaw (i.e., previous November–December, Supplementary Fig. S6), however there is no significant change of SM in the western part of the CLLH corresponding to the weak Precip anomalies there (Supplementary Fig. S6). In early spring, positive SM anomalies in the eastern part of the CLLH disappears, while the SM anomalies show overall wetness in the western part of the CLLH. The positive SM anomalies in the western part of the CLLH can be viewed as the late winter Precip signal transmission in the soil, hence benefit in situ vegetation growing (Fig. 10c, d). Besides, the IFs confirms the evident impact of the Iseesaw on the SM via Precip (Fig. 9c, e).We also calculated the surface latent heat flux anomalies, the results of which exhibit positive (weak) and significant (insignificant) anomalies in the eastern (western) part of the CLLH (Fig. 17c) corresponding to the local positive (weak) and significant (insignificant) SAT anomalies (Fig. 8b). These results imply that in the western part of the CLLH, positive Precip anomalies in later winter are stored and accumulated as wet SM anomalies to early spring, while in the eastern part of the CLLH the wet SM anomalies in late winter cannot persist to early spring due to the local enhanced evapotranspiration.

Fig. 10
figure 10

Regression patterns of SM anomalies derived from ERA5-Land (unit: m3·m−3, multiplied by 103 for better presentation) in a late winter and b early spring with respect to late winter Iseesaw during 1982–2015. c and d are the same as in (a, b), but the SM anomalies are derived from GLDAS-Noah (unit: kg·m−2). The stippled regions indicate regression coefficients significant at the 95% confidence level

For the ST anomalies, Fig. 11a, b depict that there are positive ST anomalies in the eastern part of the CLLH corresponding to positive SAT anomalies in late winter. These positive ST anomalies persist to early spring (Fig. 11c, d), favoring local vegetation growing in the eastern part of the CLLH. As illustrated in Fig. 9d, f, obvious IFs exist over the eastern part of the CLLH, demonstrating the impact of the Iseesaw on the ST via SAT. Besides, the enhanced in situ surface latent heat flux (Fig. 17) in late winter could reduce the infiltration of surface rainfall into the subsurface soil, resulting in decreased soil heat capacity and thus contributing to the maintenance of positive ST anomalies in early spring.

Fig. 11
figure 11

Regression patterns of ST anomalies derived from ERA5-Land (unit: K, multiplied by 102 for better presentation) in a late winter and b early spring with respect to late winter Iseesaw during 1982–2015. c, d Same as in a and b, but the ST anomalies are derived from GLDAS-Noah (unit: K, multiplied by 102 for better presentation). The stippled regions indicate regression coefficients significant at the 95% confidence level

In addition to the total SM and ST of top 100 cm, we also investigate at which soil depths the later winter SST seesaw can impact the SM (ST) in the western (eastern) part of the CLLH. The area-average SM (ST) anomalies of different soil layers from January to April are calculated in a rectangular blue (red) box in 23°–25.5° N (23°–29° N) and 97°–102° E (104°–108° E) as shown in Fig. 3b, which is selected as a typical region of western (eastern) part of the CLLH where in situ vegetation growth is modulated by SM (ST) as described in Sect. 3.2a. Figure 12a, b indicate that positive Pre anomaly in the western part of the CLLH penetrates to the surface soil layers in January and February, including layer L1 (0–7 cm and 0–10 cm for ERA5-Land and GLDAS-Noah, respectively) and L2 (7–28 cm and 10–40 cm for ERA5-Land and GLDAS-Noah, respectively), however the response of SM in the deep layer L3 (28–100 cm and 40–100 cm for ERA5-Land and GLDAS-Noah, respectively) is very weak. These may explain the insignificant change of total SM of top 100 cm in late winter (Fig. 10a, b). In March and April, the influence of positive Pre anomaly is no longer confined to surface soil layers but reaches deep layer L3, therefore positive total SM anomalies of top 100 cm are observed (Fig. 10c, d). In the eastern part of the CLLH, the ST anomalies exhibit similar results with the SM anomalies in the western part of the CLLH (Fig. 13a, b), which indicate similar evolution of SAT anomaly signal in the soil layers. Unlike the positive Pre anomalies, the influence of positive SAT anomalies already reaches the deep layer L3 in February, corresponding to the slightly warmer-than-normal total ST of top 100 cm in late winter (Fig. 11a, b). The above results are further verified by Liang-Kleeman IF analysis. For the IFs of SM in the western part of the CLLH (Supplementary Fig. S7), obvious IFs appear in surface layer L1 and subsurface L2 in February and reach deep layer L3 in March, corresponding to the significant positive SM anomalies through the entire soil layer. In April, the IFs in the surface layer disappear, however, they are still evident in subsurface layer and deep layer, so does the positive SM anomalies. The IFs of ST in the eastern part of the CLLH demonstrate similar features (Supplementary Fig. S8). The above results indicate the persistence of SM (ST) signal from March to April in the relatively warm-and-dry (cool-and-wet) western (eastern) part of the CLLH, and may partially explain the strengthened positive NDVI anomalies from March to April (Supplementary Fig. S1).

Fig. 12
figure 12

Regression of area-average SM anomalies over the western part of the CLLH (23°–25.5° N, 97°–102° E, rectangular blue box in Fig. 3b) a derived from ERA5-Land (unit: m3·m−3, multiplied by 103 for better presentation) and b derived from GLDAS-Noah (unit: kg·m−2) in different layers from January to April with respect to late winter Iseesaw during 1982–2015. The symbol “dot” denotes the significant point at the 95% confidence level

Fig. 13
figure 13

Regression of area-average ST anomalies over the eastern part of the CLLH (23°–29° N, 104°–108° E, rectangular red box in Fig. 3b) a derived from ERA5-Land (unit: K, multiplied by 102 for better presentation) and b derived from GLDAS-Noah (unit: K, multiplied by 102 for better presentation) in different layers from January to April with respect to late winter Iseesaw during 1982–2015. The symbol “dot” denotes the significant point at the 95% confidence level

To sum up, the above analyses imply that the vertical penetration of SM (ST) anomalies from January to April can be perceived as the transport of late winter Pre (SAT) signals from surface to deep layers in the relatively warm-and-dry (cool-and-wet) western (eastern) part of the CLLH, suggesting the restoration of later winter SST seesaw signals via the local soil memory effect. As a result, the interannual variation of early spring CLLH NDVI is significantly regulated by later winter SST seesaw of equatorial central Pacific–Philippine Sea.

3.3 The physical mechanisms responsible for the linkage between late winter SST seesaw of equatorial central Pacific–Philippine Sea and early spring CLLH NDVI: atmospheric circulation anomalies

The Pre and SAT anomalies over the CLLH in late winter are strongly affected by atmospheric circulation. In this section, we explore the circulation anomalies associated with the SST seesaw of equatorial central Pacific–Philippine Sea in late winter.

Moisture transport and upward motion play essential role in the formation of the Precip. As depicted in Fig. 15c, there is no evident vertically integrated water vapor flux into the CLLH through its boundary, suggesting the relatively weak impact of the Iseesaw-related moisture transport to the increased Precip over the CLLH in the positive phase of Iseesaw. Previous studies have elucidated that sever droughts and high-temperature events over the CLLH in winter and spring are primarily attributed to anomalous subsidence in middle and upper troposphere (Feng et al. 2014; Ding and Gao 2020; Wang and Kong 2021; Dong et al. 2023). Hence, the mechanisms responsible for the Iseesaw-related vertical motion over the CLLH are investigated through the MSE budget equation (Eq. (6)).

In the positive phase of Iseesaw, significant anomalous ascending motions dominate the troposphere over the CLLH (blue solid-dotted line in Fig. 14a), indicating \(< \omega ^{\prime} >\) is less than 0. On the other hand, the climatological MSE increases with height (red solid line in Fig. 14a), i.e., \(\frac{\partial }{\partial p}\overline{MSE}\) is less than 0. As a result, anomalous positive \(<\omega^{\prime}\frac{\partial }{\partial p}\overline{MSE}>\) appears (term A of Eq. (6) and Fig. 14b). Therefore, if the physical processes on the right-hand side of Eq. (6) increase the MSE in the column, anomalous ascending motions should be generated to increase the MSE exported out of the column to keep balance of the MSE budget, and vice visa.

Fig. 14
figure 14

Regression of a CLLH area-average vertical velocity anomalies profile (blue solid-dotted line, unit: Pa·s−1, multiplied by 103 for better presentation), b CLLH area-average MSE budget (unit: W·m−2), c CLLH area-average anomalous enthalpy advection of anomalous temperature by climatological wind profile (red solid-dotted line, unit: W·m−1, multiplied by 103 for better presentation) and d SBOB area-average (3–15° N, 75–100° E, blue box in Fig. 8a) temperature anomalies profile (red solid-dotted line, unit: K) and vertical velocity anomalies profile (blue solid-dotted line, unit: Pa·s−1, multiplied by 103 for better presentation) with respect to Iseesaw in late winter during 1982–2015. The red solid line in a denotes the CLLH area-average climatological MSE (unit: W·m−2, divided by 103 for better presentation). The big dots denote the area-average anomalies significant at the 95% confidence level

The results of MSE budget (Fig. 14b) revealed that the positive \(<\omega {\prime}\frac{\partial }{\partial p}\overline{MSE}>\) are primarily balanced by horizontal advection of anomalous moist enthalpy by climatological horizontal winds \(-<\overline{\overrightarrow{V}}\cdot \nabla h^{\prime}>\) (term C of Eq. (6) and Fig. 14b) and positive net energy flux anomalies \(F_{net}^{\prime}\)(term E of Eq. (6) and Fig. 14b) in the positive phase of Iseesaw, with the contribution of \(-<\overline{\overrightarrow{V}}\cdot \nabla h^{\prime}>\) being 7 times larger than \({F}_{net}^{\prime}\), hence the impact of positive net energy flux anomalies is rather small and will not be analyzed in this research. To investigate the detailed process of the anomalous advection of moist enthalpy by climatological winds, a further decomposition is applied (Eq. (7)). The results show that the anomalous positive enthalpy advection by climatological wind \(-<\overline{\overrightarrow{V}}\cdot {C}_{p}\nabla T^{\prime}>\) (term C1 of Eq. (7) and Fig. 14b) has the largest contribution. The vertical profile of anomalous enthalpy advection by climatological wind \(-\overline{\overrightarrow{V}}\cdot {C}_{p}\nabla T^{\prime}\) is shown in Fig. 14c, the results indicate that significant positive enthalpy advection exists in middle and upper troposphere with maximum value at 350-hPa level, corresponding to pronounced anomalous upward motion there (Fig. 14b).

In the positive phase of Iseesaw, anomalous positive enthalpy advection in middle and upper troposphere contributes to increased enthalpy entering into the atmospheric column over the CLLH, an anomalous compensating upward motion is then generated to increase the MSE exported out of the atmospheric column and thus keep balance of the MSE budget. Hence, the Iseesaw-related anomalous positive enthalpy advection by climatological wind in middle and upper troposphere is the dominant driver of the anomalous upward motion over the CLLH. Feng et al. (2014) also demonstrated that during the winter drought, the local anomalous descent over the CLLH is maintained by anomalous cold temperature advection in the middle troposphere. Therefore, now the key question is: How does the positive enthalpy advection by climatological wind in middle and upper troposphere form in the positive phase of Iseesaw?

Beyond the CLLH, in the positive phase of Iseesaw, late winter rainfall is significantly higher-than-normal over the equatorial central Pacific (Fig. 8a) associated with local SST warming (Fig. 4). Due to this enhanced heating, significant anomalous low-level flow convergence over this area is observed (Fig. 15a), with significant anomalous northerly winds to its northwest. At the same time, decreased rainfall is observed over the tropical western Pacific (Fig. 8a) associated with local negative SSTAs (Fig. 4). Consequently, an anomalous anticyclone is formed over the Philippine Sea (Fig. 15a) responds to this anomalous cooling, with significant anomalous northerly winds to its east. As a result of the combined impact of positive SSTAs over the equatorial central Pacific and negative SSTAs over the Philippine Sea, significant anomalous northerly winds prevail over the tropical western Pacific, and these anomalous northerly winds are also evident in different levels in the lower troposphere (e.g., 925-hPa level, 700-hPa level, figures not shown; 850-hPa level is shown in Fig. 16c). Via numerical simulation, Wang and Fan (1999) obtained a very similar low-level atmospheric circulation response to cold (warm) SSTAs centered at 10° N, 150° E (0° N, 140° W) in winter. Therefore, the anomalous pattern of low-level atmospheric circulation regressed onto the late winter Iseesaw is rather convincing.

Fig. 15
figure 15

Regression patterns of a wind anomalies (vectors; unit: m·s−1) and climatological specific humidity (shading; unit: g·kg−1) at 1000-hPa level, b vertically integrated WVF (vectors; unit: kg·m−1·s−1) and divergence anomalies (shading; unit: kg·m−2·s−1, multiplied by 105 for better presentation), d vertical velocity anomalies (unit: Pa·s−1, multiplied by 102 for better presentation) at 500-hPa level with respect to Iseesaw in late winter during 1982–2015. The black box indicates the CLLH (21–29° N, 97–108° E). Only the vectors exceeding the 95% confidence level are plotted in (a, b). The stippled regions indicate regression coefficients significant at the 95% confidence level in (bd). Altitude below 100 m are shaded with white in a

Fig. 16
figure 16

Regression patterns of a anomalous enthalpy advection by climatological wind (unit: W·m−1, multiplied by 103 for better presentation) at 350-hPa level, b temperature anomalies (shading; unit: K) and climatological wind (vectors; unit: m·s−1) at 350-hPa level, c wind anomalies (vectors; unit: m·s−1) and climatological temperature (shading; unit: K) at 850-hPa level, d temperature advection anomalies (unit: K·s−1, multiplied by 107 for better presentation) at 850-hPa level with respect to Iseesaw in late winter during 1982–2015. The black box indicates the CLLH (21–29° N, 97–108° E). The stippled regions indicate regression coefficients significant at the 95% confidence level in (a, b, d). Only the vectors exceeding the 95% confidence level are plotted in c. Altitude below 1500 m are shaded with white in (c, d)

The significant prevailing anomalous northerly winds over the tropical western Pacific in the positive phase of Iseesaw exert substantial impact on neighboring regions. In a climatological-mean sense, the spatial patterns of specific humidity in lower troposphere exhibit marked north–south gradient over the tropical western Pacific (Fig. 15a). Therefore, these significant prevailing anomalous northerly winds transport less water vapor to the Maritime Continent (MC) featured by evident anomalous northerly water vapor flux and divergence of water vapor flux (Fig. 15b), and hinder the formation of Precip over the MC. Besides the specific humidity, the climatological temperature in lower troposphere increases from north to south over the tropical western Pacific (Fig. 16c), hence the anomalous northerly winds also bring off-equatorial cold air into the MC, resulting in the decreased temperature in lower troposphere there (Fig. 16d). The above results indicate that the anomalous northerly winds advect low moist enthalpy (dry and cold) air into the MC, resulting in the obviously decreased moist enthalpy there (Fig. 15c). It is reported that anomalous vertical motion is strongly regulated by abnormal moist enthalpy advection in tropical regions (Wang and Zhang 2019; Yan et al. 2020; Lu et al. 2023). According to the decreased moist enthalpy, profound anomalous sinking motions dominate the MC (Fig. 15d), which also hinder the formation of Precip over this region. Due to the combined effect of reduced water vapor content and profound anomalous sinking motions induced by the significant prevailing anomalous northerly winds over the tropical western Pacific, apparent negative Precip anomalies are generated over the MC in the positive phase of Iseesaw.

The MC is regarded as the most important atmospheric heat source region of the world (Ramage 1968). By exciting atmospheric circulation anomalies via the Gill‐type response (Gill 1980), abnormal latent heat release associated with Pre anomalies over the MC significantly influences the global weather and climate (Jin and Hoskins 1995; Neale and Slingo 2003; He et al. 2017; Chen et al. 2022). In the positive phase of Iseesaw, corresponding to the pronounced negative Pre anomalies over the MC, a pair of significant anticyclonic anomalies in lower troposphere straddling the equator over the tropical north and south Indian Ocean are observed (Fig. 15a). Along with the significant anticyclonic anomalies over tropical north Indian Ocean, significant anticyclonic water vapor flux dominates the southern Bay of Bengal (SBOB, 3–15° N, 75–100° E, blue box in Fig. 8a), with apparent anomalous water vapor flux divergence (Fig. 15b). The anomalous divergence of water vapor flux over the SBOB may lead to the formation of reduced Precip and suppressed convection there (blue box in Fig. 8a). Besides, the significant surface anticyclonic anomalies may induce surface Ekman divergence and anomalous sinking motions in lower troposphere (blue solid-dotted line in Fig. 14d), which also partially contribute to the reduced Precip there. In turn, the suppressed convection over the SBOB results in strengthened sinking motion, further enhancing the negative Precip anomalies. Therefore, significant anomalous subsidence occupied the entire troposphere over the SBOB (blue solid-dotted line in Fig. 14d).

In response to the anomalous subsidence over the SBOB in the positive phase of Iseesaw, warm temperature anomalies are formed through adiabatic heating (red solid-dotted line in Fig. 14d). Figure 16b shows the anomalous temperature regressed onto the Iseesaw and climatological wind at 350-hPa level, at which maximum positive enthalpy advection exists over the CLLH (Fig. 14c). It is observed that besides the apparent positive temperature anomalies, climatological southwesterlies prevail over the Bay of Bengal. These climatological southwesterly winds advect anomalous warm air from the SBOB to the CLLH, leading to the positive enthalpy advection over the CLLH (Fig. 16a) and increasing of local MSE. As a result, an anomalous compensating upward motion is generated (blue solid-dotted line in Fig. 14a), which is beneficial to the increased Precip.

Besides the influence of Precip via abnormal vertical motion, the Iseesaw also exerts obvious impact on the CLLH SAT in late winter, especially in the eastern part. In the positive phase of Iseesaw, obvious positive sea level pressure (SLP) anomalies are seen over the MC and adjacent regions (Fig. 18c) associated with in situ negative Precip anomalies and significant descending anomalies. Climatologically, a distinct low (high) pressure system centers around the MC (Siberia) in the SLP field during boreal winter, forming a strong south-north pressure gradient between the tropical region and midlatitudes. This pressure gradient primitively drives near-surface northeasterly to the south of about 30° N in East Asia and contributes to the formation of East Asian winter monsoon (EAWM), which is one of the most active components of the climate system during boreal winter. Accordingly, the increased SLP in MC and adjacent regions induced by local suppressed convection has an essential impact on the EAWM. We have calculated the correlation coefficient between the EAWM intensity index and the Iseesaw, which is – 0.49 during 1982–2015 (significant at the 99% confidence level). The EAWM intensity index is defined following the method proposed by Wang and Chen (2014) based on a combination of normalized area-averaged SLP over Siberia, the North Pacific, and the MC. This result indicates that the EAWM intensity would be weaker (stronger) in the positive (negative) phase of Iseesaw. Being the western flank of anomalous anticyclone over Philippine Sea in lower troposphere, evident southerly wind anomalies are observed over the South China Sea in the lower troposphere (Figs. 15a, 16c). The climatological temperature in the lower troposphere exhibits clear north–south gradient there (Fig. 16c), hence the abnormal southwest winds would transport warm air into the eastern part of the CLLH. To confirm the impact of this abnormal southwest wind to the SAT variation, Fig. 16d shows horizontal temperature advection anomalies at 850 hPa regressed onto late winter Iseesaw. Obviously, anomalous warm temperature advections are located in southern China, including the eastern part of the CLLH. It is noted that the anomalous warm temperature advections over the south China correspond well with the positive SAT anomalies therein in the positive phase of late winter Iseesaw (Fig. 8a). The above results indicate that the formation of the warm SAT over the eastern part of CLLH in the positive phase of late winter Iseesaw is mainly attributed to the warm temperature advection induced by the Iseesaw-related anomalous anticyclone over Philippine Sea in lower troposphere.

4 Discussions

4.1 The physical mechanism of the later winter SST seesaw between the equatorial central Pacific and the Philippine Sea

In the previous sections, the impacts of later winter SST seesaw between the equatorial central Pacific and the Philippine Sea on early spring CLLH NDVI are analyzed. The CPI and the PSI is highly correlated, moreover, the result of Liang-Kleeman IF analysis indicates a causal relationship between the PSI and the CPI, suggesting that the CPI (PSI) exerts substantial (weak) influence on the PSI (CPI). In the following, the physical mechanism that the CPI exerts impact on the PSI is investigated.

Figure 17a, b depict the partial regression pattern of later winter SSTAs with respect to the Iseesaw independent of the CPI and PSI, respectively. After removing the PSI, apparent positive SSTAs remain unchanged over the equatorial central Pacific while the negative SSTAs disappear over the Philippine Sea. In contrast, significant SSTAs over the equatorial central Pacific and the Philippine Sea almost vanish after removing the CPI. Besides, in the positive phase of PC1, significant warm SSTAs over the equatorial central Pacific have already emerged in preceding autumn while no significant SSTAs appear over the Philippine Sea (Supplementary Fig. S3 c). These results verified the substantial impact of the CPI on the PSI.

Fig. 17
figure 17

Partial regression patterns of a SSTA (unit: °C) independent of the PSI, b SSTA (unit: °C) independent of the CPI with respect to Iseesaw in late winter during 1982–2015. Regression patterns of c surface latent heat flux anomalies (unit: W·m−2) and climatological wind (vectors; unit: m·s−1) at 1000-hPa level, d sea surface wind stress curl anomalies (unit: N·m−3, multiplied by 108 for better presentation) with respect to Iseesaw in late winter during 1982–2015. Only the vectors exceeding the 95% confidence level are plotted in (a, b). The stippled regions indicate regression coefficients significant at the 95% confidence level in. The white box in (d) indicate the mean NECB location 12–14° N and 127–130° E

Climatologically, the northeast trades prevail over the tropical western Pacific in lower troposphere (Fig. 17e). In the positive phase of Iseesaw, an anomalous anticyclone is situated over the Philippine Sea (Fig. 15a), with anomalous northerly winds to its east. Superimposed on the background northeasterly trade winds, the anomalous northerly winds increase surface wind speed and hence surface latent heat flux (Fig. 17c), favoring the formation of the negative SSTAs over the eastern part the Philippine Sea. However, over the western part the Philippine Sea, anomalous southerly winds to the west of this anomalous anticyclone are opposite to the climatic wind direction, resulting in a reduction of surface latent heat flux (Fig. 17c) and hindering formation of the negative SSTAs there. Hence, the negative SSTAs are attributed by other physical mechanisms.

The Pacific North Equatorial Current (NEC) is a westward current near the equatorial Pacific driven by climatological surface equatorial easterlies (Fig. 17e), and it plays an important role in the heat budget of the Pacific Warm Pool. In the positive phase of Iseesaw, significant anomalous surface equatorial westerlies are located over the equatorial western Pacific provoked by enhanced convection around the equatorial central Pacific (Fig. 15a). These surface equatorial westerlies are opposite to the climatic surface equatorial easterlies direction, leading to the weakened NEC. This may partially explain the negative SSTAs over the tropical western Pacific warm pool region in the positive phase of Iseesaw (Fig. 18a). After encountering the western boundary along the Philippine coast, the NEC bifurcates into the northward (southward) flowing Kuroshio (Mindanao) Current. The intensity of the Kuroshio Current and the Mindanao Current can be indicated by the latitude of the North Equatorial Current Bifurcation (NECB). When the NECB shifts poleward, the volume transport of the Kuroshio Current is reduced, and vice versa (Qiu and Lukas 1996; Qu et al. 1998; Kim et al. 2004). Qiu and Chen (2010) demonstrate that variation of the NECB can be well captured by wind‐induced local sea surface height (SSH) anomalies around the mean NECB location near the coast of the Philippines (12–14° N and 127–130° E, white box in Fig. 17d). Wang et al. (2020) indicate that the meridional shifts of the NECB and the SSH around the mean NECB location is negatively correlated.

Fig. 18
figure 18

Early spring NDVI anomalies over the CLLH obtained by regression on a the normalized INiño4, b the normalized IAO in late winter during 1982–2015. e The time series of normalized INiño4 (blue solid-dotted line) and normalized IAO (red solid-dotted line) in late winter. Stippling in (a) indicates where the NDVI anomalies are significant at the 95% confidence level

Figure 17d presents surface wind stress curl regressed onto the Iseesaw (Eq. (8)), it is shown that pronounced negative anomalous surface wind stress curl are observed over the Philippine Sea (Fig. 17d), associated with anomalous surface anticyclone there (Fig. 15a). This may lead to anomalous Ekman downwelling and decreasing of the SSH around the mean NECB location. As a result, the NECB shifts poleward and the transport of surface warm water to the Philippine Sea is weakened. The negative SSTAs over the western part of the Philippine Sea may result from the weakened Kuroshio Current in the positive phase of Iseesaw.

To sum up, the CPI exerts substantial impact on the PSI via wind-induced anomalous surface latent heat flux and wind-induced anomalous meridional shifting of the NECB. However, it is also noted that negative weak SSTAs still exist after removing the CPI (Fig. 17b). These results suggest that the CPI might be a primary factor of the PSI, but not the only factor.

4.2 The impacts of ENSO on early spring CLLH vegetation interannual variation

As shown in Fig. 4, it is noticeable that positive SSTAs are not only confined in the equatorial central Pacific but also extend to the equatorial eastern Pacific, suggesting the potential impact of late winter ENSO on the CLLH NDVI in early spring.

Statistically, the correlation coefficients between PC1 and late winter Niño indices INiño3, INiño3.4, INiño4, ICP_ENSO and IEP_ENSO are 0.22, 0.26, 0.33, 0.31and 0.17 during 1982–2015, respectively. None of the correlation coefficients reach the 95% confidence level while the correlation coefficients between PC1 and late winter INiño4 (blue solid-dotted line in Fig. 18c) as well as ICP_ENSO (figure not shown) exceed the 90% confidence level. Besides, although later winter ENSO and SST seesaw are highly correlated (the correlation coefficients between Iseesaw and INiño4 as well as ICP_ENSO are 0.93 and 0.87, respectively), the partial correlation coefficients between PC1 and Iseesaw independent of INiño4 as well as ICP_ENSO are 0.38 and 0.35, respectively, which still exceed the 95% confidence level. These results indicate that late winter SST seesaw between equatorial central Pacific and Philippine Sea, containing other signals besides the ENSO, can exert influence on the interannual variation of early spring CLLH NDVI. Meanwhile, the regression map of early spring NDVI anomalies over the CLLH onto the normalized later winter INiño4 is shown in Fig. 18a. Significant positive NDVI anomalies are limited in the eastern part of the CLLH in the positive phase of INiño4 whereas significant positive NDVI anomalies are located in most of the CLLH in the positive phase of Iseesaw. (Fig. 2b). The regression maps of early spring NDVI anomalies over the CLLH onto the normalized ICP_ENSO are similar to those onto the normalized INiño4 (figure not shown).

To investigate the mechanism of ENSO on CLLH NDVI, Fig. 19a, b depict the INiño4-regressed anomaly patterns of Pre and SAT in late winter, respectively. In the positive phase of INiño4, warm anomalies of SAT are located in the eastern part of the CLLH whereas there are weak Pre anomalies over the CLLH. The ICP_ENSO-regressed anomaly patterns of Pre and SAT exhibit similar results (figures not shown). It is noteworthy that the warm anomalies of SAT in the eastern part of the CLLH (Fig. 6c) are well coherent with the local increased NDVI in the positive phase of INiño4 (Fig. 4c). Hence, the INiño4-related warm anomalies of SAT in the eastern part of the CLLH are beneficial for local ST increasing and vegetation growing in early spring. On the contrary, the INiño4-related weak anomalies of Pre in the western part of the CLLH (Fig. 6d) are unfavorable for local SM increasing, leading to the insignificant changes of local NDVI in early spring (Fig. 4c). These results further confirm that later winter in situ SAT (Pre) plays a profound role in the early spring NDVI interannual variation in the eastern (western) part of the CLLH.

Fig. 19
figure 19

Regression patterns of a Precip anomalies (unit: mm·day−1) and b SAT anomalies (unit: °C) with respect to INiño4 in late winter during 1982–2015. c, d Same as in (a, b), but with respect to the IAO. The Precip and SAT data are derived from GPCP and CRU-TS v4.05, respectively. The black box indicates the CLLH (21–29° N, 97–108° E). The blue box indicates the southern Bay of Bengal (3–15° N, 75–100° E) as described in Sect. 3.3. The stippled regions indicate regression coefficients significant at the 95% confidence level

4.3 Contrast of impacts between the late winter SST seesaw and the ENSO on early spring CLLH vegetation interannual variation

Later winter ENSO and SST seesaw are highly correlated, the correlation coefficients between Iseesaw and INiño4 as well as ICP_ENSO are 0.93 and 0.87, respectively (exceeding the 99% confidence level). The partial correlation coefficients between PC1 and Iseesaw independent of INiño4 as well as ICP_ENSO are 0.38 and 0.35, respectively, which still exceed the 95% confidence level. These results indicate that late winter SST seesaw, containing signals of the ENSO and others besides the ENSO, can exert substantial influence on the interannual variation of early spring CLLH NDVI.

As shown in Figs. 3b and 18a, significant positive NDVI anomalies are located in the western and eastern part of the CLLH in the positive phase of the Iseesaw whereas significant positive NDVI anomalies are limited in the eastern part of the CLLH in the positive phase of INiño4. By examining the Precip and SAT anomalies (Figs. 8 and 19a, b), it is shown that increased SAT appears in the eastern part of the CLLH in the positive phase of the Iseesaw and the INiño4, however, enhanced Precip appears in the western part of the CLLH only in the positive phase of the Iseesaw. This can explain the relatively larger correlation coefficient between PC1 and Iseesaw (0.44) than that of INiño4 (0.33).

For the SSTs over the tropical Pacific, nearly equal increased SSTs over the equatorial central Pacific are observed in the positive phase of the Iseesaw and the INiño4 (Fig. 20a, b), however, the negative SSTAs over the Philippine Sea are weaker in the positive phase of the INiño4. The regressed area-average SSTA over the Philippine Sea with respect to the Iseesaw and the INiño4 is – 0.23 °C and – 0.16 °C, respectively. The area-average SLP over the Philippine Sea (10–20° N, 120–150° E, blue box in Fig. 20c, d) is used to measure the variability of the Philippine Sea anomalous anticyclone (Wang and Zhang 2002). The regressed area-average SLP over the Philippine Sea with respect to the Iseesaw and the INiño4 is 94.15 hPa and 87.63 hPa (Fig. 20c,d), respectively, indicating less influence of the INiño4 on the Philippine Sea anomalous anticyclone than the Iseesaw.

Fig. 20
figure 20

Regression patterns of a SSTA (unit: °C), c SLP (unit: hPa) with respect to INiño4 in late winter during 1982–2015. b, d Same as in (a, c), but with respect to the Iseesaw. The blue box in (c, d) indicates the region of Philippine Sea anomalous anticyclone (10–20° N, 120–150° E) as described in Sect. 4.3. The stippled regions indicate regression coefficients significant at the 95% confidence level

As analyzed in Sect. 3.3, the Philippine Sea anomalous anticyclone plays an essential role in the reduction of Precip over the SBOB. The regressed area-average Precip over the SBOB with respect to the Iseesaw and the INiño4 is – 0.51 mm·day−1 and – 0.45 mm·day−1 (Figs. 8a, 19a), respectively, and the INiño4 produces weaker subsidence over the SBOB due to the less decreased Precip. The enthalpy entering into the atmospheric column over the CLLH via the enthalpy advection in middle and upper troposphere is reduced, hence there is no significant change of the Precip over the CLLH. As a result, the NDVI anomalies in the western part of the CLLH are weak.

Therefore, Early spring CLLH NDVI is more closely related with the late winter SST seesaw of equatorial central Pacific–Philippine Sea than with the ENSO.

4.4 The impacts of AO on early spring CLLH vegetation interannual variation

Although there exists notably coherent variation between Iseesaw and PC1, with a correlation coefficient being 0.44, it is noted that the PC1 reaches its minimum value of – 2.71 in 2010 (Fig. 3d) for the period 1982–2015 while the Iseesaw is still positive (1.26) as shown in Fig. 18c, indicating an obvious out-of-phase variation between late winter SST seesaw and early spring CLLH NDVI of this year. Hence, other factors rather than the late winter SST seesaw of equatorial central Pacific–Philippine Sea played a vital role in the year 2010.

From autumn of 2009 to spring of 2010, the CLLH experienced a long-lasting drought, which was the most severe drought during the last 50 years and was regarded as a once-in-a-century drought (Yang et al. 2012). Due to this drought, vegetation productivity and terrestrial carbon cycling in the CLLH exhibited strong reductions (Zhang et al. 2012; Li et al. 2019a). Yang et al. (2012) demonstrated that this drought was featured by significant anomalous descent over the CLLH which was driven by the extreme negative phase of the AO. The normalized IAO in late winter (red solid-dotted line in Fig. 18c) reaches its minimum value of – 2.59 in 2010 for the period 1982–2015, and the correlation coefficients between PC1 and late winter IAO is 0.35 during 1982–2015, which is significant at the 95% confidence level. Fig. S9 presents the regressed geopotential height anomalies at 500-hPa level, it is clear that significant increased (decreased) geopotential height over the Arctic (mid-latitude) region. The above results indicate that late winter AO also substantially impacts the CLLH NDVI in early spring.

The regression map of early spring NDVI anomalies over the CLLH onto the normalized later winter IAO is shown in Fig. 18b. In the positive phase of IAO, significant positive (negative) NDVI anomalies are observed in the western (eastern) part of the CLLH. This dipole pattern of IAO-related NDVI anomalies is obviously different from that of PC1 and Iseesaw as well as INiño4. Figure 19c, d depict the IAO-regressed anomaly patterns of Pre and SAT in late winter, respectively. In the positive phase of IAO, increased (decreased) Pre (SAT) anomalies are located in most part of the CLLH. This abnormal wetter and colder land surface condition contributes to the increased SM in the western part of the CLLH and is favorable for local vegetation growing, corresponding to the positive NDVI anomalies there. However, this condition hinders the increasing of ST in the eastern part of the CLLH and is unfavorable for local vegetation growing, corresponding to the negative NDVI anomalies there. Furthermore, the correlation coefficients between late winter IAO and Iseesaw is – 0.16 during 1982–2015, indicating an independent relationship between late winter SST seesaw of equatorial central Pacific–Philippine Sea and the AO.

5 Conclusions

Early spring is a unique season of vegetation growth over the CLLH. Using a combination of reanalysis, land assimilation and satellite-based NDVI data sets, this study advances the knowledge of late winter tropical SSTAs’ impact on the early spring vegetation growth over the CLLH.

Firstly, the spatiotemporal characteristics of early spring CLLH NDVI are examined. It is indicated that the main mode of early spring CLLH NDVI, accounting for 33.6% of the total interannual variance, exhibits same-sign pattern. This homogenous pattern is found to be statistically linked to the late winter SST seesaw of equatorial central Pacific–Philippine Sea. In the positive phase of late winter Iseesaw, there are warm (cold) SSTAs in the equatorial central Pacific (Philippine Sea). The correlation coefficient between late winter Iseesaw and the principal component time series (PC1) corresponding to the main mode of early spring CLLH NDVI is 0.44 during 1982–2015 (exceeding the 99% confidence level). Accordingly, the NDVI in the early spring tends to increase in most of the CLLH in the positive phase of late winter Iseesaw.

The possible physical mechanisms responsible for this statistical linkage are analyzed. From a local perspective, early spring SM (ST) and late winter Pre (SAT) are identified as the key climate factors that affect early spring NDVI in the western (eastern) part of the CLLH via correlation and Liang-Kleeman IF analyses, while early spring incoming solar radiation has little impact. The positive phase of late winter Iseesaw results in stronger (higher)-than-normal Pre (SAT) in the western (eastern) part of the CLLH, the abnormal signals of late winter Pre (SAT) may be recorded by land surface and may persist to early spring owing to the SM (ST) memory. Accompanied with sufficient incoming solar radiation, wetter(warmer)-than-normal land surface condition in the western (eastern) part of the CLLH is favorable for in situ vegetation growth in early spring. For the late winter Iseesaw-related atmospheric circulation anomalies, enhanced upward motion is the critical factor that results in increased Precip over the CLLH. Through the MSE budget diagnosing, positive enthalpy advection by climatological wind in middle and upper troposphere is revealed to be the main reason leading to the enhanced upward motion in the positive phase of late winter Iseesaw. Beyond the CLLH, significant SST warming (cooling) and enhanced (reduced) Precip are located over the equatorial central Pacific (Philippine Sea) in the positive phase of late winter Iseesaw, as a result of these combined influences, an anomalous anticyclone is formed over the Philippine Sea in lower troposphere with anomalous northerly winds to its east. These anomalous northerly winds advect off-equatorial dry and cold air into the MC and result in negative Pre anomalies there. Consequently, an anomalous anticyclone over the SBOB in lower troposphere is provoked via a Gill-type Rossby wave atmospheric response. Featured by significant water vapor content reduction, this anomalous anticyclone gives rise to negative Pre anomalies over the SBOB, accompanied with evident local subsidence and significant warming through adiabatic heating. In middle and upper troposphere, climatological southwesterly winds advect anomalous warm air from the SBOB to the CLLH, leading to the positive enthalpy advection and increasing of local MSE. As a result, an anomalous compensating upward motion is generated, which is beneficial to the increased Precip. On the other hand, the positive SLP anomalies over the MC and adjacent regions as well as the southerly wind anomalies being the western flank of anomalous anticyclone over Philippine Sea over the South China Sea indicate a weaker-than-normal EAWM. The positive SAT anomalies in the eastern part of the CLLH would be a result of warm temperature advections induced by the southerly wind anomalies.

Finally, the physical mechanism of the later winter SST seesaw, the influence of ENSO and AO on the CLLH vegetation, the contrast of impacts between the late winter SST seesaw and the ENSO on CLLH vegetation is discussed. The outcomes indicate that the CPI exerts substantial impact on the PSI via wind-induced anomalous surface latent heat flux and wind-induced anomalous meridional shifting of the NECB; then both the ENSO and the AO affect the CLLH vegetation; what’s more, the late winter SST seesaw, containing signals of the ENSO and others besides the ENSO, can exert substantial influence on the interannual variation of early spring CLLH NDVI.

The results of this study contribute to a better understanding of the interannual variability of early spring vegetation over the CLLH. To our knowledge, the present results are the best that have ever been reported for the impact of remote large-scale climate factors on CLLH vegetation interannual variation. Of note, there are still some drawbacks in this study. First of all, previous studies demonstrated that human activities have significant impact on vegetation variations in the CLLH (Chen et al. 2021; Han and Song 2022; Fan et al. 2023). Although all the data, including the NDVI, were linearly detrended before analysis, the potential human influences on vegetation may still exist, which are not considered in this study and may lead to some limitations of the relevant results. Then, the influence of early spring SM and ST on CLLH NDVI is based on reanalysis data set and land surface assimilation product due to the lack of observed SM and ST data in the CLLH. Therefore, it is important to be aware of the uncertainties of SM and ST data derived from ERA5-Land and GLDAS-Noah. Finally, the possible physical mechanisms of this study were identified from the observation analysis. In the future, more in-depth studies based on model simulations will be required. Besides, we only adopted the NDVI as a proxy for vegetation in the current study. Other proxies, such as Leaf Area Index, solar-induced chlorophyll fluorescence, gross primary productivity, could be applied to reduce the uncertainty in vegetation variations in the future.