Observationally Based Analysis of Land–atmosphere Coupling

The temporal variance of soil moisture, vegetation and evapotranspiration over land has been recognized to be strongly connected to the temporal variance of precipitation. However, the feedbacks and couplings between these variables are still not well understood and quantified. Furthermore, soil moisture and vegetation processes are associated with a memory and therefore they may have important implications for predictability. In this study we apply a generalized linear method, specifically designed to assess the reciprocal forcing between connected fields, to the latest available observational data sets of global precipitation, evapotranspiration, vegetation and soil moisture content. For the first time a long global observational data set is used to investigate the spatial and temporal land variability and to characterize the relationships and feedbacks between land and precipitation. The variables considered show a significant coupling among each other. The analysis of the response of precipitation to soil moisture evidences a robust coupling between these two variables. In particular, the first two modes of variability in the precipitation forced by soil moisture appear to have a strong link with volcanic eruptions and El Niño–Southern Oscillation (ENSO) cycles, respectively, and these links are modulated by the effects of evapotranspiration and vegetation. It is suggested that vegetation state and soil moisture provide a biophysical memory of ENSO and major volcanic eruptions, revealed through delayed feedbacks on rainfall patterns. The third mode of variability reveals a trend very similar to the trend of the inter-hemispheric contrast in sea surface temperature (SST) and appears to be connected to greening/browning trends of vegetation over the last three decades.


Introduction
Soil moisture (SM) is an important variable of the climate system, playing an important role in the feedbacks between land surface and atmosphere.SM is important in determining climate variability at a wide range of temporal and spatial scales and controls hydrologic and energy cycles (Seneviratne et al., 2010;Dirmeyer, 2011).Soil moistureprecipitation feedbacks have been investigated at the global (Koster et al., 2004(Koster et al., , 2009) and the regional (Pal and Eltahir, 2003;Hohenegger et al., 2009) scale through numerical simulations.Recent observational studies have focused on local land-atmosphere coupling (Santanello et al., 2009).However, a comprehensive observational study at the global scale of the SM precipitation (PRE) coupling has never been performed.As shown by several modelling studies, it is over transition zones between wet and dry climates that a strong coupling between soil moisture and precipitation can be clearly identified, and it is over these regions that "soil moisture memory" can most probably contribute to subseasonal Published by Copernicus Publications on behalf of the European Geosciences Union.and longer climate predictions (Koster et al., 2004;Ferranti and Viterbo, 2006).The term "soil moisture memory" refers to the property of soil moisture to display persistent anomalies induced by climatic events like El Niño-Southern Oscillation (ENSO) or volcanic eruptions.Since slowly varying states of the land surface can be predicted weeks to months in advance, the response of the atmosphere to these landsurface anomalies can contribute to seasonal prediction.The large discrepancies among model results evidence the need of observational analysis of soil moisture-precipitation feedbacks (Seneviratne et al., 2010).The observational study by Alessandri and Navarra (2008) clearly identified a link between rainfall and land-surface-vegetation variability indicating an important delayed feedback of the land surface to the precipitation pattern.In this regard, a mechanism by which vegetation may provide delayed memory of El Niño and La Niña events is identified.
Predictability of climate at seasonal and longer timescales stems from the interaction of the atmosphere with slowly varying components of the climate system such as the ocean and the land surface (Shukla and Kinter, 2006).However, much of the model improvements so far have been obtained over ocean, where extensive availability of observations allowed model progresses and reliable application of assimilation techniques (Rosati et al., 1997;Alessandri et al., 2010Alessandri et al., , 2011)).In contrast, forecast performance over land is substantially weaker compared to the ocean (Wang et al., 2009;Alessandri et al., 2011).Since most of the applications of climate predictions would serve economic interests that are land-based, there is an urgent need to improve climate forecasts over land.Long-term improvements in understanding land-climate interactions and feedbacks over land must come from the enhancement of the description of the physical processes on the basis of dedicated process studies and observational databases.This can be suitably pursued firstly by analysing the newest available satellite-derived observational data sets that can lead to a better understanding and quantification of land-surface-atmosphere feedbacks.The better knowledge will then help us to conceive improved systems for the simulation of climate and for the improvement of its prediction at seasonal and possibly longer timescales.Here a global array of relevant up-to-date high-quality data sets is acquired, harmonized and analysed.The comprehensive data set is analysed to characterize the seasonal-mean interannual variability in land-surface variables and to improve understanding of the relationship and feedbacks between land and climate.The analysis method is based on the coupled manifold (CM) technique (Navarra and Tribbia, 2005), which was specifically designed to analyse covariation between fields considering both the local and remote forcing of one field to the other.The CM technique has proved to be successful for the analysis of different climate fields, like precipitation, vegetation characteristics, sea surface temperature (SST), and temperature over land (Alessandri and Navarra, 2008;Cherchi et al., 2007;Wang et al., 2011a).Recently, the CM tech-nique has been also applied to investigate the relationship between surface temperature and electricity demand in summer (De Felice et al., 2014).By taking advantage of the new global array of relevant up-to-date high-quality data sets, the present work substantially extends the analysis previously performed by Alessandri and Navarra (2008) and, for the first time, it includes SM and evapotranspiration (ET) feedbacks on PRE.
This paper is organized as follows.The observational data sets are described in Sect. 2. Section 3 describes the analysis method and gives a brief introduction of the CM technique.Section 4 presents the results.Summary and discussion of the main results of this study are given in Sect. 5.

The observational data sets
The data sets used for this study are all observationally based in order to make the analysis as independent as possible from global circulation model limitations and biases.Highquality up-to-date observational data sets of precipitation (PRE, from the Global Precipitation Climatology Project -GPCP), evapotranspiration (ET, from the University of Montana), soil moisture (SM, from the European Space Agency -ESA) and leaf area index (LAI, from Boston University) have been acquired and prepared.The selection of the data sets is based mainly on two criteria: (1) the period is covered as long as possible and (2) with global spatial coverage.The observed monthly PRE data set is described in Adler et al. (2003).ET values are satellite-based estimates from the Global Inventory Modeling and Mapping Studies (GIMMS) and MODIS (Zhang et al., 2010).The SM data set (Liu et al., 2011(Liu et al., , 2012) ) is the most complete record of this variable, based on active and passive microwave satellite sensors.The LAI data set (Zhu et al., 2013) is a long-term global data set resulting from the application of a neural network algorithm to the NDVI3g product from GIMMS satellite data.All land-surface data sets (SM, ET, LAI) are satellite products independent on the PRE data set, which is based on rain gauges.Despite that both ET and LAI products have been acquired by using the AVHRR sensor, the data sets have been produced by independent research groups which used completely different methodologies.The LAI product has been generated by applying a neural network algorithm on the NDVI satellite product, while the ET data set has been produced by using a modified Penman-Monteith approach including eddy covariance and meteorological data from the FLUXNET tower network.The time period, depending from the availability of the data sets, is 24 years  for ET and29 years (1982-2010) for the other variables.Original data sets come with various sampling frequencies, ranging from daily to monthly.See Table 1 for a summary of the characteristics of the retrieved data sets.
The data have been preprocessed and prepared for the subsequent analysis (Table 1).The preprocessing included space   (Zeng et al., 2013), the coverage reduces substantially during the winter season.The SM data set derives from blending passive and active microwave satellite retrievals.Figure 1b shows the percentage of SM missing data for each grid point.All grid points with a percentage of missing number larger than 30 % (white areas in Fig. 1b) have not been considered in the analysis.Over regions characterized by particularly dense vegetation and high canopies, both satellite products are unable to provide reliable estimates (Liu et al., 2012).Conversely, non-vegetated areas are associated with NaN values in the LAI data set.In order to evaluate the effect of major volcanic eruptions on land-atmosphere coupling, we used the stratospheric aerosol optical depth (AOD) at 550 nm, available from the NASA GISS data set (Sato et al., 1993).To evaluate the effect of ENSO, we compute the NINO3 index based on the Hadley Centre Global Sea Ice Coverage and Sea Surface Temperature (HadISST 1.1, 1870-present;Rayner et al., 2003) data set.

The analysis method
The CM technique (Navarra and Tribbia, 2005) seeks linear relations between two atmospheric fields Z and S (which in general are assumed to be rectangular matrices) of the kind The subscript "for" indicates the component of the field forced by the other variable (hereinafter forced manifold), while "free" indicates the free manifold.The free manifold contains the effects of nonlinearities.The linear operators A and B express the link between Z and S. A expresses the effect of S on Z, while B represents the effect of Z on S. In general, A and B are different.A and B are found by solving the Procrustes minimization problem: Following Navarra and Tribbia (2005), the technique is applied to the principal components of Z and S; therefore the coefficients of the linear operators A and B express the relations between the modes of the two variables.Canonical correlation analysis (CCA) scaling (data scaled by the covariance matrices) is applied to the principal components (PCs) of the variables before solving the Procrustes problem: where Ẑ and Ŝ are the CCA-scaled variables.Please refer to Navarra and Tribbia (2005) for further details of the CM technique.
As explained in Cherchi et al. (2007), after applying the CCA scaling, the elements of A and B are correlation coefficients and can be tested (with a significance test based on the Student t distribution) to reject the null hypothesis of being equal to zero.To improve the robustness of the analysis, each element of the A and B matrices has been verified to be different from zero at the 1 % significance level, following the method proposed by Cherchi et al. (2007).
The CM has two main advantages compared to other methods.The first one is that, when applied to a couple of climate fields (i.e.PRE and SM), CM is able to separate one field (i.e.PRE) into two components: the first component (forced) is the portion of PRE variability that is connected to the SM variability, whereas the second (free) is the part of PRE that is independent of SM.Therefore, the CM technique enables finding robust relations between fields in the presence of strong background noise.The second advantage is that the CM technique is able to detect both local and remote effects of the forcing variable.This is not possible with other methods such as SVD (singular value decomposition; Bretherton et al., 1992).
In the present analysis the CM technique has been applied to the seasonal-mean inter-annual anomalies.The climatological seasonal cycle has been removed and the data have been stratified using the seasons: JFM (January-February-March), AMJ (April-May-June), JAS (July-August-September) and OND (October-November-December).The JFM, AMJ, JAS, and OND stratification is used by Alessandri and Navarra (2008) in their CM study of vegetation and rainfall, which we will use to compare our results.The trends are kept for their relevance as possible indicators of climate change.
The LAI and SM data sets contain missing values, whose number and position significantly vary with time.The application of the CM algorithms requires that the number and position of the missing values is constant with time.Hence, if an NaN is present in a given grid point at any time, then it requires that grid point to be marked as NaN, thus losing a great amount of information.In order to retain as much information as possible from the data, we decided to replace the missing values with climatological values provided that their total number, considering a particular grid point, does not exceed a given threshold.We selected different thresholds for SM and LAI in order to obtain as similar as possible spatial coverage of the two variables.The chosen threshold is 10 for LAI and 30 % for SM.The results are robust with respect to a ± 10 % change in the threshold values.As shown in Fig. 1b, the areas more affected by the replacement of SM missing values (30 % of values replaced by climatology) are northeastern Europe, the east coast of central South America, eastern China and the Korean Peninsula.Since the replacement of missing values with climatology reduces time variability, the coupling in these regions may be underestimated as a consequence.We note that these gap-filled regions do not correspond to transition zones between wet and dry climates (Koster et al., 2000).Therefore, they are not expected to display a strong coupling between SM and PRE, nor significantly affect the main results of present study.
Since the main interest of the work is on the land surface, the ocean values are masked out from the PRE data set.A preliminary analysis (not shown) revealed that their inclusion resulted in a more difficult interpretation of the empirical orthogonal function (EOF) patterns (Bretherton et al., 1992), due to the interaction of phenomena on different space and timescales which are not connected to land variables.

Results
The CM technique has been applied to analyse the reciprocal forcing between PRE and the observed surface variables (SM, ET, LAI).The global-scale reciprocally forced temporal variances between PRE and the land-surface variables is reported in Table 2.In all, 19 % of the PRE variability is forced by SM.On the other hand, 17 % of the SM variance appears to be forced by PRE.In total, 18 % of the variability in PRE is forced by ET and 14 % of the variance of ET is forced by PRE.Considering the coupling between PRE and LAI, 17 % of the variance of PRE appears to be forced by LAI and 14 % of the variability in LAI is forced by PRE.All the variance ratios in Table 2 are significant at the 1 % level.The chance of coincidentally getting as high or higher ratios has been tested by means of a Monte Carlo bootstrap method (1000 repetitions).
Since SM is the most important land-surface parameter affecting seasonal to interannual variability/predictability of precipitation (Koster et al., 2000;Zhang et al., 2008), the coupling between SM and PRE will be analysed in detail in the following.

Reciprocal forcing between PRE and SM seasonal-mean anomalies
Figure 2 shows the ratio of the forced/total variance over land.The ratio of SM variance forced by PRE is in panel a, while panel b shows the ratio of PRE variability which is accounted for by the SM variability.For each grid point, the null hypothesis of coincidentally getting as high or higher variance ratios has been tested using a Monte Carlo bootstrap method (1000 repetitions).The regions where the ratio values are not significantly different from zero at the 1 % level are dotted.The observed SM variability appears to be intensely forced by PRE over the Sahel and central eastern Africa, southern Africa, the Middle East, the semi-arid region of central western Asia, the Indian subcontinent, Argentina, eastern Brazil and Australia.Note that, due to the limitations of the satellite estimates discussed in Sect.2, large areas in Russia and the Amazon Basin are not covered in the SM data set.
The larger observed effects on PRE due to SM inter-annual variability (Fig. 2b) occur in eastern Brazil, La Plata Basin, the Sahel, Asian boreal forests, the Middle East, Pakistan, Indonesia, and northern and eastern Australia.Most of these regions correspond to transition zones between dry and wet climates, where evaporation is highly sensitive to soil moisture (Koster et al., 2000).Here we refer to the transition regions between very dry and very humid environments, as individuated by Koster et al. (2000).By using the CM technique (described in Sect.3), the seasonal-mean PRE anomalies are separated into forced and free components, where forced and free refer to the influence of the SM variation.The variance explained by each mode of the PRE forced field is reported in Table 3.The EOF analysis shows that the first three components of the variability in the forced PRE field together account for 48 % of total variance.The first two PCs do not display trends, while the third PC is dominated by a clear trend, as will be discussed later.The first mode of variability in the forced PRE field explains 26 % of the total variance.The corresponding PC displays two significant peaks at years 1983 and 1992 (Fig. 3a).The PC is significantly correlated (maximum correlation coefficient equal to 0.56 at lag 0) with the stratospheric AOD.AOD peaks in correspondence to the two major eruptions of the period: 1983 (El Chichón) and 1992 (Pinatubo).The peaks in the AOD time series correspond to those of the forced PRE PC1, suggesting that this mode of variability is related to changes in the solar radiation at the ground, confirming that absorption and reflection of solar radiation by aerosol are particularly effective in reducing the hydrological cycle.The fast response of the precipitation anomalies to the radiation change induced by large tropical volcanic eruptions is in agreement with the results of the lag-correlation analysis by Gu and Adler (2011), who found 0 time lag between stratospheric aerosol signal and PRE.The lagged correlations of PC1 and AOD (Fig. 3b) show that significant (at the 5 % level) correlations endure up to about 2 years after the aerosol peak (i.e.behind the autocorrelation period of AOD itself; Fig. 3b dashed line).This result indicates that SM may provide a memory of the major volcanic eruptions for PRE.Table 4 shows the variance explained by each EOF mode of the whole original PRE field (that is, forced + free components).The link between PRE and volcanic eruption signal is also evident in the first mode of variability in the total rainfall field, as confirmed by the correlation of the cor-  4).
Figure 3c shows the spatial pattern of the first EOF of the PRE anomalies forced by the SM.A clear negative signal is present over areas characterized by a wet climate (Amazon Basin, India and Indonesia).In these regions the stratospheric aerosol emitted during the volcanic eruptions has the effect of reducing the intensity of the hydrological cycle (Alessan-  (Trenberth and Dai, 2007).In particular, according to Joseph and Zeng (2011) and Iles et al.
(2013), the negative signal over the monsoon regions may indicate a suppression of the monsoon linked to the effects of the aerosol released during major eruptions.Further, different from our results and other observational (Trenberth and Dai, 2007) and modelling (Joseph and Zeng, 2011)   the other hand, over transition zones (US Great Plains, Argentina, Middle East) the dimming effect may result in reduced evapotranspiration during the hot/dry season, which drives an increase in SM (Wild et al., 2009).During the following cool/wet season, the enhanced SM can induce a lagged increase in the portion of PRE forced by SM.That can explain the increased PRE over transition areas.On the other hand, the reduction of PRE over the southern Asia monsoon region and the enhancement of PRE over the semi-arid areas of central western Asia is consistent with the monsoondesert mechanism (Rodwell and Hoskins, 1996;Cherchi et al., 2014): the reduction of radiation caused by the stratospheric aerosol drives a reduction of convection over monsoon regions and a consequent reduction of PRE over southern Asia therefore abating Rossby wave-induced subsidence over the Middle East and the eastern Mediterranean (Cherchi et al., 2014).
The second PC of PRE forced by SM, explaining 14 % of total variance, is dominated by a large-scale oscillation (Fig. 4a).The corresponding principal component (full line) displays an high correlation coefficient of 0.60 with the NINO3 index (average of the SST in the tropical Pacific region 5 • S-5 • N, 210-270 • E; dashed line) at lag 2 (significant at the 1 % level), indicating that EOF2 represents the portion of the rainfall forced by SM that is related to the ENSO (Philander, 1989) variability.The second mode of forced PRE response due to SM variability appears to be lagged by one to several seasons with respect to the ENSO phase (Fig. 4b), with the strongest correlations with the NINO3 index two seasons after the maximum El Niño or La Niña intensity and significant correlations enduring until the lag 5 season (i.e.behind the autocorrelation period of ENSO itself; Fig. 5b dashed line).The results indicate that the effects related to ENSO in the SM may induce a delayed forcing on PRE.Therefore, SM appears to provide a biophysical memory of ENSO on the global precipitation pattern.The signal of ENSO can also be evidenced in the second mode of variability in the total rainfall field as indicated by the correlation of the corresponding PC (explaining 5 % of total PRE variance) with NINO3 (Table 4).Again, the lag at which maximum correlation is attained is the same (lag 2) as in the forced field, but the correlation coefficient is 0.60 for the forced field and 0.43 for the total PRE field.
The spatial pattern of the second EOF of the PRE anomalies forced by SM (Fig. 4c) displays the signature of the tripole pattern over South America typical of ENSO teleconnections (Ropelewski and Halpert, 1989).Similarly, negative PRE anomalies are shown over Brazil, southern Africa, northern India and Indochina, displaying the land-surface feedback to the reduced rainfall related to the positive phase of ENSO there (Trenberth et al., 1998).On the other hand, positive precipitation anomalies characterize the west and east coasts of North America, central America, the dry and semi-arid region of northern Venezuela, La Plata Basin, the Horn of Africa, the Sahel, Europe, central and eastern Asia, southern India and the east coast of Australia.Most previ- ous research showed reduced precipitation over India during ENSO years (Ropelewski and Halpert, 1989;Trenberth et al., 1998).The positive anomalies of PRE forced by SM over southern India related to the positive phase of ENSO evidence an interesting negative feedback of the land surface on the effect of ENSO on the rainfall over India.
The third PC of the PRE forced by the SM, explaining 8 % of forced variance, displays a trend (Fig. 5a) corresponding to a clear signal of increasing precipitation over the Sahel, southeastern Europe, central Asia, northeastern Asia, the Great Plains of North America, and Nordeste (Brazil) and the northern part of South America (Fig. 5b).The trend of increasing precipitation is particularly strong over the Sahel, where, according to Hagos and Cook (2008), it can be related to a warming of the northern tropical Atlantic Ocean which, through a modification of the associated cyclonic circulation, enhances moisture transport over the region.In contrast, a decrease in precipitation is evident over most of the Southern Hemisphere (SH), northwestern Russia, eastern Russia, northern India, China and the western US, showing a northsouth polarity of the precipitation trend.The above trend pattern strongly resembles the trend pattern of global rainfall annual mean anomalies described by Munemoto and Tachibana (2012, hereinafter MT12).The authors associated this northsouth polarity with a relatively larger warming of the Northern Hemisphere (NH) compared to the SH that characterized the last three decades, starting from the early 1980s.MT12 found that the trend of the SST corresponds to an increase in the specific humidity in the NH with respect to SH that enhances (reduces) precipitation in the NH (SH).Although the focus of MT12 is on the Sahel region, the authors defined a global index, the north-south SST (NS-SST) polarity index, which successfully captures the global signal of the precipitation trend.The NS-SST index is defined as the areaaveraged NH SST annual mean anomalies minus the SH SST anomalies.The NS-SST index (computed from HadISST), normalized by its standard deviation, and its trend are plotted in Fig. 5a.Note that here the NS-SST index is computed from the seasonal mean anomalies instead of the annual mean anomalies used in MT12; nonetheless, the trend is not affected.

Mediation effects of ET and LAI on the coupling between PRE and SM
To investigate how the coupling between rainfall and soil moisture is mediated by evapotranspiration and vegetation we further applied the CM technique between the component of PRE forced by SM and ET (LAI), obtaining the component of PRE forced by SM which is also forced by ET (LAI).As summarized in Table 5, 20 % of the inter-annual variability in the PRE anomalies forced by the SM is estimated to be globally forced by the ET variation.It is important to note here that 19 % × 20 % = 3.8 % represents only the ET forc-  ing on PRE mediated by SM and not the whole ET forcing on PRE, which is actually 18 % (Table 2).At the same time, 23 % of the variance of PRE forced by SM is evaluated to be also forced by LAI; therefore, the LAI forcing on PRE mediated by SM corresponds to 17 % × 23 % = 3.9 %.
Figure 6a shows the ratio of the variance of PRE forced by the SM which is also forced by the ET with respect to the total forced rainfall variance.Figure 6b shows the same plot but for the LAI.The "hotspots" in Fig. 6a are similar to those found in Fig. 2b over the Sahel, the Horn of Africa, eastern Europe, Asian boreal forests, central Asia, the west coast of the US, eastern Brazil and La Plata Basin.This indicates that in all these regions the link between PRE and SM is at least in part mediated by ET.Not surprisingly, the same regions also display a link with vegetation (Fig. 6b).Furthermore, vegetation appears to significantly affect rainfall variability over the semi-arid regions that are not dependent on ET such as central western Asia, southeastern Africa, southeastern Asia and western Australia, suggesting that in these regions the SM forcing on PRE is mediated by vegetation state (e.g.stress of vegetation will affect PRE there).
To analyse how the response of PRE forced by SM to climate events and the rainfall trend are mediated by ET (LAI), we applied the CM technique between each of the physical fields corresponding to the first three modes of variability in PRE forced by the SM and ET (LAI).Here we take the physical fields corresponding to the first three modes of variability in PRE forced by SM and further decompose them to extract the parts of each mode that are forced individually by ET and LAI.This analysis allows for solving how ET and LAI contribute to each component of PRE forced by SM, which has been identified to be linked to external climate forcing (volcanic eruptions, ENSO and trend).Overall, considering the global land, 21 % of the variance displayed by the first mode (linked to volcanic eruptions) of PRE forced by the SM is forced by the ET and 27 % by LAI (Table 6).As for the second mode (connected to ENSO), 38 % of the variance is forced by ET and 36 % by LAI.Concerning the third mode (displaying a trend), 31 % of the variance is forced by ET and 29 % is forced by LAI.Rainfall variability forced by the ET and LAI decomposed through EOF analysis is reported in Table 7.Interestingly, the third PC of the PRE forced by the ET (explaining 7 % of the forced variance) is correlated with AOD, with a maximum correlation coefficient of 0.41 at lag 6. Analogously, the second PC of the PRE anomalies forced by the LAI (explaining 10 % of the forced variance) is correlated with AOD, with a maximum correlation coefficient of 0.41 at lag 3, suggesting that both ET and vegetation contribute to providing memory of volcanic eruptions, modulating at longer scales the effect of the SM forcing on PRE.The first PC of PRE forced by ET (explaining 30 % of the forced variance) is found to be significantly correlated with the NINO3 index, with a correlation coefficient of 0.52 at lag 0. The first PC of PRE forced by LAI (explaining 27 % of the forced variance) also has a maximum correlation coefficient of 0.67 at lag 0 with the NINO3 index, indicating that vegetation acts as the mediator at longer scales of the signal between SM and PRE.This result is consistent with the relationship found by Alessandri and Navarra (2008) between precipitation forced by vegetation (NDVI) and ENSO and with the delayed vegetation response to ENSO signal found by Zeng et al. (2005).All the above correlation coefficients passed a significance test at the 1 % level.To determine the regions where the mediating effects of ET and LAI have the larger influence on the coupling with respect to the stratospheric volcanic eruptions, the first mode of variability in PRE forced by the SM has been correlated with the total components of PRE forced by the ET and LAI.The correlation coefficients are shown in Fig. 7a for PRE forced by the ET and Fig. 7b for PRE forced by the LAI.Only the regions where correlations passed a significance test at the 5 % level are shaded.Black upward (white downward) triangles denote areas with positive (negative) values of the first EOF of the PRE anomalies forced by the SM (Fig. 3c).The correlations are positive almost everywhere (i.e. the effects of both ET and LAI tend to amplify the response of rainfall to large volcanic eruptions) and the patterns are very similar for ET and LAI, indicating that the feedback of ET may Table 7. Rainfall variability forced by the ET and LAI decomposed through EOF analysis.Each line displays the EOF explained variance (column 2) and the corresponding PC correlation with relevant climatic indices (column 3).Here the maximum PC correlation is reported considering lagged correlations in the range −16 to +16.Only the correlation coefficients significant at 1 % level are reported.Plata Basin, western central Asia, the Horn of Africa, southern Africa, the Asian monsoon region, Indonesia and Australia.Over these regions, evapotranspiration and vegetation activity are radiation limited (Seneviratne et al., 2010).Nevertheless, while over some regions (southern part of North America, La Plata Basin, the Middle East, western central Asia and the Horn of Africa) ET and LAI contribute to an increase in rainfall, over other regions (northern South America, southern Africa, Indian monsoon region, Australia) they contribute to rainfall reduction.As discussed in Sect.4.1, over most of the SH (apart from La Plata Basin and the Horn of Africa) and the Asian monsoon region there is a reduction of precipitation that can be associated with the dimming effect and the consequent reduction of the hydrological cycle.In humid regions the rainfall reduction can stress vegetation and may reduce its growth, with effects lasting up to 1 year (Wang et al., 2011b).On the other hand, over most of the arid and semi-arid regions (Middle East, western central Asia), the reduced evapotranspiration during past seasons induced by the dimming effect may increase SM and therefore attenuate the stress on vegetation.This, in turn, has a positive effect on precipitation.
The point-by-point correlation coefficient between the second mode of variability (related to ENSO) in PRE forced by the SM and the total field of PRE forced by the ET is shown in Fig. 8a.The correlation between the second mode of PRE forced by the SM and the total field of PRE forced by the LAI is shown in Fig. 8b.The sign of the feedback between PRE and SM is indicated by the second EOF of PRE forced by the SM overlaid on the plot.Large positive correlations up to 0.6 are found globally over most of the land areas.ET has a positive feedback on the increase in precipitation over the west coast of the US, the dry and semi-arid region of northern Venezuela, La Plata Basin, the Sahel, northern Europe, India, central and eastern Asia, and the southeastern coast of Australia.Still, a positive feedback is present over Brazil, southern Africa and Indochina, but in this case ET leads to further reduction of PRE.A negative feedback of ET is seen over Mexico.In this region the positive ENSO phase induces wet and cool conditions (Trenberth et al., 1998) associated with an increase in PRE forced by SM that is contrasted by a reduction of ET.As for vegetation, it contributes to rainfall enhancement over east and west coasts of the US, La Plata Basin, northern Europe, the Horn of Africa, the semi-arid region of western central Asia and eastern Asia.Conversely, vegetation mediates precipitation reduction over Brazil, southern Africa and Indochina.
Figure 9 shows the point-by-point correlation coefficient between the third mode of variability in PRE forced by the SM (displaying a linear trend; see Fig. 5) and the total fields of PRE forced by the ET (Fig. 9a) and PRE forced by the LAI (Fig. 9b) with the third EOF of PRE forced by the SM overlaid on it.The feedback of ET on this mode of variability in PRE is not significant over most of the NH.A positive effect of ET is seen over the semi-arid regions of the SH, but while ET mediates an increase in rainfall over the Sahel, it leads to further reduction of PRE (Fig. 9a) over Bolivia and Australia.On the other hand, ET has a negative feedback over the humid region of Tanzania where it contrasts the reduction of PRE.The pattern of the feedback of LAI (Fig. 9b) is very different from that of ET.Overall, the vegetation has a positive feedback on the rainfall anomaly pattern forced by the SM.In particular, large correlations up to 0.6 are seen over the Sahel, the east coast of the US, western South America, eastern Europe, tropical southern Africa, western central Asia, Asian boreal forests, central and eastern Asia, the Indian monsoon region and eastern Australia.The strong signal over the Sahel is in agreement with Zeng et al. (1999)   trend similar to that of the NS-SST index, analogous to the third PC of PRE forced by SM, while no trends are found in the first five PCs of PRE forced by ET.

Conclusions
A global array of relevant up-to-date high-quality data sets (soil moisture, evapotranspiration, leaf area index and precipitation) is acquired, harmonized and analysed.For the first time a long comprehensive global observational data set is used to characterize the land variability as a function of the space and timescales and to improve understanding of the relationships and feedbacks between land and climate.By applying the CM technique on the seasonal-mean inter-annual anomalies, the relationship and the coupling between the acquired surface variables is assessed considering all the seasons.
The analysis shows a considerable degree of reciprocal forcing and coupling in the land-surface variables considered.The reciprocal forcing with precipitation is particularly strong for the soil moisture, with 19 % of the inter-annual variability in the precipitation over continental areas that are forced by the SM variation.Conversely, 17 % of the SM variance is forced by PRE.
The PRE forced by the SM is dominated by a prominent decadal-scale drying, initiated by the perturbation of the abrupt Mt Pinatubo eruption.In 1991, PC1 of the dominant forced mode of PRE shows an abrupt decrease and the negative anomaly continues increasing in the subsequent years until 1994.It is only after 1995 that the rainfall starts to slowly recover towards the pre-eruption levels.In 1997, the signal sums up with that of ENSO.It appears that the persistence of the negative SM anomalies leads to increasing stress conditions for the vegetation, thus leading to a larger ET response at longer time lags after the perturbing event.Our interpretation is that the persistence of the negative SM anomalies provides the memory of the initial perturbing event and our analysis indicate that, through this mechanism, the effect of Mt Pinatubo eruption can last for several years and its memory appears to extend and sum to the following 1997-1998 El Niño event.The second PC of the PRE forced by the SM displays a large-scale oscillation correlated to ENSO variability, with significant correlations enduring behind the autocorrelation period of ENSO itself and up to more than 1-year lag.This indicates that ENSO effects on SM induce a delayed forcing on PRE.The third PC of the PRE forced by the SM is dominated by a trend, positive over most of the NH and negative over most of the SH.This trend appears to be related to the inter-hemispheric SST contrast which corresponds to an increase in the specific humidity in the NH with respect to the SH that enhances (reduces) precipitation and SM in the NH (SH).
The combined analysis of the PRE modes related to the external climate forcings (volcanic eruptions, ENSO, SST trend) and the rainfall forced by ET and LAI evidences the role of ET and LAI as the mediators between SM forcing and rainfall.In particular, it appears that both ET and LAI tend to provide a positive feedback on PRE over most of the regions, contributing to further enhancing or reducing rainfall depending on the regions of the globe, with large differences between wet, transition and semi-arid climates.Nevertheless, the response to ENSO is characterized by a negative feedback of ET over regions where the positive ENSO phase induces wet and cool conditions (i.e.Mexico).
It is important to note that the coupling with SM revealed by the present analysis has to be considered an underestimation of the real coupling, due to the incomplete cover of the SM data set.Nevertheless, the present investigation identifies the regions characterized by a strong coupling and suggests most possible mechanisms linking the considered variables.
Since SM has been recognized as the most important landsurface parameter affecting seasonal to interannual variability in precipitation (Koster et al., 2000;Zhang et al., 2008), the present paper focused on the coupling between SM and PRE.Detailed analysis of the reciprocal forcing between ET and LAI, LAI and SM, and ET and SM will be the subject a future paper that will further address the specific coupling among land-surface variables.

Data availability
Evapotranspiration data set available from the Numerical Terradynamic Simulation Group (NTSG) of the University of Montana (http://www.ntsg.umt.edu/project/et). Leaf

Figure 1 .
Figure 1.(a) Global mean missing values in the time series (in %): LAI (full) and SM (dashed).(b) Map of the percentage of SM missing data for each grid point.All grid points with a percentage of missing number larger than 30 % (white areas in b) have not been considered in the analysis.

Figure 2 .
Figure 2. Ratio of the forced variance to the total variance.(a) The fraction of SM variance forced by PRE.(b) The fraction of PRE variance forced by the SM.Dots are placed over areas covered by the forced variable data set but where variance ratio values did not pass a significance test at the 1 % level.

Table 4 .
Total rainfall variability decomposed through EOF analysis.Each line displays the EOF explained variance (column 2) and the corresponding PC correlation with relevant climatic indices (column 3).Here the maximum PC correlation is reported considering lagged correlations in the range −16 to +16.Only the correlation coefficients significant at the 1 % level are reported.Variance explained Correlation with climate indices PC ., 2012) with a consequent reduction of SM, PRE and continental discharge

Figure 3 .
Figure 3. (a) First normalized PC of the PRE anomalies forced by the SM (full line and filled circles), after cutoff low-pass filtering at 2 year −1 frequency.Dashed line (and cross marks) stands for the normalized stratospheric aerosol optical depth (AOD).Lines stand for 5-year exponential moving average, while marks represent each single season.(b) Lagged correlations between AOD and PC1 of the forced PRE.The dashed curve is the autocorrelation function of the AOD.Marks indicate significance at the 5 % level.(c) First EOF of the forced PRE.Arbitrary units.

Figure 4 .
Figure 4. (a) Second normalized PC of the PRE anomalies forced by the SM (full line and filled circles).Dashed line (and cross marks) stands for the normalized NINO3 index.Lines stand for three-season running means, while marks represent each single season.(b) Lagged correlations between NINO3 index and PC1 of forced PRE.The dashed curve is the autocorrelation function of the NINO3 index.Marks indicate significance at the 5 % level.(c) Second EOF of the forced PRE.Arbitrary units.

Figure 5 .
Figure 5. (a) Third normalized PC of the PRE anomalies forced by the SM (full line and filled circles).Dashed line (and cross marks) stands for the normalized NS-SST index.Lines stand for three-season running means, while marks represent each single season.Coloured lines represent the trends (red for the PC, blue for the NS-SST index).(b) Third EOF of the forced PRE.Arbitrary units.

Figure 6 .
Figure 6.Ratio of the forced variance to the total variance.(a) The fraction of PRE variance forced by the SM which is also forced by the ET.(b) The fraction of PRE variance forced by the SM which is also forced by the LAI.Dots are placed over areas where variance ratio values did not pass a significance test at the 1 % level.

Table 6 .
Ratios of the global-scale forced variance over the total variance resulting from the application of the CM technique between the first three modes of PRE forced by SM and the total fields of ET and LAI.ET LAI PRE forced by SM mode 1 0.21 0.27 PRE forced by SM mode 2 0.38 0.36 PRE forced by SM mode 3 0.31 0.29

Figure 7 .
Figure 7. Point-by-point correlation of the first mode of variability in PRE forced by SM with (a) the total fields of PRE forced by ET and (b) the PRE forced by LAI.Data have been filtered using a cutoff low-pass filter at 1 year −1 frequency.Only areas where correlations passed a significance test at the 5 % level are shown.Black upward (white downward) triangles denote areas with positive > 0.01 (negative < −0.01) values of the first EOF of the PRE anomalies forced by the SM (Fig. 3c).

Figure 8 .
Figure 8. Point-by-point correlation of the second mode of variability in PRE forced by SM with (a) the total fields of PRE forced by ET and (b) the PRE forced by LAI.Data have been filtered using a cutoff low-pass filter at 1 year −1 frequency.Only areas where correlations passed a significance test at the 5 % level are shown.Black upward (white downward) triangles denote areas with positive > 0.01 (negative < −0.01) values of the second EOF of the PRE anomalies forced by the SM (Fig. 4c).

Figure 9 .
Figure 9. Point-by-point correlation of the third mode of variability in PRE forced by SM on the total field of (a) PRE forced by ET and (b) PRE forced by LAI.(c) Magnitude of the LAI change over 1982-2010, quantified using a linear model under the assumption of monotonic change.Data have been filtered using a cutoff low-pass filter at 1 year −1 frequency.Only areas where correlations (a-b) and trend (c) passed a significance test at the 5 % level are shown.Black upward (white downward) triangles denote areas with positive > 0.01 (negative < −0.01) values of the third EOF of the PRE anomalies forced by the SM (Fig. 5b).

Table 2 .
Ratios of the global-scale forced and free variance with respect to the total variance resulting from the application of the CM technique between PRE and SM, ET and LAI.

Table 5 .
Ratios of the global-scale forced and free variance with respect to the total variance resulting from the application of the CM technique between PRE forced by SM and ET as well as LAI.