Next Article in Journal
Impact of Duration of Land Abandonment on Infiltration and Surface Runoff in Acidic Sandy Soil
Next Article in Special Issue
Seasonal Changes in Phosphorus in Soils and Vegetation of Vegetated Filter Strips in Cold Climate Agricultural Systems
Previous Article in Journal
Effect of Strip-Till and Variety on Yield and Quality of Sugar Beet against Conventional Tillage
Previous Article in Special Issue
Impacts of Low Disturbance Liquid Dairy Manure Incorporation on Alfalfa Yield and Fluxes of Ammonia, Nitrous Oxide, and Methane
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Agriculture
Volume 12
Issue 2
10.3390/agriculture12020167
agriculture-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Tile Drainage Flow Partitioning and Phosphorus Export in Vermont USA

by
Ryan Ruggiero
1,
Donald Ross
1 and
Joshua W. Faulkner
1,2,3,4,*
1
Department of Plant and Soil Science, University of Vermont, Burlington, VT 05405, USA
2
Extension Center for Sustainable Agriculture, University of Vermont, Burlington, VT 05405, USA
3
Department of Civil and Environmental Engineering, University of Vermont, Burlington, VT 05405, USA
4
Gund Institute for the Environment, University of Vermont, Burlington, VT 05045, USA
*
Author to whom correspondence should be addressed.
Agriculture 2022, 12(2), 167; https://doi.org/10.3390/agriculture12020167
Submission received: 4 December 2021 / Revised: 11 January 2022 / Accepted: 23 January 2022 / Published: 25 January 2022
(This article belongs to the Special Issue Optimizing Nutrient Management in Cold Climate Agroecosystems)
Browse Figures
Review Reports Versions Notes

Abstract

:
Tile drainage (TD) has been identified as a potential non-point source of phosphorus (P) pollution and subsequent water quality issues. Three fields with TD in Vermont USA were monitored to characterize hydrology and P export. Fields were in corn silage and used minimal tillage and cover cropping practices. Preferential flow path (PFP) activity was explored by separating TD flow into flow pathway and source connectivity components using two hydrograph separation techniques, electrical conductivity end member unmixing, and hydrograph recession analysis. TD was the dominant P export pathway because of higher total discharge. Drought conditions during this study limited surface runoff, and possibly resulted in maximum PFP activity in the active clay soils. The non-growing season dominated annual P loading for two of the three study years. Peak P concentrations in TD occurred during events following manure injection in the fall, as well as in the spring post cover crop termination and post-planting. Intra-event analysis of rainfall pulses showed that TD flow and P concentrations were higher because of higher intensity pulses. This study highlights the impacts of current manure management, as well as the potential for climate change to increase P transport to TD.

1. Introduction

Reducing nutrient loss and subsequent water degradation is a challenge for agriculture as we explore the boundaries of crop and livestock yields [1]. Phosphorus (P), among other required nutrients, is applied to farmland to increase fertility; it is easily transported from the soil to runoff and eventually surface waters, resulting in non-point source P pollution [2]. Accumulated legacy P from long-term application of P fertilizers and manure increases the difficulty of managing agricultural runoff [3,4]. Despite this, dairy cropping systems continue to apply manure P since the land application is the primary economically viable method for its disposal, and nitrogen contained in the manure is needed for crop production [5]. In addition, in the face of climate change, where rainfall has and is expected to continue to increase during the period leading up to planting, tile drainage (TD) has increased as a solution to concerns around spring field accessibility and crop yields [6,7,8]. TD alters the hydrologic processes that control P transport during storm events, and in some instances is regarded as a best management practice (BMP) for reducing P in agricultural runoff [9,10,11,12]. However, TD installed in fields with poorly drained soils, consisting of finely textured clays that are prone to desiccation cracking, usually embody preferential flow pathways (PFP) [13]. PFP plays an important role in TD hydrology and nutrient export and has been highlighted as a significant transport mechanism for P [14]. PFP permits rapid movement of water, reducing resorption of P to the soil matrix, which results in high P concentrations in TD [15,16].
Previous studies have used two different hydrograph separation methods for estimating PFP contributions to TD, and until recently the two have yet to be combined to clarify matrix-PFP interactions and their impact on P transport [17]. The first method, hydrograph recession analysis (HRA), separates matrix and PFP contributions into slow flow (SF) and quick flow (QF), respectively, by examining the hydraulics of the system. Here, the soil profile is assumed to drain via multiple reservoirs, thus PFP is separated from matrix contributions based on groupings of pore sizes that drain at similar rates [18,19]. The second method, a mass balance or end-member mixing analysis (EMMA), separates matrix and PFP contributions into old and new water, respectively, via a tracer that accumulates in infiltration as contact with soil media increases [20]. This method estimates matrix and PFP contributions based on contact time with the bulk soil [21,22,23]. By combining the two methods into a four-component hydrograph model, the mixing of infiltrating water between the matrix and PFP can be assessed, thus improving estimates of P export to TD [17].
Increasing temporal and spatial contact between P and the soil matrix can result in lower P concentrations in runoff [24]. This is achieved by tillage and fertilizer incorporation, which also have the effect of disrupting PFP, thus reducing P loss to TD [23,25,26]. However, mechanical action is known to reduce soil structure and aggregate stability that can increase sediment in runoff [27], which is easily transported to TD via PFP [13]. Manure injection (MI) is a BMP that incorporates manure P with the soil matrix with minimal aggregate disturbance [28]. MI can reduce P losses in surface runoff [29], however, it is still unclear the effect of MI on P losses via TD [24,30].
While P concentrations in agricultural runoff correspond to P levels in the soil [31], peak P loading occurs during rainfall and snowmelt events because hydrology is usually the controlling factor of P export in agricultural watersheds [14]. P transport dynamics are also a function of season and antecedent moisture, which allude to the biogeochemical process controlling P solubility and export [14,15]. In the Lake Champlain Basin (LCB), 38% of P loading to the lake has been attributed to agriculture [32], and rainfall and climate trends suggest more extreme periods of wetness and drought [8]. Thus, climate change may have a direct impact on agricultural P export [33].
This study aims to improve the understanding of P transport to TD in the LCB, and to help adapt nutrient management decisions to the impacts of climate change. Extensive TD research in the USA Midwest provided an analytical foundation for characterizing P transport in the Northeast, where intensive dairy cropping agroecosystems are also prevalent, yet there is a lack of edge of field (EoF) research [34]. Year-round monitoring of rainfall and surface and TD discharge metrics was performed to (1) characterize P transport from TD in Vermont, USA (VT) and (2) assess the role of PFP and rainfall dynamics on P transport. Effects of antecedent moisture condition (AMC), seasonality, P application timing, and transport/source limitation behavior were explored to characterize P transport dynamics [35,36,37]. A dataset of high temporal resolution measurements of TD flow, TD P concentrations, and rainfall from an archetypal VT dairy agroecosystem was used for this analysis. Also, the four-component hydrograph separation model proposed by Nazari et al. [17] was used to assess PFP activity, as well as matrix-PFP interactions. A unique intra-event rainfall pulse analysis was performed using temporal rainfall data, where it was hypothesized higher intensity rainfall pulses would result in higher P concentrations in TD because of PFP activity.

2. Materials and Methods

2.1. Site Description

The study area was chosen because Lake Champlain contains eutrophic areas, and agriculture has been identified as a significant P contributor in the watershed [7]. Other agricultural watersheds have similar geographies, land management practices, and water quality issues, which makes this study valuable for management decisions elsewhere. For example, TD is common in the US Midwest and has been identified as a major P export pathway to the Great Lakes [14]. Furthermore, fields with TD in this study and the US Midwest have similar soils (high clay contents) [11], both respond rapidly to rainfall indicating highly-active PFP networks [21], and both have the potential to export the majority of annual P via TD [14].
TD was measured year-round from intensive dairy forage production fields, located within the western LCB (Addison County, VT, USA; Figure S1. Three different TD networks draining two separate field sites were monitored. The first site is drained by the TD network ‘F1’ and the second field site, ‘DC’, is drained by the two separate TD networks of DC south (DCS) and DC north (DCN). The TD network at F1 was installed in 2016 and both networks at DC were installed in 2018. TD network areas were determined from installation maps provided by the installer; F1 is 14.16 ha, DCS is 8.0 ha and DCN is 4.85 ha. At both sites, TD was installed to a depth of 1 m and is spaced 7.62 m (25 ft) apart. TD laterals at F1 are 101.6 mm (4 in) in diameter and are connected to a 203.2 mm (8 in) TD main, and at DC the TD laterals are a diameter of 101.6 mm and are connected to 152.4 mm (6 in) TD mains. The sites are 3.2 km apart and the mean annual rainfall and temperature for the area are 94 cm and 7.8 °C, respectively. Site soils were not identical but are dominated by fine-textured soils that are prone to desiccation cracking. F1 is of the Vergennes clay soil series, while soils at DC are dominated by the Covington and Panton silty clay series with a small vein of the Swanton fine sandy loam series in the east/upper parts of the field [38]. Soil P was tested late in the growing season of 2021 and levels were in the low-optimal range for field crops (6 mg/kg Modified Morgans).
Both sites were in corn silage production during the study period and since TD installation occurred, yet, F1 was in hay production before TD installation. Before TD was installed, the farmer had formed multiple broad swales extending across the entire field that was sloped to the edges to improve surface drainage. This grading occurred at DCS, and the southern portion of DCN, but not at F1. These surface features remain, and the broad swales drain to surface inlets that have individual outlets and are not linked to the TD network. During the study period and in the few years prior, both sites received light chisel tillage before planting and dairy manure was injected after the corn harvest (Table 1). In the fall of 2020, the farmer performed deep tillage along the field topographic contours to effectively create surface roughness features (i.e., ‘water bars’) throughout the field that were intended to intercept surface runoff and promote infiltration.

2.2. Field Measurements

TD flow and nutrient data were collected at F1 for the entirety of the 2019 and 2020 water years (WY) and three events in the 2021 WY. At DC, data were collected for the entirety of the 2020 and 2021 WY and throughout October of the 2022 WY. Rainfall was measured using both a tipping bucket rain gauge (Onset Computer Corp., Bourne, MA, USA) and a manual rain gauge at DC. Rainfall was estimated for the 2019 WY at F1 using NOAA NOWData for the Burlington, VT area [39].
TD flow was monitored at the EoF where TD mains were discharged. Flow from the TD outlet at F1 was directed into an H-flume (1.5 ft) and was measured using the flume along with a compound weir (Thel-Mar, LLC, Dickson, TN, USA) inserted into the pipe. Stage in the weir at the TD outlet was measured using a bubbler flow module (Teledyne ISCO 730 Bubbler Flow Module) and converted to flow rate using rating curves. Below the full capacity of the weir of 10.7 L/s, the weir rating curve was used to determine the flow. At or above this flow rate, a relationship between the water pressure at the TD outlet and the stage in the flume was used to determine the flow.
TD flow at both DCS and DCN was measured using in-line electromagnetic flowmeters (ModMag M1000, Badger Meter, Inc., Collegeville, PA, USA). The flow meters were set back several meters from the TD main outlet and installed below ground. Manholes were dug to access the TD main and the main was then cut and the flow meter placed at the end of the new outlet (Figure S2). The end of the flow meter was allowed to drain freely into the manhole, where a standpipe was installed on the inlet of the old TD outlet to control the water level in the manhole. The electromagnetic flowmeters required pipe-full conditions during measurements, thus the outlet elevation using the standpipe was set just above the top of the flow meter outlet.
Automatic water samplers (6712, Teledyne ISCO, Lincoln, NE, USA) were used to record stage or flow at five-minute intervals at each TD outlet, as well as sample the discharging TD water. At F1 samples were taken using an anchored intake just before the outlet and at DC the sample tubing was anchored at the outlet of the flow meter. Samplers contained 24 1-L bottles and a two-part program was typically used; the first set of bottles was reserved for composite baseflow sampling (part A) and the second set was reserved for discrete event flow sampling (part B). Baseflow sampling was time-based and was disabled during events. Event sampling was triggered based on the rise in flow rate, and rather than 1 h composite samples throughout the hydrograph [13,22], event sampling was discrete and both time and flow proportional sampling was used to maximize resolution during peak flows. Event pacing was determined based on prior observations of the network response and flow characteristics, while baseflow pacing occurred every six hours and used composite sampling. Minor modifications to the automatic sampler programming occurred throughout the study period to account for storm size, baseflow, etc. yet the general approach remained consistent.
Surface runoff was measured from distinct surface watersheds above the DCS and DCN TD networks as part of a paired watershed study in the calibration phase occurring at the DC site (Figure S3). Surface runoff was measured using automatic water samplers with methods like those at the tiles except for the flow measurement device. At DCS surface runoff was collected by a surface inlet set back several meters from the EoF and directed underground through a 0.3 m (12 in) pipe to the EoF, where the flow was measured using a bubbler module to determine the stage in a compound weir inserted into the pipe. At DCN surface runoff was directed using wooden wing-walls into a 0.47 m (1.5 ft) H-flume located at the EoF where a bubbler module was used to determine the stage in the flume (Figure S2). At DCS water samples were collected using an intake anchored in the outlet pipe and at DCN the intake was placed in a plastic box anchored at the flume outlet.

2.3. Water Quality Analysis

Events were anticipated from the weather forecast, however, a wireless modem was used to communicate with the automatic samplers to remotely determine if event sampling had occurred. Event samples were retrieved from sites and transported to the lab within 24 h, and grab samples were returned to the lab the same day. Samples were analyzed at the University of Vermont Agriculture and Environmental Testing Laboratory located in Burlington VT by standard methods for total P (TP) (SM 4500-P F: alkaline persulfate digestion and flow injection analysis) and soluble reactive P (SRP). Sample splits for SRP were filtered using a 0.45 um membrane filter and frozen until the analysis. Most event samples in this study were turbid and thus centrifuged until they were non-turbid and then decanted to ease with filtering. The P fraction of TP-SRP is equal to particulate P (PP) plus dissolved unreactive P, where dissolved unreactive P is organic and PP smaller than the filter pore [40]. However for simplicity, herein it is assumed that TP-SRP is PP [13,36,41]. Runoff sample and rainwater (collected from the manual rain gauge) EC was measured using a benchtop EC meter (Amber Science, Inc., Eugene, OR, USA).

2.4. Analytical Methodology

Runoff event initiation points corresponded with the initial rise in TD flow, usually from low- or no-flow conditions. Events during the dry season usually ended when TD flow returned to zero, while events during the wet season usually had long recession limbs. For these events, the last visible inflection point on the recession limb was used to determine the end of the event hydrograph [42]. Event volumes were calculated as the area under the hydrograph, time to peak was calculated as the time between the start of the event and the maximum flow rate for the event, runoff ratio (RR) was calculated as the ratio of total event runoff depth to the total rainfall depth, and hydrograph response time was calculated as the time difference between the start of rainfall and the initial hydrograph rise [22]. Event total rainfall was determined as the rainfall in the 24 h leading up to the event start time plus rainfall during the event.
Chemographs were constructed by linearly interpolating between sample bottle concentrations to achieve a continuous concentration dataset for each event at the resolution of the flow data (i.e., 5 min). Starting and ending concentrations were assumed to be entirely SRP and set to the average concentrations in baseflow samples at each site during the study period. Loadographs were constructed by multiplying the continuous flow (i.e., hydrograph) and nutrient concentration (i.e., chemograph) datasets. Event loads of TP, SRP, and PP (TP-SRP) were determined by integrating the loadographs, and event flow weighted mean concentrations (FWMC) were back-calculated by dividing the mass exports by the total event volume. At F1, for six relatively small to moderately sized events, insufficient samples were obtained to calculate P export as described above. For these, an average P concentration for the event was assumed based on the samples available and concentrations in other events temporally nearby. Also at F1, for three large events, equipment failure resulted in both incomplete hydrographs and incomplete TP samples. For these, export was estimated from either another nearby monitored tile or from previous storms of similar rainfall intensity. At DC, loadographs were constructed for events that fell between 12 May 2020 and 21 July 2021, while compositing sampling was used outside of this period [43]. At the DC site, P export during missed events was estimated using regression relationships with nearby tiles. Baseflow loading was simplified because sampling was intermittent, and P concentrations were generally very low. Thus, daily baseflow P loading was set to a constant site-specific value, where it was assumed this loading rate occurred every day there was TD flow but no event hydrograph.

2.4.1. Four Component Hydrograph Separation

This study closely follows the methods used in Nazari et al. [17] to combine the two hydrograph separations, EMMA and HRA, into one four-component model. The methods are presented in detail in that study, and modifications to them are described here. First, their approach uses 30-min data while here we used 5-min data, which maximized the resolution of the flow data due to significantly shorter recessions observed in this study. Our results suggested that daily aggregated data of flow and P concentrations would result in under-estimating P loading due to flashy TD hydrographs. We also found that there was high multicollinearity between flow components in the four-component hydrograph separation model. Thus, we chose not to pursue the multiple linear regression analysis. We performed the master recession curve analysis in RC 4.0 using the matching strip method [44]. Results from the master recession curve analysis showed that QF:SF recession constant ratios were greater than 3, thus suggesting that reservoir distinctions were valid (Table S1) [13,17]. Concerning the individual event hydrograph separations, intermediate reservoirs [45,46] of TD hydrograph recessions in log-normal space were frequently observed, making it difficult to isolate a single inflection point to represent the peak of SF [17]. Thus, to separate QF and SF for each event, the recession constant of the shallowest of the observed linear reservoirs was used as the slope of the linear increase in slow flow from the start of the event. The intersection of this straight line with the recession of the actual hydrograph represented the end of QF and the subsequent hydrograph was set to SF.
An additional modification from the Nazari et al. [17] methods was that continuous EC data were not available in this study because EC was measured in discrete samples (i.e., individual event and baseflow bottles) [47]. Thus, only the events where baseflow EC was available prior to the event were used in the source contribution analysis and subsequently the four-component hydrograph separation. This was determined since the rapid dilution of EC was observed between baseflow samples and the event’s first discrete sample, even if event sampling occurred prior to TD flow reaching 1.0 L/s. To perform the EMMA with the discrete EC values, EC was linearly interpolated between sample values, and the value for baseflow EC prior to the event was used as starting and ending values for the matrix signature (i.e., new water = 0 at the start and end of the event). In addition, rainwater EC was consistently several orders of magnitude less than surface runoff and surface runoff was not available for most events. Thus, the minimum value of surface runoff EC was used across all events as the signature for new water. This assumption represents the minimum amount of mixing of surface runoff with the layer of interaction [48] and is the most accurate signature of the new water end-member available in this study.

2.4.2. Rainfall Pulse Analysis

TD readily and rapidly responded to rainfall, thus, rainfall pulses were analyzed to better understand how pulse intensity affects TD discharge metrics associated with the timing of a rainfall pulse. Pulse data were determined from the temporal rainfall data from the tipping bucket rain gauge. The tipping bucket rain gauge records the duration of time between rainfall increments of 0.254 mm (0.01 in). The tipping bucket rain gauge data was aggregated on an hourly basis by clock hour to obtain rainfall pulses. Consecutive 30-min pulses greater than zero were lumped together into a single pulse, and pulse metrics of total volume, maximum intensity, and duration were calculated. The maximum pulse intensity of the pulse group was set equal to the highest 30-min intensity contained within the group. The maximum intensity was compared to the quartiles of the period of record of the 30-min pulse data to assign levels to pulse group (herein ‘pulse’) intensities. Pulses that fell below Q4 were assigned to ‘Low’ and pulses above Q4 was set to ‘High’. The response variable window for a pulse was the hour following the start of a pulse to the hour following the end of a pulse. The mean of the three P species loadographs, along with the maximum of the hydrograph, in the response variable window was calculated as the response variables associated with the pulse.
A two-sample t-test was used to compare means from the two groups of pulse intensities. This analysis was essentially an experiment where rainfall pulses are treatments, and the response variables are TD mean P loading rates and maximum flow rate in the time window associated with the rainfall pulse. For the assumption of independence to be met between replicates, we must assume that the P pool is inexhaustible, thus there is an equal chance for the effect of pulse intensity to generate a given TD P loading rate regardless of how much was exported previously. This was not the case for events close to MI since data suggested old water contained high soluble P prior to events, suggesting source limitation post-MI. Osterholz et al. [36] suggested removing observations six months post P application to isolate the effects of legacy P on TD concentrations. However, because of the drought and frequent P applications, very few events and rainfall pulses met this assumption in this study. As a result, we removed pulses that occurred within 100 days of MI to minimize the effect of source limitation P transport.

2.4.3. Statistical Analysis

The Kendall’s rank correlation coefficient was used to determine the relationship between event rainfall and TD discharge metrics and time series, namely day of the hydrologic year, number of days since MI, and number of days since P application, which included MI, cover crop termination, and fertilizer application during planting (Table 1) [35,36]. Pearson correlations were used to correlate TD discharge and rainfall metrics. Event rainfall metrics included total rainfall, which was determined as the rainfall in the 24 h leading up to the event plus the rainfall during the hydrograph, max hourly rainfall intensity, and 24-h and 7-day and 30-day rainfall totals. Event TD discharge metrics included P species loads and FWMC, total water discharge, peak flow rate, response time, time to peak, and RR. The Kruskal-Wallis nonparametric test was used to evaluate significant differences when data was not normally distributed. Factor groups with significant Kruskal-Wallis p-values were compared using the Dunn test post-hoc analysis to determine significant differences between individual levels [49]. Analysis was performed using R software [50].

3. Results

3.1. General Hydrology and P Transport

Rainfall was 108, 63, and 70 cm during 2019, 2020, and 2021 WY, respectively. Event rainfall metrics, namely total rainfall, maximum hourly rainfall intensity, and 24 h and 7- and 30-day antecedent rainfall, were positively correlated to the day of the hydrologic year (Figure S4). For maximum hourly rainfall intensity, this is the expected result because of high-intensity summer thunderstorms, however, correlations were highly influenced by a very wet July in the 2021 WY (Figure S4). An abnormally dry to moderate drought period occurred from June 2020 to August 2021, and during this period, TD regularly responded to rainfall without surface runoff occurring [51]. A lack of surface runoff was also attributed to high surface roughness from subsoiling in the fall of 2020. TP and SRP FWMC in this study and elsewhere in the LCB [13] were higher in the surface runoff than in TD. However, event peak TP concentrations in TD flow were found to exceed surface runoff TP FWMC on some occasions. Findings in other cold climate regions have found that the majority of surface P export occurs during the non-growing season because the majority of annual flow may stem from snowmelt [52,53]. Results in this study showed that snowmelt did not result in higher surface runoff volumes relative to summer thunderstorms. Also, P concentrations in surface runoff and TD were higher in the summer than during the snowmelt event, thus the non-growing season was not a significant P export period via surface runoff. TD was the dominant P export pathway, and results coincide with the notion that hydrology dictates P transport via TD [14]. Most of the annual TP export through TD was during event flows, where larger events generally exported more P (Table S2).
P export was highly variable between events and years. In fall 2020 there was no baseflow prior to and after event hydrographs because of the drought, suggesting the entirety of the event flows and P was transported to TD via PFP. Also, RR for these events was among the lowest and hydrograph response times were relatively high (above the interquartile range; Table S2), suggesting that a soil moisture deficit needed to be overcome to allow PFP activity [54,55]. In the 2019 and 2020 WY, most of the annual TP loading through TD was during the non-growing season, while in the 2021 WY, TP loading was similar between the growing and nongrowing seasons (Table 2). This stems from the drought that limited the number of events during the growing season in WY 2020 and during the non-growing season in WY 2021. Despite this, average event volumes and RR were still highest during the non-growing season across all sites and WY (Table 2).

3.2. Rainfall, TD Discharge, and P Export Dynamics

3.2.1. Seasonal Differences

In the summer, the event maximum rainfall intensity was the highest, and time to peak flow rate and RR were significantly lower than the other seasons (Figure 1). Event TP and SRP FWMC in TD were not significantly different between the four seasons, however, SRP loads were significantly higher in the winter (Figure 1). King et al. [56] found that while TD SRP FWMC was lowest in the winter, SRP loads were significantly higher than summer and fall values. The winter represents a period when legacy P acts as a homogeneous P source during events and thus P is more transport limited [36]. SRP FWMC may be lower in winter because it’s the longest time since P application, thus, legacy P is the main contributing source of P to TD [36]. However, since P-mineral soil binding energies decrease due to wetter soils in winter, SRP loading may be higher [14].
The mean event TD RR was 44%; F1 and DCN had mean RR of 55% and 54% respectively, while DCS had a mean of 26% (Table S3). There were 5 events where RR was greater than one, all of which corresponded to fall rain-on-snow events. Mean RR at DCN during the NG season of the 2020 WY was 94% (Table 2), suggesting that shallow groundwater was significantly augmenting event flows at this site during this period. At the field scale with extremely fine-textured clay soils, relatively short response times to rainfall and flashy hydrographs are usually attributed to PFP [13,22]. During the NG season soils stay wetter for longer and may be closer to field capacity. Thus, clay soils would be expected to swell and the PFP network would decrease [57]. However, rapid TD response to rainfall occurred throughout this study and over a large range of AMC. This suggests that PFP was highly active regardless of AMC, which agrees with others who have found that PFP activity persists during the non-growing season in fine-textured soils [17,21,58]. PFP transport to TD in this study could also be a function of the recent TD installations (2016 at F1 and 2018 at DC). Tiles backfilled incorrectly or backfilled during a drought period when soils are dry and blocky may result in large PFP that exist directly over TD laterals, which may not diminish until years of freeze-thaw cycles help settle the backfill material [59].

3.2.2. P Export Relative to P Application Periods

Days since P application was the better predictor of event TP and SRP loads, peak sample concentrations, and FWMC, as opposed to days since MI (Figure S5). Nevertheless, P concentrations in event samples and event FWMC were highest post-MI. There were four MI periods observed during this study, i.e., the autumns of each year. There was moderate and severe drought pre-MI in the falls of the 2019 and 2021 WY, respectively, while conditions pre-MI in the 2020 and 2022 WY were abnormally dry and normal, respectively. In the falls of the drier years, event TP export through TD was relatively low because of small FWMC (<1000 ugP/L) and small event discharges in the 2019 and 2021 WY, respectively. In the falls of the wetter years, event FWMC remained elevated throughout the fall (>1000 ug P/L), whereas in the 2020 WY, the result was that fall events contributed to most of the annual P export. P concentrations in TD were also elevated following the P application periods later in the WY, namely, post cover crop termination and during planting. In 2019, 2020, and 2021 WY, P concentrations in event samples post cover crop termination were elevated relative to events prior. As the glyphosate-P application rates were relatively small (Table 1), this suggests that the cover crop was the main contributing P-source during this period, which can occur during the spring thaw if their cells lyse over winter [60]. In the 2019 WY, the only year that had events post-planting and post fertilizer application, event peak P concentrations spiked slightly higher post fertilizer application than post cover crop termination.
The relationship between P transport and AMC in this study is contrasted to other studies that have shown that event P transport to TD on clay soils is reduced during wetter AMC due to the swelling of desiccation cracks [61]. While the PFP network may have been reduced during wetter AMC via the shrink-swell nature of the soils, here, we were unable to show evidence of this. This can partly be attributed to the coarse nature of measurements, i.e., we did not examine the soil profile. Nevertheless, throughout our study, TD hydrographs responded rapidly to rainfall and had flashy rising and falling limbs suggesting that most event flows were transported through highly active PFP regardless of the AMC. Granted we believe that the PFP network represented a maximum condition during this study period because of several factors, including prolonged drought, recent TD installations, and subsoiling.
Because of the relationships between P export and AMC, and that P export decreased as the time since P application increased, our findings suggest that TD P loads could be reduced if MI/P application is timed. In some crops, the ideal time might be during the growing season, when plant P uptake and soil vegetative cover are at a maximum, but this is not possible in many annual crops due to the inability to access the field without damaging the crop. Another ideal time for MI/P application would be when PFP is minimized, and matrix flow dominates TD. This could be during the non-growing season when the clay soils swell and PFPs close, however in cold climates manure application is often prohibited during the winter months to reduce surface P losses from manure application on frozen or saturated soils.

3.2.3. P Concentrations in TD during Events

There were 840 TD event samples analyzed for TP, and 467 were also analyzed for SRP. Overall, while linear trends were significant for TP, flow rate explained little of the variation in P concentrations during events (R2 = 0.05). This agrees with others who have concluded that flow rate poorly predicts P concentrations in TD due to high variability at low flow rates [17,47]. In this study, the highest TP and SRP concentrations were post-MI and coincided with initial event samples. When removing samples that occurred within the first hour of the event, the relationship between SRP concentrations and the flow rate was significant, and when removing samples taken within the first 10 h, linear trends drastically improved (Figure 2).
Post MI, the highest event P concentrations in TD corresponded to the initial hydrograph response, while in the spring and summer, the highest event sample concentrations corresponded to the highest flow rates (Figure S6). While this suggests P transport was source limited in the fall (post-MI) and transport limited in the spring (post cover crop termination/fertilizer application), events during the fall also showed transport limitation behavior apart from the samples corresponding to the initial hydrograph rise. The peak of the drought during this study occurred in the fall of 2020, where one large event on 30 September occurred 3 days post-MI, however, a large soil moisture deficit limited TD flow. There was no TD baseflow after this event and the next event produced a TD response, which was on 12/01. Because of a lack of baseflow, we hypothesize that manure P was still abundant in the surface soils prior to the events in December. Two events in succession in December 2020 are discussed below to highlight the role of AMC on SRP transport under macropore flow. Events on 1 December 2020 and 25 December 2020 occurred during the drought and there was no baseflow prior to the event. However, soil AMC was much higher for the second event due to snow on-field prior. Under this condition, peak event TP and SRP concentrations were 3–4 times higher than the peak concentrations on 1 December 2020 (Figure S6).
The use of MI in this study resulted in a highly labile source of P near the surface, which resulted in high SRP concentrations in TD when conditions were ideal for PFP transport of P-laden pre-event water. TD hydrographs may be initiated via unsaturated macropore flow, which originates at the surface and is usually associated with little to no matrix-PFP mixing. As the event continues, flow transitions to saturated macropore flow that transports both event and pre-event (i.e., matrix-displaced) water [23]. Matrix-PFP interactions in soils with high labile P (e.g., soon after P application and/or high legacy P concentrations) may contribute significantly to SRP concentrations in TD during events under saturated macropore flow because a higher proportion of pre-event-like water is transported [23,26]. However here and in Williams et al. [62], peak event SRP concentrations at the beginning of events suggest that unsaturated macropore flow has the potential to be the dominant source of TD SRP concentrations if soil AMC is high.

3.3. Flow Pathway and Source Connectivity

3.3.1. QF, SF, New and Old Water Metrics

Differences in the timing of the two hydrograph separations suggest that flow pathways were not equivalent to source contributions [17]. On average, the peak of old water occurred over 20 h prior to SF (Table 3). Also, the average timing of peak QF contributions occurred over one hour prior to peak contributions from new water, however, there were five events at DCS and ten events at DCN where QF and new water timing were equivalent (Table S4). QF and new water, as well as SF and old water, had different cumulative volumes, further supporting the notion that flow pathway and source contribution estimates are not equivalent in TD [17].
Mean contributions of QF and new water were 78% and 48% of total event TD flow, respectively (Table 3). Peak flow contributions of QF and new water were 96% and 80%, respectively. While Nazari et al. [17] reported lower peak flow contributions for QF (and similarly for new water) than in this study, mean new water (and QF) contributions were similar to Nazari et al. [17] and others also using EC-EMMA hydrograph separations [22,63]. Also, we observed similar but slightly lower mean RR across all sites as compared to Nazari et al. [17], whose study site was roughly twice the size of DCS and four times the size of DCN. Here, total, QF, and new water RR were consistently higher at DCN, the smaller of the two DC fields, suggesting field size and thus TD area may have less of an impact on event volumes and pathway/source contributions. Differences between RR may be explained by the fact that DCN has a greater topographic relief than DCS. Surface ponding at DCN will have a higher chance of finding the soil cracks because it is sloped, and in addition, the deep tillage may have resulted in substantially more PFP, which led to more QF and new water at DCN.

3.3.2. Four Component Hydrograph Separation

Correlation analysis suggests that QF, new water, and QFnew were the superior predictors of event P loads (Table S5). Previous studies have suggested that PFP control SRP transport to TD because of the strong correlation between concentrations of water-extractable P in the layer of interaction and SRP in TD [17]. Nazari et al. [17] showed that including the four-component hydrograph separation improved the prediction of SRP concentrations in TD, suggesting matrix-PFP mixing is an important indicator of SRP in TD. Thus, the four-component hydrograph separation approach could be a valuable tool for the LCB. For example, matrix-PFP interactions could improve field scale-P export models like the VT P-Index, a tool used by farmers and agricultural service providers that estimates the risk of P loss from fields based upon P source and transport factors [64]. Here, loads of all three P species were similarly correlated to the PFP hydrograph components (Table S5), suggesting that PFP activity also controls PP transport to TD [13]. PP may constitute the majority of TD TP as was shown here and elsewhere in the LCB [65], as the dispersive-like clay soils in this region are easily suspended in infiltrating waters, and sorbed PP is then transported to TD via PFP.
On average, QFnew represented a high proportion of total flow (46%), and peak contributions were as follows: 77% for QFnew, 69% for QFold, 12% for SFnew, and 41% for SFold (Table 3). Nazari et al. [17] found lower peak contributions from QFnew (66%) and QFold (27%), yet higher peak contributions from SFnew (33%) and from SFold (98%). This could be because of the differences in data collection and analysis methods, as well as differences in soils, tillage, and slope. For example, conservation tillage was employed in the Nazari et al. [17] study while here deep tillage (i.e., subsoiling) was performed to promote infiltration.
QFnew had higher mean volumes, yet QFold volumes were equal to or exceeded that of QFnew for eight of the smaller 28 site events in the analysis (Table S4). An example of QFold dominating the entirety of event flows was during back-to-back events of 2 July and 3 July of 2021. For the duration of a multi-peak event on 2 July 2021, QFnew remained higher than QFold. The event the following day, which occurred on saturated soils, QFold remained higher than QFnew. This suggests that in back-to-back events AMC can have an impact on matrix-PFP mixing. However, AMC metrics in this study, namely 24-h and 7- and 30-day rainfall totals, were not positively correlated with QFold volumes (Table S5). While higher AMC corresponds to wetter macropore walls, where the potential for matrix-PFP interactions increases [23,26,61], results here suggest that soils needed to be above field capacity for QFold to be the dominant transport pathway to TD at these sites. We observed similar correlations between QFold and QFnew volumes (r = 0.49) as in the Nazari et al. [17] study (r = 0.40), who also found that only the smaller events had QFold volumes exceeding QFnew (3 of 27 events). QFold had a much smaller standard deviation than QFnew, suggesting QFold has the potential to be the dominant transport pathway during smaller events simply because there is less QFnew.
Here and in Nazari et al. [17], AMC and event rainfall were not good predictors of QFold and the other three hydrograph components volumes. We found negative correlations between 7-day rainfall and QFnew and SFnew, and between maximum rainfall intensity and SFold and SFnew (Table S5). Also, we observed positive correlations between time to peak and the volumes for each of the four components except for QFold. Since QFold time to peak was correlated to the overall hydrograph time to peak, it is expected QFold volumes would be as well. QFold was the only component not correlated to the above-mentioned rainfall and TD discharge metrics, and had the lowest coefficient of variation, suggesting this pathway is controlled by different mechanisms relative to the other three components.
The peak of SFnew occurred at an average time of 27.4 h into the event, which is slightly quicker than in Nazari et al. [17] (32 h). Nazari et al. [17] found a negative relationship between the timing of SFnew and 10-day rainfall, while here we found a semi-strong negative relationship between SFnew with 24 h and 30-day rainfall, and a very strong negative relationship with maximum rainfall intensity (Table S5). The SFnew pathway represents event water transported to deeper horizons via PFP that absconds back into the matrix prior to the tile [17]. This pathway represents groundwater recharge via PFP, and frequent recharge of groundwater via PFP suggests P from the surface could be transported and accumulate as legacy P in much deeper soil horizons [66].

3.4. Rainfall Pulse Analysis

There was an uneven number of pulses identified (Table S6), and pulse data were non-normal, thus TD discharge metrics associated with rainfall pulses were compared using a non-parametric comparison of means. P loading rates in TD during periods associated with high-intensity rainfall were significantly greater than during periods associated with low-intensity rainfall (Figure 3). Also, high-intensity rainfall pulses resulted in significantly higher flow rates as compared to the lower intensity rainfall pulses. Others have found relationships between rainfall metrics and P concentrations and loading in TD at the inter-event scale, however, this is the first known study to examine these relationships at the intra-event scale. Vidon and Cuadra [45] found that P concentrations were more variable during larger rainfall events as compared to smaller ones, where we found similar results with higher variation associated with higher pulse intensities.
Despite soil P levels being relatively low for agricultural fields that are considered hot-spots for legacy-P export (6 mg/kg Modified Morgans) [31,36,67], here, annual TD P loading was still relatively high for the LCB [10,68]. Because of several factors during this study, the PFP networks in these fields possibly represented a maximum condition (see Section 3.2.2). As time proceeds and/or under different AMC, the PFP network may diminish, resulting in less rapid infiltration and possibly different rainfall-driven P transport dynamics. This may also shift the P export-dynamic between the surface and subsurface, where the surface may have higher runoff volumes and thus P export, since less water is being channeled to TD via PFP. Moreover, climate change is expected to continue to increase the intensity of both drought and rainfall events in the LCB [8], and results here suggest this has implications for P transport via TD. In the LCB, P in the lake bottom sediment is easily mobilized to the overlying water column, thus PP exported at any time throughout the year may contribute to cyanobacteria blooms [69]. Legacy P will continue to be a challenge for this region, which persists in the watercourse, acting as a source of P even after BMP reduces P export from agricultural fields [4].

4. Conclusions

TD was the dominant P export pathway at two sites, and annual and seasonal loads varied and were affected by the drought. The results of the two- and four-component hydrograph separations aligned with previous TD studies. Here, PFP was the dominant transport pathway during events. Results agree with the notion that annual P export from TD has the potential to be the greatest during the non-growing season, since it is the time of year with the highest TD discharge [14]. Spikes in TD P concentrations during events coincided with periods of P influxes, namely MI, herbicide application/cover crop termination, and fertilizer application at planting.
MI showed the potential for very high P concentrations in TD and AMC appeared to affect manure P export in the fall. We showed different P export dynamics between wet and dry years, as well as between events with wet and dry soil AMC. With the hydrograph separation techniques used here, we can monitor and quantify PFP activity, and thus it is possible to manage manure application during periods of high PFP transport, which would likely be the most effective way to reduce P transport to TD. It is still unclear if MI should be a BMP for field P losses as more farms use TD to adapt to climate change. Future work should include more event data to be able to quantify the differences between events following MI on wet and dry soils, as well as if MI is an improvement over the surface application for subsurface P export [24,30]. Also, while subsoiling likely reduced surface runoff, its effect on TD P export is still unknown. The intra-event rainfall analysis showed P concentrations were higher in TD following higher intensity rainfall pulses, yet more data are needed to confirm this, as climate change is expected to have an impact on both AMC and rainfall characteristics.

Supplementary Materials

The following are available online at https://www.mdpi.com/article/10.3390/agriculture12020167/s1, Table S1: Master recession curve constants for each site, Figure S1: Map of study location in the northeast United States. The state of Vermont in Orange and Addison County, Vermont in Black. Source: Jenny Bower, UVM PSS Department. Base layers were obtained from the Vermont Center for Geographic Information and the U.S. Census Bureau., Table S2: Rainfall and tile discharge metrics for each site-event, Figure S2: Images of sampling infrastructure, Table S3: Summary statistics for event tile discharge metrics for each site, Table S4: Summary of flow pathway and source connectivity models for events included in the four-component hydrograph separation analysis. Time to peak in hours in parenthesis and all flow components are in mm, Table S5: Pearson correlation matrix for rainfall and tile discharge metrics for events included in the four-component hydrograph separation. Units for the hydrograph separation components are ‘mm’ and all other units are the same as in Table S1, Table S6: Summary statistics for High- and Low-intensity rainfall pulses during events, Figure S3: Field site schematics, Figure S4: Kendall correlations between event rainfall metrics and day of the hydrologic year (1 October–30 September), including all data (top) and removing the 2021 WY Growing season (bottom), Figure S5: Kendall correlations between days since P/days since manure application and event tile P discharge metrics, Figure S6: Four component hydrograph separation (left) and P concentrations in event samples as a function of tile flow rate (right) for events included in the four-component hydrograph separation analysis at DCS (top) and DCN (bottom). F1 was omitted from this figure because only 1 event was included in the four-component hydrograph separation analysis. Note that the scale for DCN on 25 December 2020 is extended to show the high P concentrations.

Author Contributions

Conceptualization, R.R., D.R. and J.W.F.; Methodology, R.R., D.R. and J.W.F.; Software, R.R.; Validation, R.R., D.R. and J.W.F.; Formal Analysis, R.R.; Investigation, R.R. and D.R.; Resources, R.R. and J.W.F.; Data Curation, R.R. and D.R.; Writing—Original Draft Preparation, R.R.; Writing—Review & Editing, R.R., D.R. and J.W.F.; Visualization, R.R.; Supervision, J.W.F.; Project Administration, J.W.F.; Funding Acquisition, J.W.F. All authors have read and agreed to the published version of the manuscript.

Funding

This article was prepared by the University of Vermont authors using Federal funds under NA180AR4170099 from the National Oceanic and Atmospheric Administration National Sea Grant College Program, U.S. Department of Commerce. The statements, findings, conclusions, and recommendations are those of the author(s) and do not necessarily reflect the views of Sea Grant, NOAA, or the U.S. Department of Commerce.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data available in a publicly accessible repository that does not issue DOIs. Publicly available datasets were analyzed in this study. These data and related items of information have not been formally disseminated by NOAA and do not represent any agency determination, view, or policy. This data can be found here: https://github.com/raruggie/Tile-Drainage-Flow-and-P-in-VT.

Acknowledgments

The authors would like to thank the landowner for cooperation, Dan Needham at the UVM-AETL for assistance with sample analysis, and Jenny Bower for supplying the location map on such short notice.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

P (phosphorus), TD (tile drainage), BMP (best management practice), PFP (preferential flow pathway), HRA (hydrograph recession analysis), SF (slow flow), QF (quick flow), EMMA (end-member mixing analysis), MI (manure injection), LCB (Lake Champlain Basin, EoF (edge of the field) VT (Vermont, USA), AMC (antecedent moisture condition), WY (water year), TP (total P), SRP (soluble reactive P), PP (particulate P), RR (runoff ratio), FWMC (flow weighted mean concentration).

References

  1. Sharpley, A.N. Agriculture, Nutrient Management and Water Quality. In Reference Module in Life Sciences; Elsevier: Amsterdam, The Netherlands, 2018. [Google Scholar] [CrossRef]
  2. Simpson, R.J.; Oberson, A.; Culvenor, R.A.; Ryan, M.H.; Veneklaas, E.J.; Lambers, H.; Lynch, J.P.; Ryan, P.R.; Delhaize, E.; Smith, F.A.; et al. Strategies and agronomic interventions to improve the phosphorus-use efficiency of farming systems. Plant Soil 2011, 349, 89–120. [Google Scholar] [CrossRef]
  3. Vadas, P.A.; Fiorellino, N.M.; Coale, F.J.; Kratochvil, R.; Mulkey, A.S.; McGrath, J.M. Estimating Legacy Soil Phosphorus Impacts on Phosphorus Loss in the Chesapeake Bay Watershed. J. Environ. Qual. 2018, 47, 480–486. [Google Scholar] [CrossRef] [PubMed]
  4. Sharpley, A.; Jarvie, H.P.; Buda, A.; May, L.; Spears, B.; Kleinman, P. Phosphorus Legacy: Overcoming the Effects of Past Management Practices to Mitigate Future Water Quality Impairment. J. Environ. Qual. 2013, 42, 1308–1326. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  5. Sharpley, A.N.; Chapra, S.C.; Wedepohl, R.; Sims, J.T.; Daniel, T.C.; Reddy, K.R. Managing Agricultural Phosphorus for Protection of Surface Waters: Issues and Options. J. Environ. Qual. 1994, 23, 437–451. [Google Scholar] [CrossRef]
  6. Batey, T. Soil compaction and soil management—A review. Soil Use Manag. 2009, 25, 335–345. [Google Scholar] [CrossRef]
  7. Moore, J. Literature Review: Tile Drainage and Phosphorus Losses from Agricultural Land; Stone, Environmental, Technical Report 83; Lake Champlain Basin Program: Grand Isle, VT, USA, 2016. [Google Scholar]
  8. Wolfe, D.W.; DeGaetano, A.T.; Peck, G.; Carey, M.; Ziska, L.H.; Lea-Cox, J.; Kemanian, A.; Hoffmann, M.P.; Hollinger, D.Y. Unique challenges and opportunities for northeastern US crop production in a changing climate. Clim. Change 2018, 146, 231–245. [Google Scholar] [CrossRef] [Green Version]
  9. Klaiber, L.B. Edge-Of-Field Water and Phosphorus Losses. In Surface and Subsurface Agricultural Runoff; University of Vermont: Burlington, VT, USA, 2016. [Google Scholar]
  10. Klaiber, L.B.; Kramer, S.R.; Young, E.O. Impacts of Tile Drainage on Phosphorus Losses from Edge-of-Field Plots in the Lake Champlain Basin of New York. Water 2020, 12, 328. [Google Scholar] [CrossRef] [Green Version]
  11. Madison, A.M.; Ruark, M.D.; Stuntebeck, T.D.; Komiskey, M.J.; Good, L.W.; Drummy, N.; Cooley, E.T. Characterizing phosphorus dynamics in tile-drained agricultural fields of eastern Wisconsin. J. Hydrol. 2014, 519, 892–901. [Google Scholar] [CrossRef]
  12. Kokulan, V.; Macrae, M.; Ali, G.; Lobb, D.; Morison, M.; Brooks, B. Temporal variability in water and nutrient movement through vertisols into agricultural tile drains in the northern Great Plains. J. Soil Water Conserv. 2021, 76, 317–328. [Google Scholar] [CrossRef]
  13. Nazari, S.; Ford, W.I.; King, K.W. Impacts of preferential flow and agroecosystem management on subsurface particulate phosphorus loadings in tile-drained landscapes. J. Environ. Qual. 2020, 49, 1370–1383. [Google Scholar] [CrossRef]
  14. King, K.W.; Williams, M.R.; Macrae, M.; Fausey, N.R.; Frankenberger, J.; Smith, D.R.; Kleinman, P.J.A.; Brown, L.C. Phosphorus Transport in Agricultural Subsurface Drainage: A Review. J. Environ. Qual. 2015, 44, 467–485. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  15. Dils, R.; Heathwaite, A. The controversial role of tile drainage in phosphorus export from agricultural land. Water Sci. Technol. 1999, 39, 55–61. [Google Scholar] [CrossRef]
  16. Beven, K.; Germann, P. Macropores and water flow in soils revisited. Water Resour. Res. 2013, 49, 3071–3092. [Google Scholar] [CrossRef] [Green Version]
  17. Nazari, S.; Ford, W.I.; King, K.W. Quantifying hydrologic pathway and source connectivity dynamics in tile drainage: Implications for phosphorus concentrations. Vadose Zone J. 2021, 20, e20154. [Google Scholar] [CrossRef]
  18. Fiorillo, F. The Recession of Spring Hydrographs, Focused on Karst Aquifers. Water Resour. Manag. 2014, 28, 1781–1805. [Google Scholar] [CrossRef]
  19. Schilling, K.E.; Jones, C.S. Hydrograph separation of subsurface tile discharge. Environ. Monit. Assess. 2019, 191, 231. [Google Scholar] [CrossRef]
  20. Amado, A.A.; Schilling, K.E.; Jones, C.S.; Thomas, N.; Weber, L.J. Estimation of tile drainage contribution to streamflow and nutrient loads at the watershed scale based on continuously monitored data. Environ. Monit. Assess. 2017, 189, 426. [Google Scholar] [CrossRef]
  21. Smith, E.A.; Capel, P.D. Specific Conductance as a Tracer of Preferential Flow in a Subsurface-Drained Field. Vadose Zone J. 2018, 17, 1–13. [Google Scholar] [CrossRef] [Green Version]
  22. Vidon, P.; Cuadra, P.E. Impact of precipitation characteristics on soil hydrology in tile-drained landscapes. Hydrol. Process. 2010, 24, 1821–1833. [Google Scholar] [CrossRef]
  23. Williams, M.R.; King, K.W.; Ford, W.; Buda, A.R.; Kennedy, C.D. Effect of tillage on macropore flow and phosphorus transport to tile drains. Water Resour. Res. 2016, 52, 2868–2882. [Google Scholar] [CrossRef] [Green Version]
  24. Jahanzad, E.; Saporito, L.S.; Karsten, H.D.; Kleinman, P.J.A. Varying Influence of Dairy Manure Injection on Phosphorus Loss in Runoff over Four Years. J. Environ. Qual. 2019, 48, 450–458. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  25. Gaynor, J.D.; Findlay, W.I. Soil and Phosphorus Loss from Conservation and Conventional Tillage in Corn Production. J. Environ. Qual. 1995, 24, 734–741. [Google Scholar] [CrossRef]
  26. Geohring, L.D.; McHugh, O.V.; Walter, M.T.; Steenhuis, T.S.; Akhtar, M.S. Phosphorus Transport into Subsurface Drains by Macropores after Manure Applications: Implications for Best Manure Management Practices. Soil Sci. 2001, 166, 896–909. [Google Scholar] [CrossRef]
  27. Torbert, H.A.; Potter, K.N.; Morrison, J.E., Jr. Management Effects on Nitrogen and Phosphorus Losses in Runoff on Expansive Clay Soils. Trans. ASAE 1996, 39, 161–166. [Google Scholar] [CrossRef]
  28. Chen, Y.; Ren, X. High performance tool for liquid manure injection. Soil Tillage Res. 2002, 67, 75–83. [Google Scholar] [CrossRef]
  29. Uusi-Kämppä, J.; Heinonen-Tanski, H. Evaluating Slurry Broadcasting and Injection to Ley for Phosphorus Losses and Fecal Microorganisms in Surface Runoff. J. Environ. Qual. 2008, 37, 2339–2350. [Google Scholar] [CrossRef]
  30. Coelho, B.B.; Murray, R.; Lapen, D.; Topp, E.; Bruin, A. Phosphorus and sediment loading to surface waters from liquid swine manure application under different drainage and tillage practices. Agric. Water Manag. 2012, 104, 51–61. [Google Scholar] [CrossRef]
  31. Sharpley, A.N. Dependence of Runoff Phosphorus on Extractable Soil Phosphorus. J. Environ. Qual. 1995, 24, 920–926. [Google Scholar] [CrossRef] [Green Version]
  32. Phosphorus Sources—Lake Champlain Basin Program. 2021. Available online: https://www.lcbp.org/our-goals/clean-water/nutrients-and-cyanobacteria/phosphorus-sources/ (accessed on 21 May 2021).
  33. Miller, S.A.; Lyon, S.W. Tile Drainage Increases Total Runoff and Phosphorus Export During Wet Years in the Western Lake Erie Basin. Front. Water 2021, 3, 757106. [Google Scholar] [CrossRef]
  34. Willard, A.; Harris, K.; Kahler, E.; Claro, J.; Danly, S.; Warren, B. Vermont Agriculture and Food System Plan: 2020. p. 110. Available online: https://agriculture.vermont.gov/sites/agriculture/files/doc_library/Vermont%20Agriculture%20and%20Food%20System%20Plan%202020.pdf (accessed on 20 March 2021).
  35. Hanrahan, B.R.; King, K.W.; Williams, M.R. Controls on subsurface nitrate and dissolved reactive phosphorus losses from agricultural fields during precipitation-driven events. Sci. Total Environ. 2021, 754, 142047. [Google Scholar] [CrossRef]
  36. Osterholz, W.R.; Hanrahan, B.R.; King, K.W. Legacy phosphorus concentration–discharge relationships in surface runoff and tile drainage from Ohio crop fields. J. Environ. Qual. 2020, 49, 675–687. [Google Scholar] [CrossRef] [PubMed]
  37. Williams, M.R.; King, K.W.; Baker, D.B.; Johnson, L.T.; Smith, D.R.; Fausey, N.R. Hydrologic and biogeochemical controls on phosphorus export from Western Lake Erie tributaries. J. Great Lakes Res. 2016, 42, 1403–1411. [Google Scholar] [CrossRef] [Green Version]
  38. Soil Survey Staff. Web Soil Survey; United States Department of Agriculture, Natrual Resources Conservation Service: Washington, DC, USA, 2021.
  39. National Oceanic and Atmospheric Administration, United States Department of Commerce. NOWData. Available online: https://www.weather.gov/climateservices/nowdatafaq (accessed on 17 November 2021).
  40. Pierzynski, G.M. Methods of Phosphorus Analysis for Soils, Sediments, Residuals, and Waters, 2nd ed.; North Carolina State University: Raleigh, NC, USA, 2009; Available online: http://www.soil.ncsu.edu/sera17/publications/sera17-2/abstract.htm (accessed on 29 May 2021).
  41. Williams, M.R.; Livingston, S.J.; Penn, C.J.; Smith, D.R.; King, K.W.; Huang, C.-H. Controls of event-based nutrient transport within nested headwater agricultural watersheds of the western Lake Erie basin. J. Hydrol. 2018, 559, 749–761. [Google Scholar] [CrossRef]
  42. Dingman, S.L. Physical Hydrology, 3rd ed.; Waveland Press: Long Grove, IL, USA, 2015. [Google Scholar]
  43. Twombly, C.; Faulkner, J.; Hurley, S. The effects of soil aeration prior to dairy manure application on edge-of-field hydrology and nutrient fluxes in cold climate hayland agroecosystems. J. Soil Water Conserv. 2021, 76, 1–13. [Google Scholar] [CrossRef]
  44. Gregor, M.; Malík, P. HydroOffice. Available online: http://hydrooffice.org/Files/UM%20RC.pdf (accessed on 20 March 2021).
  45. Husic, A.; Fox, J.; Adams, E.; Backus, J.; Pollock, E.; Ford, W.; Agouridis, C. Inland impacts of atmospheric river and tropical cyclone extremes on nitrate transport and stable isotope measurements. Environ. Earth Sci. 2019, 78, 36. [Google Scholar] [CrossRef]
  46. Schilling, K.E.; Helmers, M. Tile drainage as karst: Conduit flow and diffuse flow in a tile-drained watershed. J. Hydrol. 2008, 349, 291–301. [Google Scholar] [CrossRef]
  47. Vidon, P.; Cuadra, P. Phosphorus dynamics in tile-drain flow during storms in the US Midwest. Agric. Water Manag. 2011, 98, 532–540. [Google Scholar] [CrossRef]
  48. Sharpley, A.N. Depth of Surface Soil-runoff Interaction as Affected by Rainfall, Soil Slope, and Management. Soil Sci. Soc. Am. J. 1985, 49, 1010–1015. [Google Scholar] [CrossRef] [Green Version]
  49. Mangiafico, S.S. An R Companion for the Handbook of Biological Statistics Version 1.3.3. 2015. Available online: rcompanion.org/documents/RCompanionBioStatistics.pdf (accessed on 5 July 2021).
  50. R Core Team. R: A Language and Environment for Statistical Computing; R Core Team: Vienna, Austria, 2021; Available online: https://www.R-project.org/ (accessed on 1 February 2020).
  51. Macrae, M.L.; Ali, G.A.; King, K.W.; Plach, J.M.; Pluer, W.T.; Williams, M.; Morison, M.Q.; Tang, W. Evaluating Hydrologic Response in Tile-Drained Landscapes: Implications for Phosphorus Transport. J. Environ. Qual. 2019, 48, 1347–1355. [Google Scholar] [CrossRef] [Green Version]
  52. Van Esbroeck, C.J.; Macrae, M.L.; Brunke, R.I.; McKague, K. Annual and seasonal phosphorus export in surface runoff and tile drainage from agricultural fields with cold temperate climates. J. Great Lakes Res. 2016, 42, 1271–1280. [Google Scholar] [CrossRef]
  53. Tiessen, K.H.D.; Elliott, J.A.; Yarotski, J.; Lobb, D.A.; Flaten, D.N.; Glozier, N.E. Conventional and Conservation Tillage: Influence on Seasonal Runoff, Sediment, and Nutrient Losses in the Canadian Prairies. J. Environ. Qual. 2010, 39, 964–980. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  54. Vidon, P.; Hubbard, L.; Soyeux, E. Seasonal solute dynamics across land uses during storms in glaciated landscape of the US Midwest. J. Hydrol. 2009, 376, 34–47. [Google Scholar] [CrossRef]
  55. Williams, M.R.; Livingston, S.J.; Heathman, G.C.; McAfee, S.J. Thresholds for run-off generation in a drained closed depression. Hydrol. Process. 2019, 33, 2408–2421. [Google Scholar] [CrossRef]
  56. King, K.W.; Williams, M.R.; Fausey, N.R. Contributions of Systematic Tile Drainage to Watershed-Scale Phosphorus Transport. J. Environ. Qual. 2015, 44, 486–494. [Google Scholar] [CrossRef]
  57. Chen, C.; Roseberg, R.J.; Selker, J.S. Using microsprinkler irrigation to reduce leaching in a shrink/swell clay soil. Agric. Water Manag. 2002, 54, 159–171. [Google Scholar] [CrossRef]
  58. Poon, D.; Whalen, J.K.; Michaud, A.R. Re-conceptualizing the Soil and Water Assessment Tool to Predict Subsurface Water Flow through Macroporous Soils. Front. Water 2021, 3, 704291. [Google Scholar] [CrossRef]
  59. Messing, I.; Wesström, I. Efficiency of old tile drain systems in soils with high clay content: Differences in the trench backfill zone versus the zone midway between trenches. Irrig. Drain. 2006, 55, 523–531. [Google Scholar] [CrossRef]
  60. Øgaard, A.F. Freezing and thawing effects on phosphorus release from grass and cover crop species. Acta Agric. Scand. Sect. B-Plant Soil Sci. 2015, 65, 529–536. [Google Scholar] [CrossRef]
  61. Grant, K.N.; Macrae, M.L.; Ali, G.A. Differences in preferential flow with antecedent moisture conditions and soil texture: Implications for subsurface P transport. Hydrol. Process. 2019, 33, 2068–2079. [Google Scholar] [CrossRef]
  62. Ford, W.; Williams, M.R.; Young, M.B.; King, K.W.; Fischer, E. Assessing Intra-Event Phosphorus Dynamics in Drainage Water Using Phosphate Stable Oxygen Isotopes. Trans. ASABE 2018, 61, 1379–1392. [Google Scholar] [CrossRef] [Green Version]
  63. Stone, W.W.; Wilson, J. Preferential Flow Estimates to an Agricultural Tile Drain with Implications for Glyphosate Transport. J. Environ. Qual. 2006, 35, 1825–1835. [Google Scholar] [CrossRef] [PubMed]
  64. NRCS VT and UVM Extension. Vermont Phosphorus Index (VT P-Index); USDA: Washington, DC, USA, 2020.
  65. Beauchemin, S.; Simard, R.R.; Cluis, D. Forms and Concentration of Phosphorus in Drainage Water of Twenty-Seven Tile-Drained Soils. J. Environ. Qual. 1998, 27, 721–728. [Google Scholar] [CrossRef]
  66. Tian, J.; Boitt, G.; Black, A.; Wakelin, S.; Chen, L.; Cai, K.; Condron, L. Mass balance assessment of phosphorus dynamics in a fertilizer trial with 57 years of superphosphate application under irrigated grazed pasture. Nutr. Cycl. Agroecosyst. 2019, 114, 33–44. [Google Scholar] [CrossRef]
  67. Ketterings, Q.M.; Czymmek, K.J.; Reid, W.S.; Wildman, R.F. Conversion of modified morgan and mehlich-III soil tests to morgan soil test values. Soil Sci. 2002, 167, 830–837. [Google Scholar] [CrossRef]
  68. Chikhaoui, M.; Madramootoo, C.A.; Eastman, M.; Michaud, A. Estimating Preferential Flow to Agricultural Tile Drains. In Proceedings of the 2008 Providence, Providence, RI, USA, 29 June–2 July 2008. [Google Scholar] [CrossRef]
  69. Smith, L.; Watzin, M.C.; Druschel, G. Relating sediment phosphorus mobility to seasonal and diel redox fluctuations at the sediment-water interface in a eutrophic freshwater lake. Limnol. Oceanogr. 2011, 56, 2251–2264. [Google Scholar] [CrossRef]
Agriculture 12 00167 g001 550
Figure 1. Boxplots of TD discharge and rainfall metrics for four hydrologic seasons (Fall: October–December, Winter: January–March, Spring: April–June, Summer: July–September). Letters denote significant differences between seasons. Jitter is used for plotting only and does not correspond to dates.
Figure 1. Boxplots of TD discharge and rainfall metrics for four hydrologic seasons (Fall: October–December, Winter: January–March, Spring: April–June, Summer: July–September). Letters denote significant differences between seasons. Jitter is used for plotting only and does not correspond to dates.
Agriculture 12 00167 g001
Agriculture 12 00167 g002 550
Figure 2. Event TP and SRP concentrations as a function of TD flow rate. Regressions were performed after removing the points in green, i.e., samples taken within the first 10 h of events.
Figure 2. Event TP and SRP concentrations as a function of TD flow rate. Regressions were performed after removing the points in green, i.e., samples taken within the first 10 h of events.
Agriculture 12 00167 g002
Agriculture 12 00167 g003 550
Figure 3. Boxplots showing the difference between pulse and TD discharge metrics for high and low-intensity pulses.
Figure 3. Boxplots showing the difference between pulse and TD discharge metrics for high and low-intensity pulses.
Agriculture 12 00167 g003
Table
Table 1. Site management practices and timing. P.A.E. is P application equivalent. Nitrogen side dress was entirely urea ammonium nitrate (UAN), and the cover crop was Winter Rye (Secale cereale).
Table 1. Site management practices and timing. P.A.E. is P application equivalent. Nitrogen side dress was entirely urea ammonium nitrate (UAN), and the cover crop was Winter Rye (Secale cereale).
F1 DC
DateActionP.A.E (g/ha)DateActionP.A.E (g/ha)
10/6/18Manure injected152,4359/29/19Cover crop planted
5/19/19Pop-up22,94910/12/19Manure injected152,435
Corn planted 4/27/20Cover crop terminated154
Starter15,3004/30/20Light chisel till
Cover crop terminated (A standard application rate of 0.34 kg (0.75 lbs) of glyphosate (acid equivalent) per acre assumed)154 Pop-up22,949
7/2/19Nitrogen side dress Corn planted
9/25/19Corn harvested Starter15,300
9/29/19Cover crop planted 6/23/20Nitrogen side dress
10/6/19Manure injected152,4359/12/20Corn harvested
5/5/20Light chisel till 9/19/20Cover crop planted
Pop-up 9/28/20Manure injected121,948
Corn planted 11/13/20Subsoiled
Starter15,3004/15/21Cover crop terminated154
Cover crop terminated1545/21/21Light chisel till
6/20/20Nitrogen side dress 5/24/21Pop-up22,949
9/18/20Corn harvested Corn planted
9/22/20Cover crop planted Starter15,300
10/6/20Manure injected121,9486/20/21Nitrogen side dress
Table
Table 2. WY-annual estimates of TD TP export (g/ha). Volume is mean event volume (mm), RR is mean event runoff ratio, G is growing season, and NG is a non-growing season.
Table 2. WY-annual estimates of TD TP export (g/ha). Volume is mean event volume (mm), RR is mean event runoff ratio, G is growing season, and NG is a non-growing season.
2019 2020 2021
BaseflowEvent BaseflowEvent BaseflowEvent
SiteMetricGNGGNGWY TotalGNGGNGWY TotalGNGGNGWY Total
F1TP.Load3.34.5536852.31396.1 a1.44.36414841554.1 b
Volume 8.9717.65 10.41
RR 0.530.6 0.51
DCNTP.Load 6.721.26432683359.9 c15.416742.2762.81536.4 d
Volume 2.7317.02 9.1823.2
RR 0.090.92 0.350.59
DCSTP.Load 1.65.12113131340.6 e3.73.9533.1416.5957.2 f
Volume 1.035.18 7.5520.15
RR 0.030.26 0.290.41
Note sampling errors resulted in the following number of missing events and load estimates: a Growing: 1 event, 171 g/ha, Non-Growing: 4 events, 261 g/ha; b Growing: 3 events, 64 g/ha, Non-Growing: 7 events, 630 g/ha; c Non-Growing: 14 events, 2770 g/ha; d Growing: 1 event, 2 g/ha; e Non-Growing: 14 events, 555 g/ha; f Growing: 1 event, 3 g/ha.
Table
Table 3. Summary statistics for TD flow metrics, fraction of total flow (FTF), and time to peak (TTP) for total TD event flow and hydrograph separation components from the flow pathway (QF, SF), source connectivity (new, old), and four-component models.
Table 3. Summary statistics for TD flow metrics, fraction of total flow (FTF), and time to peak (TTP) for total TD event flow and hydrograph separation components from the flow pathway (QF, SF), source connectivity (new, old), and four-component models.
Flow MetricFTFTTP
(N = 25)(N = 25)(N = 25)
Flow Volume (mm)
Mean (SD)15.8 (12.7)NA (NA)12.8 (13.3)
Median [Min, Max]13.4 [1.23, 55.0]NA [NA, NA]8.83 [0.917, 49.6]
Missing0 (0%)25 (100%)0 (0%)
QF
Mean (SD)13.2 (10.6)0.821 (0.108)12.0 (11.5)
Median [Min, Max]11.3 [0.893, 45.0]0.828 [0.579, 0.973]8.83 [0.917, 49.6]
SF
Mean (SD)2.61 (2.51)0.179 (0.108)34.3 (21.2)
Median [Min, Max]2.23 [0.0955, 10.1]0.172 [0.0273, 0.421]28.5 [10.0, 89.7]
New
Mean (SD)9.23 (9.61)0.493 (0.191)13.2 (13.5)
Median [Min, Max]7.25 [0.325, 41.2]0.542 [0.102, 0.803]9.75 [1.33, 52.8]
Old
Mean (SD)6.58 (3.81)0.507 (0.191)11.3 (9.57)
Median [Min, Max]6.60 [0.900, 13.8]0.458 [0.197, 0.898]8.17 [0.833, 32.3]
QF old
Mean (SD)4.30 (2.50)0.344 (0.150)11.3 (9.57)
Median [Min, Max]4.11 [0.575, 9.63]0.321 [0.0845, 0.696]8.17 [0.833, 32.3]
QF new
Mean (SD)8.89 (9.00)0.477 (0.179)13.2 (13.5)
Median [Min, Max]7.09 [0.318, 38.2]0.529 [0.100, 0.776]9.75 [1.33, 52.8]
SF old
Mean (SD)2.28 (1.92)0.163 (0.105)32.0 (17.8)
Median [Min, Max]1.95 [0.0955, 7.25]0.133 [0.0273, 0.419]27.1 [9.92, 66.7]
SF new
Mean (SD)0.337 (0.706)0.0157 (0.0270)28.0 (26.3)
Median [Min, Max]0.0262 [0, 3.08]0.00279 [0, 0.121]28.2 [0, 89.8]
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Ruggiero, R.; Ross, D.; Faulkner, J.W. Tile Drainage Flow Partitioning and Phosphorus Export in Vermont USA. Agriculture 2022, 12, 167. https://doi.org/10.3390/agriculture12020167

AMA Style

Ruggiero R, Ross D, Faulkner JW. Tile Drainage Flow Partitioning and Phosphorus Export in Vermont USA. Agriculture. 2022; 12(2):167. https://doi.org/10.3390/agriculture12020167

Chicago/Turabian Style

Ruggiero, Ryan, Donald Ross, and Joshua W. Faulkner. 2022. "Tile Drainage Flow Partitioning and Phosphorus Export in Vermont USA" Agriculture 12, no. 2: 167. https://doi.org/10.3390/agriculture12020167

APA Style

Ruggiero, R., Ross, D., & Faulkner, J. W. (2022). Tile Drainage Flow Partitioning and Phosphorus Export in Vermont USA. Agriculture, 12(2), 167. https://doi.org/10.3390/agriculture12020167

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop