Loblolly Pine Productivity and Water Relations in Response to Throughfall Reduction and Fertilizer Application on a Poorly Drained Site in Northern Florida

Abstract

Loblolly pine (Pinus taeda L.) forests are of great ecological and economic value in the southeastern United States, where nutrient availability frequently limits productivity. The impact of fertilizer application on the growth and water relations of loblolly pine has been investigated by numerous studies; however, few field experiments have examined the effects of drought. Drought is of particular interest due to the potential for climate change to alter soil water availability. In this study, we investigated the impact of fertilizer application and a 30% reduction in throughfall on loblolly pine productivity, transpiration, hydraulic conductance, and stomatal conductance. The study was installed in a ten-year-old loblolly pine plantation on a somewhat poorly drained site in northern Florida. Throughfall reduction did not impact tree productivity or water relations of the trees. This lack of response was attributed to abundant rainfall and the ability of trees to access the shallow water table at this site. Fertilizer application increased basal area production by 20% and maximum leaf area index by 0.5 m²∙m⁻², but it did not affect whole-tree hydraulic conductance or the sensitivity of stomatal conductance to vapor pressure deficit. During the spring, when leaf area and vapor pressure deficit were high, the fertilizer-only treatment increased monthly transpiration by 17% when compared to the control. This relationship, however, was not significant during the rest of the year.

1. Introduction

The concentration of carbon dioxide (CO₂) and other greenhouse gases in the Earth’s atmosphere have increased markedly due to the combustion of fossil fuels and land use changes. The concentration of atmospheric CO₂ has increased from 280 parts per million (ppm) to over 400 ppm since the start of the Industrial Revolution [1,2]. Climate change simulations using the assumptions of the Intergovernmental Panel on Climate Change A2 scenario, which represents a largely business-as-usual future, suggest that the southeastern U.S. will experience an increase in average temperature of 2.5 to 4.7 °C by 2085 [3]. Estimates for changes in precipitation for this region are much less certain and range from a 22% reduction to a 9% increase in summertime (June, July, August) precipitation by 2070–2099 [3]. Increased average temperature and potential decreases in precipitation would decrease soil water availability and increase vapor pressure deficits (D). These changes could influence the carbon gain and water use of loblolly pine (Pinus taeda L.) plantations due to the sensitivity of stomatal conductance to these variables.

Loblolly pine accounts for 84% of all seedlings planted in the southeast United States (U.S.) and is the nation’s leading timber species [4,5]. Within its native range, which extends from Texas eastward to Florida and northward to Delaware, there are approximately 12 million hectares of loblolly pine plantations [4,6]. These plantations have tremendous ecological and economic value. Southeastern forests sequester 76 Tg of carbon each year, which represents 13% of regional greenhouse gas emissions [7], and they play an important role in regulating water quality and yield [8,9]. Forestry is vital to the economy of the South, providing 5.5% of the jobs and 7.5% of the industrial output [10].

The productivity of loblolly pine plantations is often limited by the availability of essential nutritional elements, especially phosphorus and nitrogen [11]. As a result, the application of fertilizer to these systems has become a common practice and has been shown to increase stemwood production and leaf area index (LAI) on a diversity of sites [11,12,13]. The increases in LAI associated with fertilizer application have also been shown to increase stand water use [14,15,16]. As the climate of the southeastern U.S. begins to change the impact of fertilizer application on the water relations of loblolly pine plantations will be of particular interest due to their extent and role in regulating watershed yield [9,17].

Several studies have investigated the impact of genetics, fertilizer application, and irrigation on the growth and water relations of loblolly pine [15,16,18,19,20], but few have investigated the influence of reduced water availability on these factors [21]. Given the economic and ecologic importance of loblolly pine forests, it is important to understand how potential changes in precipitation will impact their growth and water relations and how these responses interact with soil nutrient availability. The overall objective of this study was to investigate the influence of throughfall reduction and fertilizer application on loblolly pine productivity, transpiration, whole-tree hydraulic conductance, and sensitivity of canopy-level stomatal conductance to D. To meet these objectives a replicated throughfall reduction by fertilizer application experiment was installed in a loblolly pine plantation in northern Florida.

2. Materials and Methods

2.1. Study Design

This study took place at the Pine Integrated Network: Education, Mitigation, and Adaptation Project (PINEMAP) throughfall × fertilizer application study in Taylor County, Florida (30°12′22″ N, 83°52′12″ W). The study design is replicated at sites in Georgia, Oklahoma, and Virginia (see Will et al. [22] for details). The site is located on a Melvina-Moriah-Lutterloh soil complex which is characterized as a somewhat poorly drained fine sand [23]. The site was planted in the winter of 2003–2004 with a seed orchard mix of St. Joe Timber Company open-pollinated, genetically improved stock at a spacing of approximately 1.8 × 3.0 m.

The study design is a randomized complete block containing two levels of fertilizer application and throughfall reduction in a 2 × 2 factorial arrangement, with four replicate blocks. Blocks were selected by minimizing the variability of plot-level basal area. Plots treated with fertilizer received 224N, 27P, 56K (kg·ha⁻¹) and 22.4 kg·ha⁻¹ of a micronutrient blend consisting of 6% sulfur, 5% boron, 2% copper, 6% manganese, and 5% zinc [22]. Fertilizer was applied by hand in April 2012. Plots treated with throughfall reduction received a 30% exclusion of incoming throughfall. A 30% reduction in throughfall was selected as this represented the extreme for predicted changes in drought severity within the United States at the time the study was established [24]. Throughfall exclusion construction was completed in April 2012. The factorial design results in four individual treatment types. The treatment types were:

(1)

Control (C): ambient throughfall, no fertilizer application

(2)

Throughfall Reduction (TR): 30% throughfall reduction, no fertilizer application

(3)

Fertilizer application (F): ambient throughfall, fertilizer application

(4)

Fertilizer application and Throughfall Reduction (TR+F): 30% throughfall reduction, fertilizer application.

All measurements were taken in 14.6 × 16.8 m measurement plots. The mean number of trees within measurement plots was 47. Measurement plots were surrounded by a 6 m treated buffer on all sides, followed by an additional 9.1 m buffer around the treatment plot. The gross area of the measurement and treated plots were 0.025 ha and 0.047 ha, respectively. All plots received a broadcast application of 140 mL∙ha⁻¹ of metsulfuron methyl and 4 L∙ha⁻¹ of glyphosate to remove competing vegetation prior to treatments being applied. An inventory of all trees in measurement plots was conducted in the dormant season of 2011–2012 to assess baseline conditions. There were no significant differences in the basal area, mean diameter at breast height (DBH, 1.37 m), height, or trees per hectare (p > 0.05) among treatments at this time. Additional inventories were conducted in the dormant season of 2012–2013 and 2013–2014 to assess the impact of treatments on stand growth. Inventories included measurements of DBH and total tree height. The sapwood area of all trees was estimated using an allometric equation relating sapwood area to DBH as reported in Gonzalez-Benecke and Martin [19]. In the original publication the equation was misreported as the sign of the intercept was incorrect. The corrected equation is:

$$\ln \left({\mathrm{A}}_{\mathrm{S}}\right)=-15.347+2.185\ln \left(\mathrm{DBH}\right)$$

(1)

where A_S is stem sapwood area in m² and DBH is expressed in mm.

2.2. Throughfall Exclusion Design

Throughfall exclusion troughs were installed parallel to pine rows and were constructed using pressure-treated lumber and 0.3 mm extrusion laminate plastic sheeting with 2 layers of U.V. stabilized coextruded polyethylene and a high strength cord grid (Poly Scrim 12, Global Plastic Sheeting Inc., San Diego, CA, USA). Each throughfall exclusion structure contained two 0.5 m wide troughs separated by a 0.4 m opening. The opening between the exclusion troughs was designed to allow some throughfall to access the area beneath each trough in order to reduce soil moisture heterogeneity. Exclusion structures were approximately 0.9 m in height at the low end and 1.4 m on the high end. PVC pipes and polyethylene drainage pipes were installed at the low end to remove the excluded throughfall from the gross plot area. The structures were installed such that 30% of the ground area was covered by troughs.

2.3. Meteorological and Soil Measurements

A weather station was installed on a tower above the stand in August of 2012 to measure photosynthetic photon flux density (PPFD) (PQS1- par quantum, Kipp and Zonen USA Inc., Bohemia, NY, USA), temperature and relative humidity (HMP45C, Campbell Scientific Inc., Logan, UT, USA) and precipitation (TE525MM-L, Texas Electronics Inc., Dallas, TX, USA). All sensors collected measurements every minute, and half-hour averages were stored on data loggers (CR10x, Campbell Scientific Inc.).

Soil volumetric water content was measured monthly from November 2013 to April 2014 using a Tektronix 1502B Metallic Time-Domain Reflectometer. Measurements were taken on tungsten inert gas (TIG) welding rods of varying lengths installed vertically into the soil. In each plot three pairs of TIG rods (30 cm, 60 cm, and 84 cm lengths) were installed at two locations. One set of TIG rods was installed in the tree row and another set between tree rows. In plots containing throughfall exclusions the set of TIG rods between the tree rows was placed directly under an exclusion trough. The location of the TIG rods in the tree row was randomly selected and the set of TIG rods between the trees rows were installed adjacent to this location. Measurements were converted to volumetric water content using the equation reported by Topp et al. [25]. The volumetric water content of the 30–60 cm (θ_30–60) and 60–84 cm (θ_60–84) soil regions was calculated as the difference of each reading from the reading of the next shallowest soil depth, assuming that the original reading represented the weighted average of the entire measured soil profile:

$${\mathsf{\theta }}_{30\text{–}60}=\frac{{\mathsf{\theta }}_{0\text{–}60}-(0.5\times {\mathsf{\theta }}_{0\text{–}30})}{0.5}$$

(2)

$${\mathsf{\theta }}_{60\text{–}84}=\frac{{\mathsf{\theta }}_{0\text{–}84}-(0.7143\times {\mathsf{\theta }}_{0\text{–}60})}{0.2857}$$

(3)

where θ_0–60 is the volumetric water content of the 0–60 cm soil region and θ_0–84 is the volumetric water content of the 0–84 cm soil region.

The volumetric water content of the top 20 cm of soil (θ₂₀) was measured monthly from April 2013 to December 2013 using a Campbell Scientific CS620 water content sensor. Measurements were taken in the tree row, a quarter of the distance between tree rows, and half of the distance between tree rows. For plots with throughfall exclusion structures this represented readings in the tree row, at the edge of an exclusion trough and directly under a trough.

Soil bulk density and texture was measured in the top 20 cm of soil in the winter of 2012–2013 using a 150 cm³ soil bulk density sampler. Soil samples were collected for the 0–10 cm and 10–20 cm depths at three locations in each of the 16 plots. One sample was taken in the tree row, one was taken a quarter of the distance between tree rows, and one was taken half of the distance between tree rows. Soil samples were weighed after being oven-dried at 105 °C for 48 h. Bulk density was then calculated as soil dry weight divided by soil volume. To measure soil texture, samples were collected from eight random locations in each plot for each depth (0–10 cm and 10–20 cm). The eight samples were then randomly split into two groups for each depth and mixed, resulting in two samples per depth per plot. All samples where then homogenized, sieved to 2 mm, and roots and organic matter set aside. Forty grams of soil was taken from each sample for texture analysis using the hydrometer method outlined by Gee and Bauder [26].

2.4. Leaf Area Index, Specific Needle Area, and ANPP

Needlefall was collected approximately every four weeks from June 2012 to April 2015 from 12 0.5-m² needlefall traps randomly located in each plot. After collection any non-needle material was removed from the sample and the remaining needles were oven dried for at least one week at 70 °C and weighed. Projected leaf area index (LAI) for each collection period was calculated from needlefall as described in Martin and Jokela [27]. This method also calculated daily LAI estimates.

Specific needle area was measured on 5 trees per plot. Ten fascicles were collected from the top third of the crown of each sample tree on 13 March 2014. Needles were cut to 10 cm and the radius of each needle in each of the 10 fascicles per tree was measured. For each individual needle the radius was measured on the two flat sides and then averaged. All needles were then dried at 65 °C for four days and weighed. Specific needle area was calculated as the ratio of the surface area of all needles measured for a given tree divided by the dry weight of all the measured needles.

Aboveground net primary productivity (ANPP) for the 2013 growing season was calculated as in Martin and Jokela [27]. Aboveground woody biomass increment was calculated using yearly inventory measurements and the allometric equations reported by Gonzalez-Benecke et al. [28]. Foliage biomass production was determined from needlefall data.

2.5. Sap Flow and Transpiration

Sap flux density (J_S, g∙m⁻²∙s⁻¹) was measured in the stem xylem on a subset of five trees per plot (80 trees in total) throughout 2013 using Granier-style heat dissipation sensors [29,30]. The temperature difference (ΔT) between these probes was measured every minute, and half hour averages were stored on data loggers. Values of ΔT were converted to J_S using the empirical calibration developed by Granier [29] and confirmed by Clearwater et al. [31]:

$${\mathrm{J}}_{\mathrm{S}}=119\times {\left(\frac{{\Delta \mathrm{T}}_{\mathrm{m}}-\Delta \mathrm{T}}{\Delta \mathrm{T}}\right)}^{1.231}$$

(4)

where ΔT is the instantaneous temperature difference between the reference and heated probes and ΔT_m is the maximum temperature difference. Values of ΔT_m were determined on nights when D was below 0.2 kPa for at least two hours, as there was little to no sap flow during these periods.

Measurement trees were selected in the spring of 2012 using the “quantiles of total” technique [32]. A summary of the minimum, maximum, and mean DBH of sap flow measurement trees at the start of 2013 is given in

Table 1. Minimum, maximum, and mean DBH of sap flow measurement trees during the winter of 2012–2013 for the control (C), throughfall reduction (TR), fertilizer application (F), and throughfall reduction and fertilizer application (TR+F) treatments.

	DBH (cm)
Treatment	Minimum	Maximum	Mean
C	10	19.1	14.9
TR	8.6	19.5	15.0
F	9.6	19.3	14.8
TR+F	9	19.6	15.3

. All sensors were inserted into the north side of the tree and covered with reflective insulation to minimize the influence of solar radiation on the temperature of the stemwood. To account for radial variation in J_s, the measurement trees closest to the quadratic mean diameter of 8 plots (2 per treatment) were outfitted with specialized probes to measure J_s at the 2–4 cm and 4–6 cm depths from February 2014 to April 2014. These probes were constructed using longer needles of the same gauge and material as the needles used in 20 mm probes.

An equation was developed to correct estimates of J_S measured in the outer 2 cm of sapwood to account for radial variation using the Q_T method described in Gonzalez-Benecke and Martin [19]. For trees outfitted with radial variation probes, corrected half-hourly estimates of transpiration (E_RAD) were calculated by multiplying the J_S measured in each region of sapwood (0–2 cm, 2–4 cm, and >4 cm) by the corresponding sapwood area of each region. Uncorrected transpiration estimates (E₂₀) were then calculated by multiplying the entire sapwood area by the value of J_S measured in the outer 2 cm of sapwood. A linear regression model was then developed to predict corrected transpiration estimates (E_COR) from E₂₀, D, and PPFD:

E_COR = −0.00065 + 0.88497 × E₂₀ + 0.0087 × D − 0.00001 × PPFD

(5)

This equation accounted for 99% of the variation in E_COR. The variables D and PPFD were included in this model as they have been shown to affect the radial pattern of sap flow in loblolly pine [19].

After accounting for radial variation, transpiration for each measurement tree was calculated by multiplying the corrected J_S estimates by total sapwood area. Plot-level transpiration estimates were calculated by using the proportion of total sapwood area of measurement trees in a plot to the total sapwood area of all trees within a plot.

Gaps that occurred in the J_S dataset due to sensor error and power failure were filled using a variety of models. Regression models were constructed for each individual measurement tree relating the J_S of each tree to the overall mean J_S of all trees within a given treatment. The R² of these models ranged from 0.87 to 0.99 with a mean of 0.96. If a data gap existed for any individual tree while sensors in other trees of the same treatment were functioning properly, values of J_S were filled using these regression equations. For the entire data set of all 80 trees, 21% of the J_S data points were filled in this manner. If all of the sensors at the site were inoperable, due to power failure, then values of J_S were filled from a nonlinear function which predicted J_S from meteorological data. The R² of these models ranged from 0.86 to 0.95 with a mean of 0.92. For the entire data set of all 80 trees, 14% of the J_S data points were filled in this manner. Power failure occurred at the site during the entire month of December and values of J_S were not calculated. For this month daily estimates of transpiration (E_Day, mm) were filled using a model to predict E_Day from meteorological data. The mean R² of these models was 0.82.

2.6. Canopy-Level Stomatal Conductance

Canopy stomatal conductance (G_C, mmol∙m⁻²∙s⁻¹) was calculated using the simplified version of the Penman-Monteith equation reported in Samuelson et al. [33]:

$${\mathrm{G}}_{\mathrm{C}}=\frac{{\mathsf{\lambda }\mathsf{\gamma }\mathrm{E}}_{\mathrm{L}}}{{\mathrm{D}\mathsf{\rho }\mathrm{c}}_{\mathrm{p}}}$$

(6)

where λ is the latent heat of vaporization (2465 J∙g⁻¹), γ is the psychometric constant (65.5 Pa∙K⁻¹), ρ is the density of air (1225 g∙m⁻³), c_p is the specific heat of air (1.01 J∙g⁻¹∙K⁻¹), E_L is transpiration per unit leaf area (mm∙m⁻²∙s⁻¹), and D is vapor pressure deficit (Pa). Transpiration per unit leaf area was calculated by dividing canopy transpiration (mm∙s⁻¹) by projected leaf area index. Values of G_C were calculated only when D > 0.6 kPa due to errors associated with calculations of G_C at low D [34].

The sensitivity of G_C to changes in D was described by fitting two different models previously developed to describe this relationship. The first model is reported in Ewers et al. [35]:

$${\mathrm{G}}_{\mathrm{C}}{=\mathrm{G}}_{\mathrm{Cref}}-\mathsf{\delta }\mathrm{lnD}$$

(7)

where G_Cref is the maximum G_C in mmol∙m⁻²∙s⁻¹ at D = 1 kPa and δ is the slope of the boundary line of the relationship between G_C and ln D (kPa). The second model is used in the Physiological Principles in Predicting Growth model (3-PG) to describe the sensitivity of G_C to changes in D:

$$\ln \left({\mathrm{G}}_{\mathrm{C}}\right)=\ln \left(\mathrm{MaxCond}\right)-\mathrm{CoeffCond}\times \mathrm{D}$$

(8)

where MaxCond is the maximum G_C in m∙s⁻¹ at D = 0 and CoeffCond defines stomatal response to D in mb⁻¹ [36,37]. Both of the models were fit using the quantile regression procedure in SAS 9.3 with a quantile threshold of 0.98 (SAS Inc., Cary, NC, USA) [37].

Whole-tree hydraulic conductance (Ks, mol∙m⁻²∙s⁻¹∙MPa⁻¹) was calculated for individual trees as the slope of the linear regression of J_S and leaf water potential (Ψ_L) [38]. Measurements of Ψ_L were taken on shoot tips every two to three months using a portable pressure chamber (PMS 1000, MPS Instrument Co., Corvallis, OR, USA). Two of the four blocks were sampled on each measurement date due to time limitations. Samples were taken from four sap flow trees per plot from pre-dawn to late afternoon on each measurement day. Measurements were taken at approximately 2 h intervals resulting in four to five measurements being completed for each tree on each day. Shoot tips were collected using a pole pruner and were stored in plastic bags with wet paper towels in order to minimize desiccation. All measurements were completed within 3 min of shoot excision. The hydraulic conductance of each tree was then calculated and used to produce plot level averages for statistical analysis.

2.7. Statistical Analysis

Analysis of variance was used to analyze the effects of treatments on growth and water relations traits, including Bonferroni adjustments. (PROC MIXED; SAS Institute Inc., Cary, NC, USA). The linear model for the analysis was:

$${\mathrm{Y}}_{\mathit{ijk}}=\mathsf{\mu }+{\mathrm{F}}_i+{\mathrm{TR}}_j+{\mathrm{FTR}}_{\mathit{ij}}+{\mathrm{B}}_k+{\mathsf{\epsilon }}_{\mathit{ijk}}$$

(9)

where Y_ijk is the parameter value of the plot with the ith level of fertilizer application (F) (0, 1), the jth level of throughfall reduction (TR) (0, 1), in the kth block (B) (1, 2, 3, 4) and μ is the population mean and ε_ijk is the error term ~NID (0, σ²_ε) (NID, normally and independently distributed). For a given main effect a value of 0 represented not receiving a treatment while a value of 1 represented receiving the treatment. Repeated measures analysis of variance was used to analyze time series data with an autoregressive covariance structure. Values of whole-tree hydraulic conductance were scaled to the plot level by averaging the values for the four trees measured in each plot. A significance level of p < 0.05 was used in all analyses.

3. Results

3.1. Soil Bulk Density, Texture, and Moisture

The bulk density and texture of the top 20 cm of soil was relatively consistent among treatments. The mean soil texture was 95% sand, 4% silt, and 2% clay. There were no differences in sand (p = 0.232) or silt (p = 0.198) content among treatments, but there were differences in clay content such that the TR+F treatment had significantly more clay content than the control (0.8% vs. 3.4%)(p = 0.025). There were no other differences among treatments in clay content (p > 0.11). Soil bulk density did not differ among treatments for the 0–10 cm (p = 0.220) or 10–20 cm (p = 0.429) depths and averaged 1.12 and 1.2 Mg∙m⁻³ respectively.

Throughfall reduction lowered the volumetric water content of the upper 20 cm of soil (θ_0–20), but this effect was not apparent at deeper depths. Values of θ_0–20 were reduced by a mean of 0.014 cm³_H2O∙cm⁻³_soil in response to throughfall reduction (p = 0.001). This estimate was calculated as the weighted mean of three θ_0–20 measurements taken at different positions on each measurement date. When analyzed individually, throughfall reduction significantly decreased θ_0–20 in the tree line by 0.01 cm³_H2O∙cm⁻³_soil (p = 0.030), at the edge of exclusion troughs by 0.023 cm³_H2O∙cm⁻³_soil (p = 0.002), and directly under exclusion troughs by 0.027 cm³_H2O∙cm⁻³_soil (p < 0.001). Throughfall reduction also lowered the volumetric water content of the soil region directly under exclusion troughs by 0.048 cm³_H2O∙cm⁻³_soil in the 0–30 cm soil depth (θ_0–30, p ≤ 0.001) and by 0.032 cm³_H2O∙cm⁻³_soil in the 30–60 cm (θ_30–60, p = 0.022) soil depth. There was no effect of throughfall exclusion on the soil moisture of the 60-84 depth (θ_60-84, p = 0.084). When measured in the tree row, plots treated with throughfall exclusion had 0.032 cm³_H2O∙cm⁻³_soil less volumetric water content in the 30–60 cm depth (p = 0.022) while θ_0–30 and θ_60–84 were unaffected (p > 0.134).

3.2. LAI, Needlefall, and Growth

Projected leaf area index (LAI, m²·m⁻²) was greater in plots treated with fertilizer application in both the 2012 (p = 0.037) and 2013 (p < 0.001) growing seasons. Throughfall reduction did not significantly affect LAI. In both years LAI was at a maximum at the end of July and at a minimum at the beginning of March (Figure 1). Maximum values of LAI for fertilized and unfertilized plots were 3.7 and 3.4 in 2012 and 4.2 and 3.6 in 2013, respectively (

Table 2. Mean aboveground net primary productivity (ANPP), basal area growth (BA_G), projected leaf area index (LAI), and needlefall for control (C), throughfall reduction (TR), fertilizer application (F), and throughfall reduction and fertilizer application (TR+F) treatments with associated p values. Values for LAI represent monthly minimums and maximums within a given year. For the LAI time series data, there were no treatment by month interactions. For significant F×TR interactions, significantly different (p < 0.05) values are indicated by different small letters.

		Treatment			p-Value
	C	TR	F	TR+F	F	TR	F×TR
ANPP 2013 (Mg∙ha−1∙year−1)	26.5 a	28.1 ac	32.3 b	30.2 bc	<0.001	0.736	0.049
BAG 2012 (m2∙ha−1∙year−1)	4.2	4.0	5.0	4.7	0.007	0.230	0.789
BAG 2013 (m2∙ha−1∙year−1)	3.4	3.5	4.3	3.9	0.033	0.576	0.314
LAI 2012 (Projected)	2.0–3.4	1.9–3.3	1.9–3.6	1.9–3.7	0.037	0.938	0.308
LAI 2013 (Projected)	2.3–3.6	2.3–3.6	2.7–4.1	2.7–4.2	<0.001	0.452	0.478
Needlefall 2013 (Mg∙ha−1∙year−1)	5.5	5.6	6.4	6.6	<0.001	0.383	0.536

). Total annual needlefall in 2013 was significantly greater (p < 0.001) in plots treated with fertilizer (6.5 Mg∙ha⁻¹∙year⁻¹) compared to untreated plots (5.5 Mg∙ha⁻¹∙year⁻¹, Table 2). Specific needle area was unaffected by both fertilizer application (p = 0.465) and throughfall reduction (p = 0.693). The mean specific needle area was 104.8 cm²∙g⁻¹.

Basal area increment was greater in plots treated with fertilizer in both 2012 and 2013 (p = 0.007, p = 0.033) and was unaffected by throughfall reduction (p = 0.230, p = 0.576). Growth rates were higher in 2012 than in 2013 and ranged from 4.0 to 5.0 m²∙ha⁻¹∙year⁻¹ and 3.4 to 4.3 m²∙ha⁻¹∙year⁻¹ for each year, respectively (Table 2). The TR+F treatment had a lower absolute growth rate than the F treatment in both 2012 and 2013, but this difference was not significant (p = 0.298, p = 0.273). Height growth was unaffected by the treatments and averaged 1.4 m∙year⁻¹.

Values of aboveground net primary production (ANPP) were 26.5, 28.1, 32.3, and 30.2 Mg∙ha⁻¹∙year⁻¹ for the C, TR, F and TR+F treatments respectively. There was a significant F×TR interaction for ANPP (p = 0.049) such that the ANPP of the F and TR+F treatments were both greater than the control (p < 0.001, p = 0.012), but only the F treatment was greater than the TR treatment (p = 0.006, Table 2). There was no difference between the F and TR+F treatments (p = 0.096).

3.3. Whole-Tree Hydraulic Conductance

Whole-tree hydraulic conductance was not affected by fertilizer application (p = 0.276) or throughfall reduction (p = 0.425). There was, however, a significant effect of measurement date on Ks (p < 0.001, Figure 2). Throughout the measurement period Ks ranged from 1.63 mol·m⁻²·s⁻¹·MPa⁻¹ in August 2013 to 2.21 mol·m⁻²·s⁻¹·MPa⁻¹ in October 2012. Similar to Ks, pre-dawn water potential (Ψ_p) was not significantly affected by fertilizer application (p = 0.372) or throughfall reduction (p = 0.093), but did vary through time (p < 0.001). Values of Ψ_p ranged from −0.34 MPa to −0.53 MPa.

3.4. Radial Variation in Sap Flow

To examine the relative rates of sap flow at different depths in the xylem, we compared daily summed sap flow per unit sapwood area (J_s) at 2–4 cm and 4–6 cm depths with daily sums in the outer 0–2 cm depth. Daily sums of J_s decreased significantly with depth (p < 0.001). The J_S of the 2–4 cm depth was on average 87% of the 0–2 cm J_s in all treatments (Figure 3). The J_s of the 4–6 cm depth was 58% of the 0–2 cm J_s for the C, TR and TR+F treatments, but 76% for the F treatment. Daily J_S of the 0–2 cm depth did not differ among treatments (p = 0.598).

3.5. Transpiration

Total transpiration for 2013 (E_Year) was not affected by fertilizer application (p = 0.349) or throughfall reduction (p = 0.544). Values of E_Year were 535, 547, 596, and 546 mm∙year⁻¹ for the C, TR, F, and TR+F treatments, respectively (

Table 3. Water use and stomatal sensitivity to vapor pressure deficit reported as mean annual transpiration (E_Year), maximum daily transpiration (E_Max), sensitivity of stomatal conductance to vapor pressure deficit (-δ), stomatal conductance at vapor pressure deficit = 1 kPa (G_Cref), 3-PG parameter defining maximum stomatal conductance (MaxCond), and 3-PG parameter defining stomatal sensitivity to vapor pressure deficit (CoeffCond) for control (C), throughfall reduction (TR), fertilizer application (F), and throughfall reduction and fertilizer application (TR+F) treatments with associated p values.

		Treatment			p Value
	C	TR	F	TR+F	F	TR	FxTR
EYear (mm∙year−1)	535	547	596	546	0.349	0.545	0.334
EMax (mm∙day−1)	2.60	2.70	3.08	2.76	0.125	0.517	0.229
−δ (mmol∙m−2∙s−1∙ln(kPa)−1)	81.37	79.93	77.56	73.49	0.635	0.687	0.793
GCref (mmol∙m−2∙s−1)	138.38	139.9	137.19	126.2	0.587	0.663	0.573
MaxCond (m∙s−1)	0.0065	0.0063	0.0060	0.0057	0.178	0.630	0.897
CoeffCond (mbar−1)	0.061	0.058	0.055	0.063	0.998	0.502	0.130

). Total rainfall for 2013 was 1520 mm, resulting in annual transpiration representing 35.2%, 39.2%, 36.0%, and 35.9% of annual rainfall for the C, TR, F, and TR+F treatments, respectively [39]. Monthly transpiration rates (E_Month) varied over time, with higher transpiration rates occurring during the spring and summer months when LAI and D were high (Figure 1). There was a significant interaction between month, fertilizer application and throughfall reduction (p = 0.033) such that in May, when D was at a maximum, E_Month was significantly higher in the F treatment than the C (p = 0.013), TR (p = 0.034), and TR+F (p = 0.029) treatments. There was no effect of treatment on leaf level transpiration (E_Leaf) during this month (p = 0.401), and E_Leaf averaged 0.4 mmol∙m⁻²∙s⁻¹ for May. In March and April the E_Month of the F treatment was greater than the control (p = 0.042, p = 0.029), but not the TR (p = 0.105, p = 0.063) or TR+F (p = 0.113, p = 0.079) treatments.

3.6. Canopy Conductance

The response of G_C to D was not significantly affected by any of the treatments (Table 3). The mean sensitivity of G_C to D, as defined by the parameters δ and G_Cref, was 78.1 mmol∙m⁻²∙s⁻¹∙ln(kPa)⁻¹ and 135.4 mmol∙m⁻²∙s⁻¹, respectively. This response was consistent across all treatments. The 3-PG model uses parameters describing maximum G_C (MaxCond, m∙s⁻¹) and stomatal response to D (CoeffCond, mbar⁻¹) to calculate values of G_C under optimal conditions. Estimates of MaxCond and CoeffCond were 0.006 m∙s⁻¹ and 0.059 mbar⁻¹, respectively. These parameters were not significantly affected by treatment (p = 0.522, p = 0.41).

4. Discussion

In this study, a 30% reduction in throughfall had little effect on loblolly pine water relations and growth on a poorly drained site in northern Florida after one year of treatment. Although throughfall reduction decreased the volumetric water content (θ) of the upper 20 cm of soil by 0.014 cm³∙cm⁻³, there were no differences among treatments in the θ of the 60–84 cm soil region, even when measured directly under exclusion troughs. Values of predawn water potential (Ψ_P), which in theory reflect soil water potential, did not differ among treatments and were similar to values reported for irrigated loblolly pine stands in other studies in similar soils [16,19]. The lack of Ψ_P response to throughfall exclusion despite a reduction in θ_0–20 and magnitude of Ψ_P values measured suggests that the trees at the study site were able to access soil water reserves deeper in the soil profile. It is likely that the trees at the study site were able to access the water table as numerous studies have shown that the water table typically fluctuates between 2m and the surface on somewhat poorly to poorly drained sites in the lower Coastal Plain flatwoods [40,41,42,43,44]. The rooting depth of loblolly pine has also consistently been shown to extend below 150 cm [45,46,47]. In addition to this, annual rainfall was greater than the long term average for the region during 2013 (1520 mm vs. 1404 mm) [39]. It is likely the combined effects of abundant rainfall and tree roots’ access to a high water table limited responses to throughfall reduction. Others have also noted a limited response to drought in planted pine on poorly drained Coastal Plain sites [42,48].

4.1. Whole-Tree Hydraulic Conductance

Treatments did not affect whole-tree hydraulic conductance (K_S), although K_S did vary over time. While we did not measure organ-level K_S, the literature suggests that differences in K_S between measurement dates may have been due to changes in root hydraulic conductivity which, in loblolly pine, has been shown to respond dynamically to changes in soil moisture [49]. The mean K_S across all treatments and measurements dates was 1.84 mol·m⁻²·s⁻¹·MPa⁻¹. Bartkowiak et al. [50] found that throughfall reduction in combination with fertilizer application reduced K_S from 3.83 mol·m⁻²·s⁻¹·MPa⁻¹ to 2.76 mol·m⁻²·s⁻¹·MPa⁻¹, but throughfall reduction in absence of fertilizer application had no effect on K_S. For stands of a similar age and genotype, Gonzalez-Benecke and Martin [19] reported values of K_S ranging from 0.84 mol·m⁻²·s⁻¹·MPa⁻¹ in the control treatment to 2.92 mol·m⁻²·s⁻¹·MPa⁻¹ in the irrigation treatment. The failure of throughfall reduction to affect K_S in this experiment suggests that soil water availability was not reduced to the extent that xylem cavitation was increased. This result, in addition to the lack of response of Ψ_P, suggests that water availability was not limiting in any treatment.

4.2. LAI and ANPP

Fertilizer application increased leaf area index (LAI) and aboveground net primary productivity (ANPP). Increases in LAI are generally associated with increased productivity, and a strong correlation between foliage biomass and stem wood production has been shown in loblolly pine [51]. Minimum LAI was similar across treatments prior to study establishment, but once treatments were applied the LAI of plots treated with fertilizer application was consistently higher than untreated plots with the magnitude of difference increasing through time (Table 2). Increases in LAI in response to fertilizer application have been observed in loblolly pine plantations growing on diverse sites throughout the southeastern U.S. [12,14,27,50]. In contrast to this, Ward et al. [52] found that fertilizer application had no effect on LAI one year after application in a ten-year-old loblolly pine stand growing on a well-drained silt loam soil. This lack of response was thought to be due to stand development not having progressed to the point at which tree nutrient demand exceeded nutrient supply at the study site. In the current study, plots treated with fertilizer application had significantly greater ANPP (p = 0.001) due to increased foliage production and increased woody biomass increment. Throughfall reduction had no effect on ANPP (p = 0.734).

4.3. Sap Flow and Transpiration

Sap flow per unit sapwood area (J_S) decreased toward the pith. Others have also observed reductions in J_S with increasing depth in loblolly pine, but the magnitude and pattern of this reduction varies widely among studies. Phillips et al. [53] reported that for twelve-year old loblolly pine the J_S of the 2–4 cm sapwood region was 41% of the outer 0–2 cm J_S. For mature loblolly pine, Ford et al. [54] observed that J_S was almost always greatest in the outer 2 cm of sapwood and tended to decrease with increasing depth. In contrast, Gonzalez-Benecke and Martin [19] did not observe a consistent pattern of radial variation and often observed maximum rates of J_S in the inner sapwood regions of loblolly pine. The reason for these varying results is not readily apparent. Observed reductions in J_S with increasing sapwood depth are likely the result of many contributing factors. The K_S of sapwood has been shown to decrease with increasing depth in conifers [55,56]. This decrease in K_S is likely the result of changes in the torus-margo structure that increase resistance to flow and permeability to air. Mark and Crews [57] showed that the margo of Pinus engelmannii Carr. and Pinus contorta Dougl. var. latifolia Engelm. were most open in the outer sapwood and became increasingly blocked and encrusted with increasing sapwood depth. The pit membrane structure of older xylem elements has also been shown to exhibit an increased permeability to air making them more vulnerable to cavitation [58]. The inner sapwood may also not be as well connected to the branches and leaves of the upper canopy, where transpiration rate is greatest, due to this region originally developing to supply branches and leaves that have since died or become shaded [59,60].

After accounting for radial variation in sap flow, daily maximum stand-level transpiration (E_Max) did not significantly differ among treatments and ranged from 2.6 to 3.1 mm for the C and F treatments, respectively. Values of E_Max were within the range reported for loblolly pine stands of various structures, ages, and genetics [16,19,61]. During the spring, when D and LAI were high, fertilizer application in absence of throughfall reduction increased monthly transpiration rate. The monthly transpiration rate (E_Month) of the F treatment was significantly higher than all other treatments in May, when D was at a maximum, and was also significantly greater than the C treatment during March and April (Figure 1). This observed increase in transpiration can be attributed to the higher levels of LAI in the F plots, as there was no difference in leaf-level transpiration (E_L, mmol∙m⁻²∙s⁻¹) during this time. Ewers et al. [14] also reported increased canopy transpiration (E_C) in response to fertilizer application, but when E_C was normalized by leaf area (E_L), plots receiving only fertilizer had significantly lower E_L than the control while plots receiving both fertilizer and irrigation were not different. Samuelson and Stokes [62] reported similar values of E_L for stands treated with either irrigation or irrigation and fertilizer application. These results suggest that under conditions of abundant soil moisture E_L is similar in fertilized and unfertilized stands, and observed differences in E_C are driven by differences in LAI. Under conditions of water limitation; however, fertilizer application tends to reduce E_L [19,63]. For example, Bartkowiak et al. [50] found that throughfall reduction in combination with fertilizer application reduced E_L, but fertilizer application in absence of throughfall reduction had no effect on E_L. Ewers et al. [14] hypothesized that reductions in E_L in response to fertilizer application is the result of needles developing under conditions of relative soil drought. This relative soil drought was attributed to increases in LAI without corresponding increases in root growth and increased interception losses.

Total annual transpiration (E_Year) did not differ between treatments for 2013, despite observed differences in E_Month during the spring. Values of E_Year ranged from 535 mm·year⁻¹ in the C treatment to 596 mm·year⁻¹ in the F treatment, which corresponded to 35.2% and 39.2% of total rainfall, respectively. Bartkowiak et al. [50] also found that E_Year was unaffected by throughfall reduction or fertilizer application and reported values of E_Year ranging from 577 to 652 mm·year⁻¹, which corresponded to 41%–46% of annual precipitation. Values of E_Year as high as 930 mm·year⁻¹ have been reported for irrigated loblolly pine stands receiving the equivalent of 3000 mm·year⁻¹ of annual rainfall [19]. Samuelson and Stokes [62] also observed significant increases in E_C in response to irrigation. These results suggest that the annual transpiration of loblolly pine forests is closely tied to rainfall inputs. In contrast, Ewers et al. [19] and Ewers et al. [34] did not observe increases in E_C with irrigation, but when fertilizer and irrigation treatments were applied E_C was significantly greater than if either treatment was applied individually. This suggests that under conditions of nutrient limitation additional water inputs may not impact E_C. Fertilizer application, in absence of irrigation, has also been shown to increase E_C despite reductions in E_L [14,63]. This effect is the result of increases in LAI being significant enough to increase E_C despite reduction in E_L.

4.4. Canopy Conductance

The sensitivity of G_C to D was not affected by throughfall reduction or fertilizer application. In contrast to this, Bartkowiak et al. [50] found that throughfall reduction decreased G_Cref from 171.9 to 148.1 mmol∙m⁻²∙s⁻¹ when compared to the control. This effect was even more pronounced in the throughfall reduction with fertilizer application treatment, which decreased G_Cref from 171.9 to 119.0 mmol∙m⁻²∙s⁻¹. The contrast between the results of this study and those of Bartkowiak et al. [50] may be due to differences in site conditions and/or differences in genetics, as the genetics of loblolly pine has been shown to affect G_Cref and δ [19,20]. For example, Gonzalez-Benecke and Martin [19] found that irrigation increased G_Cref and δ in a South Carolina loblolly pine seed source while a Florida seed source was unaffected. This suggested that loblolly pine originating from Florida was less conservative in regulating water loss in response to changes in soil moisture. The lack of response of the diverse northern Florida seed sources in the current study may have resulted from adaptions to abundant water supply due to a high water table and high annual precipitation.

In loblolly pine, reductions in G_Cref have been reported in response to fertilizer addition [15,63]. This reduction in G_Cref is likely due to the same factors as those hypothesized to reduce E_L [14]. Although a reduction in G_Cref was observed in response to fertilizer application in Ewers et al. [35], when both irrigation and fertilizer treatments were applied, G_Cref was found to be greater than the control or irrigation treatments alone. Samuelson et al. [16] also observed an increase in G_Cref with fertilizer application during a growing season having above average rainfall. The results of these studies suggest that the interaction between genetics, soil nutrition, and soil water availability play an important role in the response of stomatal conductance to changes in D.

The parameters describing the sensitivity of stomatal conductance to changes in D used in the 3-PG model (MaxCond and CoeffCond) were unaffected by either fertilizer application or throughfall reduction. The value of MaxCond observed in this study (0.006 m∙s⁻¹) was identical to the number reported by Bryars et al. [64] for the parameterization of 3-PG for loblolly pine. In contrast to this Gonzalez-Benecke et al. [65] reported MaxCond to be 0.0118 m∙s⁻¹ in their parametrization. This difference may be due to Gonzalez-Benecke et al. [65] using eddy covariance data to parameterize the model as eddy covariance measures total evapotranspiration and has been shown to produce higher estimates than the sap flow technique [66]. The value of CoeffCond measured in this study (0.059 mbar⁻¹), however, was similar to the value of 0.041 mbar⁻¹ reported in Gonzalez-Benecke et al. [65], but was more than twice the value reported by Bryars et al. [64] (0.025 mbar⁻¹). This result suggests that the sensitivity of G_C to changes in D is higher than is currently reflected in the parameterization of the Physiological Principles in Predicting Growth Model (3-PG) for loblolly pine [64,65]. This may have consequences, as the model is sensitive to this relationship when calculating stand transpiration and productivity.

5. Conclusions

In this study, a 30% reduction in throughfall did not significantly affect the growth or water use of a loblolly pine plantation on a poorly drained site in the lower Coastal Plain. This lack of response was likely due to abundant rainfall during the study period as well as the ability of trees to access the shallow water table at the study site. Although throughfall reduction decreased θ₂₀, water use was not affected at the canopy or leaf-level. If a 30% reduction in precipitation was applied to an entire watershed, however, it is possible that the water use and growth dynamics of loblolly pine plantations would respond differently due to a lowering of the water table. It is also likely that a 30% reduction in throughfall over a longer period of time would produce different results. Fertilizer application significantly increased basal area production, ANPP, and LAI. This response has been observed by many others and is considered to be typical [11,12,13,14,27]. Fertilizer application did not affect the sensitivity of G_C to D, but did result in increased E_Month during the spring, when D was high. During this time the monthly transpiration rate of the F treatment was 17% higher than that of the control. This increase in canopy transpiration was due to the fertilizer application only treatment having a greater LAI than the control treatment as there was no difference in E_Leaf during this time. These results suggest that a reduction in precipitation will not affect the growth or water use of loblolly pine on poorly drained sites in the lower Coastal Plain as long as the water table remains high enough to be accessed by tree roots. If reductions in precipitation are significant enough to lower the water table below the rooting zone, however, soil water resources will likely impact the productivity of loblolly pine plantations in the region. Under such a scenario, fertilizer application would likely impact the water relations of these forests due to increases in LAI, also increasing the rate at which soil water resources are depleted.