Interannual Variation of Transpiration and Its Modeling of a Larch Plantation in Semiarid Northwest China

: Quantifying the variation of forest transpiration (T) is important not only for understanding the water and energy budget of forest ecosystems but also for the prediction, evaluation, and management of hydrological effects as well as many other ecosystem services of forests under the changes of climate, vegetation, and anthropological impacts. The accurate prediction of T, a key component of water used by forests, requires mechanism-based models describing the T response to environmental and canopy conditions. The daily T of a larch ( Larix principis-rupprechtti ) plantation was measured through monitoring the sap flow in the growing season (from May to September) of a dry year (2010), a normal year (2012), and a wet year (2014) at a shady slope in the semi-arid area of Liupan Mountains in northwest China. Meanwhile, the meteorological conditions, soil moisture, and forest canopy leaf area index (LAI) were monitored. simple applicable T model, three the evapotranspiration evapotranspiration (PET), soil water supply ability relative extractable transpiration The T model was established as a continuous multiplication of the T response equations to individual factors, which were determined using the upper boundary lines of measured data. The effect of each factor on the T in a dry year (2010) or normal year (2012) was assessed by comparing the measured T in the baseline of the wet year (2014) and the model predicted T, which was calculated through inputting the actual data of the factor (i.e., PET) to be assessed in the dry or normal year and the measured data of other two factors (i.e., REW, LAI) in the baseline of the wet year. The results showed that the mean daily T was 0.92, 1.05, and 1.02 mm; and the maximum daily T was 1.78, 1.92, and 1.89 mm in 2010, 2012, and 2014, respectively. The T response follows a parabolic equation to PET, but a saturated exponential equation to REW and LAI. The T model parameters were calibrated using measured data in 2010 and 2012 (R 2 = 0.89, Nash coefficient = 0.88) and validated using measured data in 2014 satisfactorily (R 2 = 0.89, Nash coefficient = 0.79). It showed a T limitation in the dry year 2010 for all factors (18.5 mm by PET, 11.5 mm by REW, and 17.8 mm by LAI); while a promotion for PET (1.4 mm) and a limitation for REW (4.2 mm) and LAI (14.3 mm) in the normal year 2012. The daily T model established in this study can be helpful to assess the individual factor impact on T and improve the daily T prediction under changing


Introduction
To increase the various ecosystem services supplied by forests, such as erosion control, timber production, carbon sequestration, and water regulation, afforestation has been encouraged worldwide, including in dryland regions. However, afforestation has increased the consumption of water, which is inherently scarce in the vast dryland regions [1][2][3][4][5]. The water used by forests in some sites and in some periods can even exceed annual precipitation [3,6]. Many studies have confirmed soil moisture decrease and even dried soil layers [6] and remarkable water yield reduction due to large-scale afforestation in the watersheds on the Loess Plateau in Northwest China [7][8][9][10][11][12]. Thus, increasing attention has been given to forest water consumption [13] and balanced forest-water management.
Quantifying plant transpiration (T) variation is important for not only the estimation of a water budget but also the elucidation of the role of T in an energy budget, the effects of hydrological flows on vegetation production, and ecosystem responses to climate change [14][15][16][17]. Tree T, which utilizes a large proportion of forest evapotranspiration [3], is a key factor affecting the stability and services of forests, and the regional water balance, or even water supply security [13,[18][19][20]. Therefore, it is very necessary to fully understand and accurately predict T dynamics [17,21,22], especially in vast dryland regions with variable climates and severe water shortages. To simplify the study of forest T variation, the numerous influencing parameters can be grouped into three factors: the weather-related evaporative driver, the soil moisture-related extractable soil water, and the vegetation-related canopy T capacity. Establishing a daily T model that can couple all three factors and be applied easily is the basic requirement for forest T predictions under changing environmental and canopy conditions. However, this kind of study is still very rare, although there have been many studies on the effects of individual factors.
Forest T is firstly determined and usually limited by the atmospheric potential evapotranspiration (PET), which is an integrated effect of many weather parameters, such as temperature, radiation, and vapor pressure deficit (VPD) [23,24]. PET could explain over 60% of the T variation in a black locust (Robinia pseudoacacia) plantation on the Loess Plateau of China [25]. The study showed a mean ratio of 58% (range 50-64%) of daily T to daily PET [26]. Many PET-T relation studies of various tree species in different regions [26][27][28][29] have been conducted usually using a logistic, linear, or quadratic equation to describe the T response to PET but without considering other influencing factors.
Quantifying the relationship between T and relative extractable soil water content (REW) has been an interesting topic in forest hydrology. The use of REW, rather than the absolute soil moisture, allows a direct comparison among the studies with different soil physical properties. Some of these studies were conducted on forests for a number of tree species in various counties [26][27][28][30][31][32][33][34]. These studies showed a common T variation trend with rising soil moisture, which initially increases quickly and almost linearly but thereafter gradually stops increasing and approaches a maximum when the REW reaches a threshold usually in the range of 0.2-0.5. Nonlinear (asymptotic) models, such as quadratic, logistic, logarithmic, negative exponent, and piecewise equations, are often used to depict the T response to rising REW.
Soil moisture is the key factor limiting tree growth and forest T in dryland regions. Li et al. [35] found that the reduction in forest T by REW accounted for 16% of the total potential T in a growing season for the same larch plantation as that in this study in a semiarid climate. However, this result was based on a simplified T-REW relation, leaving out the effect of evaporative demand. In fact, evaporative demand most often controls the T variation in the field, and the phenomenon that T is limited by soil moisture appears only during periods of a high soil water deficit [36,37]. For example, Ungar et al. observed that soil moisture acted as the primary driver for forest T only when the soil moisture was below 0.15 cm 3 ·cm −3 [28]. The forest T and its response to REW were often controlled or affected by varying weather conditions [38]. A study in a Cyclobalanopsis glauca forest in South China showed that the T responses to soil moisture were more intense under conditions of higher VPD [39]. To obtain an ideal T-REW relation, some studies have been conducted with seedlings in chambers under artificially manipulated weather conditions [25] or in a field where the data on cloudy, overcast, and rainy days or days with extremely large or small VPD are excluded [32]. It is possible that such an analysis is unable to properly represent the actual T behavior in the field.
In addition to PET and REW, forest canopy leaf area index (LAI) is also a key factor affecting T, as has been confirmed by many studies, such as those of Vertessy et al., Zimmermann et al.,. Based on a review of studies on 20 forests and 12 tree species under different climate types, Granier et al. derived a generic response of canopy conductance to a rising LAI within the range of 0-11 [43]. The T first increases with rising LAI but at a declining rate, and then, T stops increasing when the LAI exceeds a threshold. However, the T response to LAI is expressed differently, often as a linear, power, or single asymptote equation or piecewise equations. Even some individual studies have concluded that T is not correlated with LAI [41]. The main reason to explain the diverse T-LAI relations could possibly be the difference in the LAI range among the above studies, such as 2.9-4.8 [40], 1-11 [43], 0.6-1.6 [41], and 1-6 [42]. Only the LAI varies within a big enough range (e.g., >6), the whole T response to LAI including the LAI threshold, from that T stops increasing with LAI, can be observed.
Forest T is codetermined by weather conditions, soil moisture, and vegetation characteristics, but direct studies on the simultaneous response of T to these factors are rare. The vast majority of the existing studies are focused on the response of forest canopy conductance to environmental factors. A multiplication approach has often been used to couple the effects of individual factors (e.g., PET, REW, and LAI) [44,45]. Granier et al. used such an approach to describe the T response to multiple meteorological factors (radiation, VPD, and temperature), soil moisture, and LAI, based on the measured data of forests of 21 broadleaved or coniferous tree species under different climates (temperate, mountain, tropical, and boreal) [43]. A similar study was conducted on the forest canopy conductance of oak (Quercus robur) in the Netherlands [46,47], and maritime pine (Pinus pinaster Ait.) in the Landes de Gascogne Forest of southwestern France [23,48]. However, such studies on forest T are still rare.
In this study, the daily T variation in a semiarid larch plantation was continuously measured in the growing season of 2010-2014. Due to the power supply failure, some T data were seriously missed; especially due to the instrument fault of a plant canopy analyzer (LAI-2000 type, LI-COR, Lincoln, Nebraska, USA), the LAI could not be measured in 2011 and 2013. Therefore, only the data from the dry year of 2010, normal year of 2012, and wet year of 2014 were used. The aims of this paper were to (1) establish a simple but solid model for predicting the daily T of forests under field conditions, and (2) clarify the importance of LAI, REW, and PET to the T variation at the study site. This study will promote a better understanding of the daily T response to varying environmental and canopy conditions and supply a basis for the integrated forest-water management.

Site Description
This study was performed in the small watershed of Diediegou, located in the semiarid northern part (106°4′55″ E-106°9′15″ E, 35°54′12″ N-35°58′33″ N) of the Liupan Mountains, in the midwestern part of the Loess Plateau in Northwest China ( Figure 1). The elevation range is 1973-2615 m a.s.l. The climate is semiarid continental monsoon, with a mean annual air temperature of 6-7 °C, a frost-free period of 130 days, and a mean annual precipitation of 425 mm that mainly occurs from June to September.
This small watershed lies in the transition zone between the loess hilly area and the lithosol mountain area. The main soil type is gray cinnamon soil (partly corresponding to the Calciustoll, Argiustoll and Haplumbrept in the USA soil classification system). The native vegetation in this watershed is meadow and broadleaved deciduous forests, which had been completely destroyed in recent history [49]. Afforestation has been encouraged since the 1980s, mainly to stop severe soil erosion and produce timber. The dominant forest is now pure plantations of Larix principis-rupprechtii, which distribute only on the shady (north-facing) or semi-shady (northwestand northeast-facing) slopes in this semiarid study region [50].

Plot Setup
This study was carried out in a permanent larch plantation plot with an area of 30×30 m 2 at the foot of a shady slope at the lower reach of the small watershed. This plot has an elevation of 2055 m a.s.l., a slope gradient of 11°, and a soil thickness of >2 m and with granular sandy loam. Since the establishment of this plantation in 1986 through afforestation, no thinning has been implemented due to a strict logging ban policy. This is a young and even-aged plantation with one canopy layer and less tree size differentiation. From 2010 to 2014, the plot was characterized by a significant (p < 0.01) increase of canopy density from 0.80 to 0.85, mean tree height from 9.41 to 10.65 m, mean diameter at breast height (DBH) from 10.13 to 11.02 cm, total sapwood area from 8480 to 9910 cm 2 , and maximum LAI (LAImax) from 4.92 to 5.63 in midsummer (Table 1).

LAI and Sapwood Area Measurement
The canopy LAI was measured every 10 days during the study period from May to October with a plant canopy analyzer (LAI-2000 type, LI-COR, Lincoln, NE, USA) at 12 fixed sites in the plot.
After the measurement of DBH (cm), the sapwood area (As, cm 2 ) of each tree was calculated by Equation (1) [51] and then summed as the plot estimation: (1)

Weather Parameter Measurements
A weather station (Weatherhawk, Logan, Utah, USA) was placed in an open area 100 m away from the plot to automatically monitor the precipitation (P, mm), air temperature (Ta, °C), relative air humidity (RH, %), solar radiation (Rs, w·m −2 ), and wind speed (U, m·s −1 ) to calculate the PET (mm·day −1 ) with the Food and Agriculture Organization (FAO) Penman-Monteith Equation [52] every 5 min.

Soil Moisture Measurements
The soil water potential (ψ, MPa) for the layers of 0-10, 10-20, 20-40, and 40-60 cm were monitored with equilibrium tensiometers (EQ15; Ecomatik, Dachau, Germany) at a representative site in terms of canopy density and site conditions but at least 1 m away from tree trunks to avoid the influence of stem flow. Data were collected every 5 min by a logger (DL6; Delta-T Devices, Cambridge, UK), converted into volumetric soil moisture (θ, vol.%) with a series of equations relating ψ and θ determined at the same site [53], and then calculated into the relative extractable soil water content (REW) [54] with Equation (2). The reason for converting ψ to θ and REW, not directly using them, is to make the study results applicable for sites with different soil hydrological properties: where θact is the actual volumetric soil moisture in the main root zone of the 0-60 cm soil layer [55]; θWP is the θ of the wilting point (ψ = −1.5 MPa) with a value of 18.5%; and θFC is the θ of the field capacity (ψ = −0.01 MPa) with a value of 45.3% [36].

Sap Flow Measurements and T Calculations
The characteristics of tree height, DBH, clean bole length, and canopy diameter of all individual trees in the plot were measured in the early stage of the growing season of every year. Based on these measurements, 4-5 sample trees (Table 2) with straight trunks and healthy crowns were selected for sap flow measurements. There was no significant difference (p > 0.05) in the features of sample trees among different years. The sap flow density was measured using thermal diffusion probes (SF-L; Ecomatik, Dachau, Germany). Each set of probes consisting of four 20-mm-long sensors (S-1, a heated sensor powered by a constant current in 12 v voltage; S-0, S-2, and S-3 as reference sensors) was mounted at the breast height of the north side of the trunk after the outer bark was peeled off and then covered with aluminum foil to avoid physical damage and thermal influences from solar radiation. Data were recorded by a data logger (DL2; Delta-T Devices, Cambridge, UK) every 5 min.
The temperature difference between sensors was calculated using Equation (3). The Baseline 3.0 Software (Dr. Yavor Parashkevov, University of Duke, Durham, North Carolina, USA) was used to determine the max T ∆ . Then, the temperature data were transferred to sap flow density ( s J , where act T ∆ is the difference in instantaneous temperature (°C) between the heated sensor and reference sensors; the temperature differences between S-1 and S-0, S-1 and S-2, and S-1 and S-3, respectively. The daily T (T, mm·day −1 ) of the forest plot was calculated using Equation (5): The SF-L sensor is a product after improving the Granier sensor (SF-G) by subtracting the separately measured natural vertical background temperature gradient within sapwood (expressed as ( (3)) from the temperature difference between the two needles (expressed as 0 ΔT in Equation (3)) of the Granier sensor. Therefore, the SF-L sensor becomes significantly more accurate and the important reference point of zero sap flow (ΔTmax) remains stable over long periods. The exact determination of the point of zero sap flow requires additional continuous monitoring of an electronic dendrometer. Therefore, the SF-L sensor is not suitable for small trees with DBH < 8 cm. This led to a fact that the mean height and DBH of sample trees were more or less higher than the plot average. However, this should not be critical because the calculation procedure was the same among all years, the contribution of small trees to the daily T of the whole plot was mostly reflected by their sapwood area (see Equation (5)), and the mean sapwood area ratio of the small trees to the whole plot was just 12% within the study period. Therefore, we assumed that the systematic error in T estimation caused by the sample tree size higher than plot average has no effect on our analysis.

Derivation of the Daily T Model
The stomatal conductance control contains two aspects: plant physiological and hydraulic characteristics [57]. Physiological control is mainly concerned with stomatal response to photosynthetic drive and VPD under the water-sufficient condition [58]; while the hydraulic characteristics are mainly in the process of root water uptake-catheter transport-water loss from stomatal. Serious water stress can lead to xylem cavitation and embolism. Then, the hydraulic transmission structure of soil-root-xylem is destroyed, and the water conductivity will decline rapidly [59]. Based on above, we assume that no synergistic interactions exist among the influencing variables. The basic form of the daily T model in this study is a multiplication coupling the T response equations to PET, REW and LAI, as shown in Equation (6): where Tmax is the maximum T (mm·day −1 ) when PET, REW, and LAI show no limit to T at the same time; ƒ(i) is the equation describing the T response to factor i (PET, REW, or LAI).
The daily T variation is simultaneously influenced by many factors in the field; thus, it is difficult to derive the individual factor effect using raw field data. To see how T responded to each factor and following what function type, the disturbance from other factors not assessed was excluded through the fitting of the upper boundary line [60] based on the selected data points, which were above the values of the mean plus one standard deviation in each segment of the factor. The needles of larch trees gradually lost their vitality in October due to the decreasing temperature, but they remained at the crown and were still measured as LAI. Therefore, only the data from May to September were used in the analysis below.
A boundary line analysis (Ram Oren's H2O Ecology Group from Duke University) was performed to determine the function type of f(i), and then, the software SAS 9.4 (SAS Institute Inc., Cary, North Carolina, USA) was used to fit the parameters in Equation (6), with field-measured data of 2010 and 2012. Then, the model was validated with field-measured data of 2014. The model fit quality was assessed using both the absolute error (A) (Equation (7)) and relative error (B) (Equation (8)): where MTi is the field-measured T (mm·day −1 ); CTi is the T (mm·day −1 ) calculated by the model developed in this study; and n is the number of valid measurements. Additionally, the model fit quality was assessed by the non-dimensional efficiency criterion of Nash and Sutcliffe (E), as shown in Equation (9) [61]: where MTm is the mean of all measured T during the whole study period. E varies from minus infinity to 1, and the latter corresponds to a perfect fit.

Uncertainy and Sensitivity Analysis
The modelling approach contains the potential uncertainty source from PET, REW, and LAI. The quantification of the effect of these uncertainties on the T calculation requires the knowledge of the aforementioned model variables and of their statistical variability and correlation structure. In this study, the uncertainties were evaluated by computing a defined relative sensitivity (RS) [62], as shown in Equation (10): where CP is the relative change of a given variable (PET, REW, LAI), AS is the corresponding relative change of output (T), Ps and Pb are the postperturbed and original values of a given variable, respectively; Cb is the measured T, and Cs is the computed T using Equation (6) by inputting the given variable values with perturbation (i.e., PET increased or decreased by 10%) and the original values of other variables (i.e., REW, LAI). The perturbation magnitude (CP) was ±10% with respect to the original data.

Assessment of the Individual Factor Effects on Daily T
Based on the measured data of PET, REW, and LAI, missing daily T due to the failure of power supply or other reasons was calculated using Equation (6) to form a complete data set of daily T during the entire growing season. Then, the monthly T and total T of the growing season were calculated by summing daily T. The common reference to assess an individual factor effect on T (at the daily, monthly, and growing season scales) in year i (here the dry year 2010 and normal year 2012) is the corresponding T values (Tj) observed in year j (here the wet year 2014) under the reference conditions of PETj, REWj, and LAIj. The T in 2014 was least restricted by PET, REW, and LAI, since the REW and LAI in 2014 were the largest and the mean PET (3.04 mm·day −1 ) was slightly lower than its maximum (3.06 mm·day −1 ). The T limitation by only one factor can be calculated by inputting the actual measured value of this factor into Equation (6), while other factors remain at their reference conditions. For example, the T limitation by only PET (T(PETi)) is calculated by inputting the actual PETi and the reference values of other factors (REWj and LAIj) into Equation (6). Then, the difference between T(PETi) and Tj, i.e., ∆T(PETi), is viewed as the effect of only PET on T. The same process was used to calculate ∆T(REWi) and ∆T(LAIi). A positive or negative value represents a promotion or limitation effect. The relative effect on T was assessed by the ratios of ΔT(PETi)/Tj, ΔT(REWi)/Tj, and ΔT(LAIi)/Tj.

Environmental Conditions
The total precipitation in the growing season was 426, 521, and 592 mm in 2010, 2012, and 2014, distributed on 49, 35, and 39 rainy days, respectively ( Figure 2). The events with a precipitation ≤10 mm accounted for 73.5%, 54.3%, and 46.2% of the total precipitation, and 31.9%, 12.1%, and 9.0% of the number of precipitation days in 2010, 2012, and 2014, respectively. The precipitation events with a depth > 10 mm explained most of the interannual variation in precipitation.
In general, the PET showed the same behavior in the growing season of all study years when not considering the data of the rainy or rain-affected days (Figure 2). From early May, the PET increased rapidly until reaching its maximum in early or middle June, and then it started to decrease. The PET average (and range) was 2.63 (0.24-8.00), 3.06 (0.22-6.11), and 3.02 (0.14-7.17) mm·day −1 for the growing seasons of 2010, 2012, and 2014, respectively.
The volumetric soil water content (SWC) varied greatly within each growing season (Figure 2), mainly due to the uneven distribution of rainfall events > 10 mm. There was a two-month-long drought period in 2010, even in a month with an SWC ≤ 25%. In contrast, in 2014, the SWC was essentially always above 30%. The SWC in 2012 was at a moderate level. The mean SWC was 30.17%, 35.39%, and 37.59% in 2010, 2012, and 2014, respectively.

Variation in the LAI
The variation in forest canopy LAI in the growing season of all study years is shown in Figure 3. The LAI increased with rising temperature (and day of year (DOY)) until its maximum in the middle growing season, and then it decreased with the lowering temperature. However, the LAI variation pattern was somewhat different among individual years, e.g., the LAI decrease was much slower in 2014 mainly due to the delayed appearance of <5 °C air temperatures in late autumn. The mean LAI in the growing season of 2014 was 4.80, much higher than that in 2010 (4.05) and 2012 (4.09). The maximum LAI usually occurred on approximately the 205th day of each year. Likely due to the rainfall distribution and soil droughts in May 2010 and June 2012, the calculated maximum LAI using the regression relations in Figure 3

Variation in Daily T
The daily T variation in the growing season is shown in Figure 4. When considering only at the variation tendency, i.e., ignoring the relatively low T data for rainy or heavy cloudy days, the daily T increased from the beginning of growing season until its peak in late July (the 30th in 2010 and 2014 and the 27th in 2012) and then decreased sharply with lowering temperature (Figure 3)

T Response to Each Individual Factor
The T responded to PET in a parabolic way (Figure 5a), conforming to the binomial equation. First, the T increased with rising PET rapidly and almost linearly when PET was <4 mm·day −1 ; then, T increased slowly when the PET varied within 4-5.4 mm·day −1 ; and T ceased to increase or even began to decrease when PET was >5.4 mm·day −1 , probably due to the stomatal closure caused by excessive dry weather conditions [27]. The upper boundary line was derived as Equation (11) The T responded to REW following a saturated exponential function (Figure 5b), as expressed by Equation (12). First, T rapidly increased with rising REW when REW was <0.4; then, T increased much more slowly if REW was >0.4; and finally, T tended to saturate at 1.7847 when REW was close to 1: The T responded to LAI also following a saturated exponential function (Equation (13)), although it looks like a linear relation in Figure 5c due to the narrow variation range of LAI. The T increased gradually with rising LAI, and did not reach the saturated value in this study. However, when LAI tended to infinity, the T was close to its saturated value of 3.3407:

T Response to Multiple Factors
To depict the T response to the simultaneously varying factors of PET, REW, and LAI, a multiplication model (Equation (6)) was built based on the T response functions to individual factors (Equations (11)-(13)), and the parameters of this model (Equation (14)) were well fitted (R 2 = 0.89, Figure 6a,b) using the measured data in the growing seasons of 2010 and 2012: The sum of the daily T calculated using Equation (14) was 112.2 and 110.1 mm in the growing seasons of 2010 and 2012, respectively, and these values were only 2.09% lower and 5.96% higher than the measured sum of 114.6 and 103.9 mm, respectively. The corresponding Nash coefficients were 0.88 and 0.87, respectively.
According to this study (Table 3), the mean difference between calculated and measured daily T was <0.06 mm in the period from May to September in 2010. In 2012, the daily T difference in all months except September (0.15 mm) was essentially <0.09 mm. This indicates that the model has high precision and can be successfully used to explain and predict the daily T response of larch plantation to the varying environmental and canopy conditions.

Model Validation
The measured daily T data in 2014 were used to validate the T model (Figure 6c). The result indicated that the sum of modeled daily T is 17.28 mm (15.39%) higher than the sum of measured daily T, with an average daily T error of 0.16 mm, a Nash coefficient of 0.79, and a determinant coefficient (R 2 ) of 0.89. This suggests that the model is reliable and accurate for predicting daily T under changing PET, REW, and LAI. The positive difference was concentrated in the autumn of 2014 due to the unusually higher LAI resulting from the delayed fall of yellow needles.

Sensitivity Analyses of the T Model
The change of T caused by the perturbation (10%) of model variables (PET, REW, and LAI) varied among different years, and the variable that caused the greatest change of T was not always the same in different years (Table 4). If the PET, REW, and LAI was increased by 10%, the order of the sensitivity coefficients of the variables was PET > LAI > REW in each of the three years, and the largest sensitivity coefficient of PET appears in the wet year. This means that in the case of increasing perturbation of the three variables, the most sensitive variable to the T model is PET, especially in the wet year. If the PET, REW, and LAI decreased by 10%, the variable that caused the greatest change of T was PET, REW, and REW in the dry year 2010, normal year 2012, and wet year 2014, respectively. This meant that in the case of decreasing perturbation of the three variables, the most sensitive variable to the T model was REW in the normal year and wet year but PET in the dry year; moreover, the sensitivity coefficient of the decreasing REW to the T model increased with rising precipitation within the three years. The mean sensitivity coefficient showed that T was the most sensitive to PET, followed by LAI and REW.

The Effects of Individual Factors on T
Based on the complete data set of daily T from the entire growing season after the missing T data were filled up using the measured data of PET, REW, LAI, and Equation (14) (Table 5), as an integrated result of various factors.

Behavior of T Response to Individual Factors
For a given forest with fixed site and LAI conditions, the complex influence of many meteorological factors on daily T can be simplified as the effect of PET [25]. In this study, the daily T responded to PET in a parabolic way, and the increase rate of T declined with rising PET. This is an integrated effect of comprehensively interacting factors. On the one hand, the increase in PET (also VPD) with rising solar radiation and temperature promoted the T; on the other hand, the excessive VPD led to stomatal closure, which resulted in decreased stomatal conductance and T [27]. The variation trend of T with PET and the related PET thresholds in this study are in line with previous studies, such as the T of a mature oak (Quercus petraea) forest near Nancy in northeastern France [27], and the T and stomatal conductance of two rainforest species (Simarouba amara and Goupia glabra) growing in plantations in French Guyana [63]. This study detected a PET threshold of 4 mm·day −1 , above which the increase rate of T with PET was strongly reduced. A similar PET threshold was observed in the control experiments of black locust seedling T [25], with PET thresholds of 3.2 and 4.0 mm·day −1 for sandy loam and loamy clay, respectively. However, daily T always increased linearly with rising PET, without an apparent threshold, in studies such as those on the T of Pinus pinaster in south France and of a tropical rain forest in French Guiana [24,64] with T increase rates (T/PET) of 0.55 and 0.75, respectively, probably due to the very high air humidity in those sites where the VPD cannot reach a level leading to stomata closure.
In many previous studies, the relation between T and soil moisture or REW had been described by a segmented equation, logarithmic equation, or logistic equation [32,65,66]. In this study, T responded to REW as a saturated exponential function. This result is consistent with previous studies, such as that by Sadras and Milroy [67], Granier [68], and Granier et al. [43] on the forests of 21 coniferous or broadleaved tree species, showing that the ratio of canopy stomatal conductance to its maximum increased with rising REW until a threshold of 0.4, but thereafter, the ratio leveled off.
Studies on the T response to LAI within a large LAI variation range are still rare. Within the LAI range of <3, Bucci et al. found that the T response can be well described by a single asymptotic equation [69]. In a study on the forest T of six willow clones, T still increased with rising LAI in the range of 3-6 [70]. In studies on forest T response to thinning, pruning, and nitrogen fertilizer application [42], T increased linearly with rising LAI within a range of 1-6, but T was higher in the leaf expansion stage than in the leaf fall stage [71]. However, Granier et al. found that the canopy conductance increased linearly with rising LAI within a range of <6 but stopped increasing thereafter [43]. Based on the few case studies with high LAI, forest T increases with rising LAI until a threshold (e.g., 6.0) but at a declining rate; after this threshold is reached, T stops increasing. Because the maximum canopy LAI was lower than 6 in every year of this study, a gradually declining increase in T with rising LAI until an LAI of 5.3 occurred, but the effect of the LAI threshold was not clearly observed.

Daily T Model Coupling the Effects of Individual Factors
Forest T is codetermined by meteorological conditions, soil moisture, and vegetation characteristics. The variation in stomatal conductance of illuminated leaves in forest canopy was expressed as a continuous multiplication equation of photon flux density, temperature, VPD, leaf water potential, and ambient CO2 concentration [44]. A surface conductance model was established by multiplying the equation of LAI variation, the maximum surface conductance, and the nonlinear equations of environmental factors (solar radiation, specific humidity deficit, temperature, and soil moisture deficit) with values between zero and infinity [45]. Now, such an approach has been widely used in studies relating canopy conductance with influencing factors, such as for forests of oak (Quercus robur) in the Netherlands [46,47], maritime pine (Pinus pinaster) in southwestern France [23,48], sessile oak (Quercus petraea) in northeastern France [71], and 21 tree species under different climates (temperate, mountain, tropical, and boreal) [43]. The vast majority of such studies focused on the canopy conductance response to environmental factors, or on a certain aspect of the relations between transpiration and weather parameters, soil moisture, and vegetation structure. However, there are still relatively few studies directly relating the simultaneous responses of forest T to all main influencing factors. In this study, a daily forest T model was developed to couple the effects of PET, REW, and LAI on T, and this model can well describe the T response to these factors. A common feature of the above studies is the adopted hypothesis that no synergistic interactions exist among the influencing variables. Further experiments are needed to test whether this hypothesis is adequate. If it is not true, the synergistic interactions among variables should be studied.

Amplitude of Individual Factor Effect on T
It is difficult to directly separate and accurately evaluate the effect on T by individual factors based on field studies, since all factors vary simultaneously and the relative importance of each factor changes. In the studies of Zhao et al. [72] and Ungar et al. [28], the primary driver of T converted from SWC under lower soil moisture to PET under higher soil moisture. Čermák et al. [73] observed a T reduction of only 10% due to the REW deficit in a mixed forest of pine and spruce under moist conditions. The difference in the T reduction caused by soil drought might be influenced by numerous factors, including the actual REW and PET and the differences in their maximums, tree species, and site conditions [35]. In a 14-year experiment of throughfall exclusion (TFE) in the Quercus ilex forest in southern France, a 29% reduction in net precipitation input to soil decreased the annual T by 23% [74]. In a 14-year (2002-2015) TFE experiment of a tropical forest in the eastern Amazon, the biomass was lowered by approximately 40% due to the 14-year soil drought through the TFE treatment. The T on the plot with 50% TFE was decreased by 30% in 2015 compared with that on the control plot with normal rainfall. However, there was only a 5% difference between the observed T on the TFE plot in 2015 and the calculated T of the control plot if reflecting the 40% biomass reduction and related basal area [75]. In our study, the model calculated T limitation due to LAI was greater than that due to REW, even in dry or normal years. In summary, the study cases mentioned above imply that the T may be mainly affected by the soil moisture deficit, either directly (through REW reduction) or indirectly, e.g., through leaves biomass (LAI) reduction.
In the dry year of 2010 in this study (Table 5), except the weak PET promotion on T in June (absolute value of 1.0 mm, or relative value of 2.60%), all three factors showed a limitation on T in every month. The LAI limitation was relatively stable compared with the strongly varying effects of PET and REW, which were quickly affected by the randomly changing weather conditions, including the depth and temporal distribution of rainfall. Due to the stored soil water from rainfall and snow melt in the non-growing season, the soil moisture in the early growing season of 2010 was high; however, the soil moisture presented a decreasing tendency in the growing season ( Figure 2) due to the accumulated evapotranspiration, low rainfall, and less heavy rainy days. In the growing season of 2010, there were 52 drought (REW < 0.4) days, with 4, 7, 1, 13, and 27 days in each month from May to September. When the drought stress had a long duration, the soil moisture limitation played a dominant role in the T reduction. Bréda and Granier [71] found that the T/PET ratio decreased sharply and immediately as REW dropped below 0.4 in a sessile oak forest in France. In short, the main limiting factors of T varied temporally in the dry year of 2010, they were PET and LAI in the period of May-Jul., and REW in the later growing season (August-September) with increasing soil drought stress.
In the normal year of 2012 in this study, if considering only the amplitude of the monthly T deviation from the baseline 2014, the first driving factor was PET in the first half of the growing season (May to July) in 2012, but LAI in the second half of the growing season (August and September). It seems that the amplitude of the PET effect on T decreased with rising month order in 2012, while the amplitude of the LAI effect on T increased. Water stress is assumed to occur when REW drops below a threshold of 0.4, under which the T gradually decreases due to stomatal closure [68]. In the growing season of 2012, there were only 12 days under soil drought stress and they appeared only in June (4 days) and August (8 days). Due to the relatively high level of rainfall in 2012, there was no obvious soil drought stress compared with that in 2014, leading to a low and stable REW effect on T within the entire growing season of 2012. The PET presented both a promotion and limitation effect on the monthly T, while the LAI always presented a limitation. Therefore, the limitation on the monthly T in the growing season was mainly from LAI in 2012 with normal precipitation. There was a significant interannual and seasonal variation in the T deficit due to PET, REW, and LAI. Bréda and Granier noted that the seasonal fluctuations in T were related to LAI and soil water availability [71]; LAI was the main factor governing the T/PET ratio when soil moisture was not limiting (0.4 < REW < 1). T was limited by both annual and seasonal variation of precipitation [76].
The absolute amplitude of the relative effect on the T (Table 5) showed that from May to July, the early growing season with sufficient soil water, the main influencing factor for T was always the PET regardless of whether it was a dry (2010) or normal (2012) year. However, from August to September, it was REW in the dry (2010) year and LAI in the normal (2012) year, respectively. The difference between years was probably caused by the precipitation difference in this period.

Limitations of This Study and Suggestions for Further Research
The daily T model developed in this study was precise for 2010 and 2012 but less precise for 2014, mainly due to overestimation in the late growing season of 2014, when the air temperature was high, given the delay in air temperatures <5 °C, leading to a phenomenon in which leaf vitality had decreased (become yellow) and T capacity had notably reduced [77], but the leaves still remained on the crown. On the other hand, the current T model considered only the LAI quantity, without the biological vitality of leaves and its temporal variation. This can well explain the T overestimation in the late growing season of 2014. In addition, the needles were damaged by some diseases and insect pests (such as Adelges laricis Vall) in this study, especially in dry years. The resulting removal of sap from the germs and needles led to not only shorter and thinner new branches [78], lower LAI, and different LAI variation patterns but also needles covered by insect bodies and greasy remains, which consequently reduced the T capacity of leaves and altered the T-LAI relation. If the variation in LAI vitality can be estimated or even monitored by new instruments, then the T model can be updated to reflect the leaf vitality variation by adding an item reflecting the varying T ability of leaves based on physiological studies, and the model prediction accuracy will be greatly improved.
In the analysis and T model development in this study, only the mean REW of the main root zone (0-60 cm) was used [55], while the difference in the contribution to T by the variation of both REW and fine root quantity along soil layers was not considered. This is an important cause of the reduced sensitivity of the T model to soil moisture variation. In future studies, the model structure should be further improved to reflect the varying contribution to T by the changing soil water availability in different soil layers. If the tree roots can uptake the groundwater, then this water resource should also be considered in the future models. However, no groundwater exists on the slopes of this study, due to the very high percolation rate of the loose bedrocks in Liupan Mountains. In addition, the model is firmly bound to a pure plantation of larch, a deciduous tree species, so it may not apply to other tree species, especially the evergreen tree species and mixed forests.
Because of the constraint of the sap flow sensors used in this study, the sap flow of small trees with a DBH below 8 cm was not measured. In general, the sap flow density increases with rising DBH till to the medium DBH of trees and then decreases with further rising DBH of large trees [17]. In the plantation plot of this study, the number ratio of small trees (DBH < 8 cm) was as high as 36%, so that the absence of small sample trees surely led to an overestimation of the T. On the other hand, the sapwood area ratio of small trees was just 12%, and the T contribution of small trees was partly reflected by the sapwood area of small trees when calculating the T using Equation (5), so that the overestimation of T should be not obviously high. However, this source of error should be excluded by a combined use of different sap flow sensors in future studies.
The T response to environmental factors varies along site conditions (e.g., slope aspect, slope position, slope gradient, and elevation). Due to the lack of observation at different sites, we cannot directly test the applicability and accuracy of the T model when used at other sites. However, the model variables of PET, REW, and LAI will vary accordingly with site conditions. Therefore, we believe that the model should be applicable to predict the T at different sites. Anyway, more observations at different sites should be encouraged in the future to test the fitness of this model when applied in a broader context, and to improve this T model further.
Climate change will increasingly influence forest T in both direct and indirect ways. The direct ways include not only variation in numerous meteorological factors (e.g., precipitation, temperature, wind speed, and solar radiation) and their temporal distribution but also increased frequency and intensity of extreme weather events and related damages, such as extreme high/low temperatures, drought and rainstorms, sleet and snowstorms, and late/early frost. The indirect ways include a strengthened change in the biologic rhythm of leaf development, greater variation in LAI and its vitality, and increased amounts of disease and insect pests. All these influences will alter the quantity and variation regime of PET, REW, and LAI. Quantifying the response of PET, REW, and LAI to these climate change-related stresses, or even enhanced influences from human activities, based on multi-and trans-disciplinary studies and coupling the response relations into the model will improve the prediction ability of the daily T model under changing environmental and canopy conditions.

Conclusions
The complex influences of many environmental and vegetation factors on forest T can be simplified into three aspects, i.e., PET, REW, and LAI. PET reflects the atmospheric evapotranspiration driver mainly determined by meteorological factors; REW reflects the soil water supply capacity mainly determined by soil moisture and soil physical properties; LAI reflects the forest capacity of water consumption mainly determined by the leaves amount. The T response to rising PET follows a binomial equation, describing the T increase with a declining rate until the PET threshold of 5.4 mm·day −1 ; the T response to both REW and LAI follows a saturated exponential equation, describing the rapid T increase until the threshold of REW (0.4) and LAI (4). These T response equations to individual factors were coupled to form a mechanism-based and multiplicative model of the daily T of the studied larch plantation. Both the calibration and validation of this model showed a satisfactory accuracy. Therefore, it can predict the T response to changing environmental and canopy conditions, and separate the contribution of individual factors or their combinations to the T in different time scales. Using this T model and taking the baseline of the wet year of 2014 with the highest T, the contribution of individual factors to the T change was assessed in the monthly and annual (growing season) scales. For the total T in the growing season of the dry year of 2010, all individual factors had a limitation effect, with the order of PET > LAI > REW; while it showed a limitation effect of LAI and REW and a weak promotion effect of PET in the normal year of 2012, with the limitation order of LAI > REW > PET. It seems that the limitation effect of all three factors is bigger in the drier year than in the wetter year. However, the effect of individual factors on the monthly T varied among observation months and years. In the dry year of 2010, the main limiting factors to T changed from PET and LAI in the first half growing season to REW in the later growing season. In the normal year of 2012, the PET effect on T changed from limitation in May and July to promotion in June, August, and September, and the effect amplitude decreased as the month went on; the REW effect of T limitation was low and stable; while the LAI effect on T was always limitation and its amplitude increased. All these led to an integrated result that the limitation to monthly T was mainly contributed by canopy LAI.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: