Effects of Drainage Water Management in a Corn–Soy Rotation on Soil N2O and CH4 Fluxes

Drainage water management (DWM), also known as controlled drainage, is a best management practice (BMP) deployed on drainage ditches with demonstrated success at reducing dissolved nitrogen export from agricultural fields. By slowing discharge from agricultural ditches, subsequent anaerobic soil conditions provide an environment for nitrate to be reduced via denitrification. Despite this success, incomplete denitrification might increase nitrous oxide (N2O) emissions and more reducing conditions might increase methanogenesis, resulting in increased methane (CH4) emissions. These two gases, N2O and CH4, are potent greenhouse gases (GHG) and N2O also depletes stratospheric ozone. This potential pollution swapping of nitrate reduction for GHG production could negatively impact the desirability of this BMP. We conducted three years of static chamber measurements of GHG emissions from the soil surface in farm plots with and without DWM in a corn–soybean rotation on the Delmarva Peninsula. We found that DWM raised the water table at the drainage ditch edge, but had no statistically significant effect on water-filled pore space in the field soil surface. Nor did we find a significant effect of DWM on GHG emissions. These findings are encouraging and suggest that, at least for this farm site, DWM can be used to remove nitrate without a significant tradeoff of increased GHG emissions.


Introduction
Excess agricultural nitrogen (N) and phosphorus (P) addition to the environment has been a long-term problem in the United States, contributing to algal blooms and eutrophication [1]. Within the Chesapeake Bay watershed agricultural production on the Delmarva peninsula accounted for~90% of the nitrogen inputs [2][3][4], much of which is derived from synthetic nitrogen fertilizer [5]. Despite the need for nitrogen to increase crop yields, more than half of the N applied can be lost to the environment in a variety of pathways such as leaching and gas emissions [6,7]. Surface ditches and tile drains installed to remove water from wet fields act as conduits for export of dissolved nutrients. Drainage water management (DWM) is a best management practice (BMP) designed to slow this movement and thereby reduce reactive N loss from agricultural fields that have drainage ditches [8,9].
Drainage water management (DWM) functions by placing drainage control structures (DCS) at the end of drainage ditches with the goal of slowing water discharge and raising water levels in ditches. By keeping water in the field for longer, the water table upstream of the ditch rises, creating more widespread anaerobic soil conditions [10]. When N, as nitrate, is exposed to more reducing conditions microbially driven denitrification is more likely [11,12]. Some studies have demonstrated that DWM has been effective at reducing nitrate loading into downstream waters [13][14][15], including preliminary results from the farm of present study [16]. However, incomplete denitrification can result in nitrous oxide (N 2 O) production, a potent greenhouse gas (GHG) and a reactant in stratospheric ozone destruction [17,18]. In addition, nitrification has also been documented to produce N 2 O as a byproduct [11]. Agricultural soil management, driven by anthropogenic nitrogen inputs, is responsible for 75% of the total N 2 O emissions in the US and 45% of total N 2 O emissions in North America [19,20]. In addition to N 2 O, reducing conditions can also create higher methane (CH 4 ) emissions via methanogenesis [21]. Methane is the second most important GHG, with potency 25-fold greater than CO 2 in a 100 year time frame [22].
Broadly, soil N 2 O fluxes are controlled primarily by inorganic-N content, carbon availability, soil moisture, microbial community abundances, soil temperature, soil pH, and soil texture [11,12]. Fertilizer N can be rapidly converted to nitrate, through nitrification, and nitrate can then become substrate for denitrification. Both nitrification and denitrification can contribute to large N 2 O emission spikes [23], especially when high soil nitrate concentrations co-occur with high soil moisture conditions [24,25]. The combination of high N input with DWM that retains water in the field could result in lower hydrologic losses of N as nitrate through the ditch but with the tradeoff of higher losses of N 2 O.
Previous work examining the impact of water table changes and DWM on GHG emissions in agriculture fields has yielded contradictory results. Some studies found that undrained fields had higher N 2 O emissions than drained fields [26][27][28][29] while other studies found that there was no difference in N 2 O emissions between drained and undrained fields [30][31][32][33]. For CH 4 , one study [26] found higher emissions in undrained fields while others found no drainage effect on CH 4 emissions [30][31][32]. Some of these studies did not directly evaluate DWM manipulations, creating higher and lower water tables at specific times of the year, rather they explored drainage impact in general [26,27,32], or drainage impact with subsurface irrigation [29,31]. The variability in conclusions and paucity of direct assessments of DWM mean that further research is needed to understand the impacts of this particular BMP on the potential for increased GHGs. Here, we compare measured soil emissions of N 2 O and CH 4 in corn and soybean fields with controlled drainage (DWM) and free drainage (non-DWM) at a farm in eastern Maryland, on the Delmarva Peninsula. We hypothesized that DWM would increase soil GHG emissions.

Study Site Location
The research site was located on a private farm within the Choptank watershed on the coastal plain of the Delmarva Peninsula on the Eastern Shore of Maryland (Figure 1, N 39.11973, W 75.8069). Historically a wet landscape, the land was made arable by installing drainage ditches. Harvested crops mostly contribute feed to nearby concentrated poultry production. Precipitation averages 110 cm per year. This study site is an operational farm with a collaborating farmer willing to accommodate our research needs while working within the timing of standard agricultural management. A DWM system was already in place and a direct evaluation of that system was possible in a relatively small area. The ditches split the site into four similarly sized fields ranging from 3 to 6 hectares each, with an average length of 220 m. The dominant soil type is a Fallsington series on 0-2 percent slopes classified as fine-loamy, mixed, active, mesic Typic Endoaquults.

Farm Management
The crop management plan was a no-till, corn-soybean rotation with a winter radish cover crop. The typical growing season ranged from April to October. Approximately 213 kg ha −1 of synthetic urea-ammonium nitrate (UAN) was applied during corn years and no N fertilizer was applied during soybean years. Approximately 22 kgN ha −1 was knifed-in during corn planting, and the other 191 kgN ha −1 was broadcast-applied in late May or early June. During corn planting, the fields received a sulfur addition of 15 kgS ha −1 . Liquid or dry potash was applied during corn and soybean years. Herbicide was applied prior to or soon after planting and sometimes in June to suppress undesired vegetative growth. Prior to harvest, daikon radish seed was applied by air to serve as the winter cover crop. Manure and lime were applied every 5 years, which did not occur during this project's duration. The site had pivot irrigation that was used to supplement natural precipitation events in the summer months. All the research fields had the same crop rotation schedule, N fertilization rate, and irrigation rate from 2017 to 2019.

DWM Management
The DWM strategy consisted of inserting and removing boards within the drainage control structures (DCS) to slow the ditch discharge while accommodating standard farm practices ( Figure S1). For farm activities such as planting, N fertilization, and harvest, the DCS boards were removed and fields were drained to create a drier condition for equipment use. This management timeline was dependent on precipitation events, groundwater levels, and evapotranspiration; therefore the strategy was modified to adapt to changing environmental conditions. We refer to the time period when the boards were in place in the DWM DCS as "active" and the period when they were removed from the DWM DCS as "inactive". Figure 1 shows a site map depicting the static chamber transect locations and borders of study plots. Non-DWM, or free drainage, was implemented on two fields (A and B) and controlled drainage, or DWM, implemented in the other two fields (C and D). In each field the soil chamber sampling design consisted of a transect of eight chambers installed perpendicular to the ditch. Each transect was approximately 12 m in length and the soil chamber collars were equally spaced. The replicate plots were adjacent to each other in order to increase the horizontal fetch of micrometeorological measurements of a companion part of this study that will be reported elsewhere.

Experimental Design
Soil and soil gas flux measurements were taken approximately monthly, depending on climatic conditions, but more temporally intensive sampling was done around planting, irrigation, N fertilization, and harvest to better characterize trace gas dynamics during those periods. Because preliminary studies showed some diel variation ( Figure S2), the order of chamber sampling was randomized to ensure that the same field was not measured at the same time of day during every sampling date, thus minimizing potential sampling biases due to diel variation. Fluxes were sampled starting at approximately 8:00 EST and extended until approximately 17:00 EST.
Two slotted PVC (5 cm diameter) piezometers, used to measured groundwater levels, were installed approximately 1.5-2 m deep in the four research fields. One piezometer was installed within 1 m of the perimeter of the ditch while another was placed in the field 15 m perpendicular to the ditch piezometer ( Figure 1). The piezometers were located within 15 m of the chamber transect. Solinst water level data loggers were placed in each piezometer which recorded pressure and temperature every 30 min. The pressure was then converted to depth of water below the surface by subtracting barometric pressure recorded separately by a Solinst barologger. Thirty-minute water depths were aggregated to daily mean for this analysis.
There were four primary controls on the statistical power of this experimental design: (1) plot replication, (2) number of chamber measurements in each plot, (3) number of days of sampling over three years and (4) number of gas concentration measurements used to calculate each flux. This experimental design was limited to two replicate plots per treatment because four was the maximum number of flux towers that could be established in the companion micrometeorological study. The within-plot statistical power of this study was provided by eight chambers per field, which was needed for measurements known to have high spatial variability, as is the case with N 2 O and CH 4 [34]. This study had excellent statistical strength in the number of gas concentration measurements per flux calculation, as there was a single chamber headspace concentration recorded every second (1 Hz), which increased the flux estimation accuracy, lowered the flux detection limit, reduced the time needed for each flux measurement, and afforded improved statistical analysis [35].

Flux Measurements
A Picarro G2308 cavity ringdown spectrometer and a static soil chamber system were used to measure N 2 O and CH 4 soil gas fluxes. The instrument specifications state that the raw precision (1σ) of the instrument was <25 ppb + 0.05% of the reading for N 2 O and <10 ppb + 0.05% of the reading for CH 4 . The chambers consisted of a PVC chamber top with a curled vent, 12 cm length and 1 mm ID, and a matching PVC chamber collar. The chamber was approximately 25 cm in diameter. The chamber collar height was 5 cm and was inserted 3-4 cm into the soil. The chamber top height was 2.5 cm and had a rubber gasket around the perimeter to ensure a tight fit onto the collar. The approximate total volume of the chamber, dependent on insertion depth into the soil, was 4 L. The flux measurement system used a recirculation technique through inlet and outlet tubing connected to the instrument, air pump, and the chamber top. Twenty-meter lengths of 6 mm OD polyethylene tubing connected the gas analyzer to the chamber top. Given the average flow rate of the pump (~230 mL min −1 ), the tubing diameter, and length, it took approximately 45 s for gas to reach the analyzer from the chamber. Because of the high frequency (1 Hz) of reporting concentration by the instrument and its high sensitivity to change in concentration, the chamber top was left on each collar for only~5 min, which was sufficient time to detect a non-zero slope of increasing or decreasing headspace concentrations, while minimizing a potential chamber artifact of altering the concentration gradient between soil and headspace gas. The slope was used to calculate a flux using Equation (1) where dConc/dt = slope (µL L −1 min −1 ), P = air pressure (atm), R = gas law constant (L atm mol −1 K −1 ), T = air temperature (K), V = chamber volume (L), and A = chamber area (m 2 ). Flux units of µmol m −2 min −1 were used for statistical analyses, but when flux interpolations and annual estimates were conducted, the flux units were converted to µg N 2 O-N m −2 min −1 .

Soil Chemical and Physical Characteristics
Volumetric water content (VWC) of the soil was measured using a manual Decagon 5TM soil moisture probe. The probe was inserted vertically into the soil surface to a depth of 5 cm. At each chamber location (n = 32), three VWC measurements were taken outside the chamber area, and the mean calculated. The probe calibration was checked in the laboratory prior to field usage and the VWC was obtained using the Topp equation [36].
Soil samples used for gravimetric soil moisture, nitrate and ammonium concentration analyses were taken using a 2.2 cm diameter soil push probe during each gas flux measurement. At each chamber (n = 32) three, 5 cm deep soil plugs surrounding the chamber were collected and placed together into a labeled bag. The bag was placed in a refrigerated cooler, transported back to the laboratory and processed within 48 h. At the lab, the soil within each bag was homogenized and sieved using a 10 mesh (2 mm) pan. Six grams of the field-moist sieved soil were extracted in 30 mL of 2M potassium chloride (KCL), placed on an oscillating shaker table for at least 2 h, and vacuum filtered through a Whatman No. 42 filter paper. Extracts were then analyzed for nitrate, nitrite, and ammonium concentration using flow injection analysis Lachat autoanalyzer [37]. If samples were not processed in the analyzer within 12 h, they were stored at 4 • C until analysis within 30 days of soil collection. Soil push probe samples from five sampling dates selected over the three-year study were also analyzed for total soil carbon and total soil nitrogen using a Carlo Erba NC2500 elemental analyzer.
The surface bulk density of the soil was measured using a 5 cm length, 5 cm diameter PVC core. Air-dried soil was sieved with 10 mesh (2 mm) pan and weighed. The bulk density was calculated as the <2 mm air-dried mass divided by the volume. Two bulk density surface cores were taken in each field for two years of this study and averaged by field and year. Surface bulk density measurements were used to convert VWC to waterfilled pore space (WFPS) by calculating the total porosity assuming a particle density of 2.65 g cm −3 [38]. Soil texture was assessed twice using hand augered soil samples at 10 cm intervals to an approximate depth of 1 m, twice, in each field. The soils were sieved using a 10 mesh (2 mm diameter) pan, air-dried for at least 72 h, and measured for soil texture using the hydrometer method [39].

Non-Parametric DWM Impact Model
The purpose of this model was to determine if the dependent variables (gas fluxes, piezometer depth, soil moisture) were impacted by DWM. Because of non-normality and heteroskedasticity of the dependent variable data, a non-parametric Wilcoxon test was used. The piezometer water depth, soil moisture and fluxes were assessed separately as dependent variables of interest. The independent variables were treatment (DWM/controlled drainage; versus non-DWM/free drainage) and board period (active-DCS boards inserted in the DWM plots; versus inactive-DCS boards removed for the farmer's management activities). A pairwise post hoc test was performed to estimate p-values assessing the difference between board periods under both treatment conditions.

Mixed-Effects Model
A generalized linear mixed-effects model (GLMER) in R program was used to assess the impact of treatment and other predictor variables on N 2 O and CH 4 fluxes. Both N 2 O and CH 4 fluxes were non-normal and right skewed, so a natural log transformation was applied. A small positive offset was added to make all values positive before log transformation [40]. A gamma log link distribution was applied in the GLMER model [41].
For N 2 O, the independent variables of treatment, fertilization period, year, and board period were categorical variables assessed as fixed effects. Chamber location (1-8, with 1 being closest to the ditch and 8 being furthest) was a continuous variable serving as a measure of distance from ditch edge. After distance from the ditch was found to not be a significant fixed effect (p > 0.05), chamber number was then considered a random effect. Sampling dates were considered part of the fertilization period if they occurred within a window of 4 weeks after a fertilization event. Several studies have shown a peak of N 2 O emissions within days of fertilization, followed by a gradual decline over 2-4 weeks [24,42,43]. Board period was the same categorical variable as described in Section 2.7. Interaction effects between all fixed effects were also evaluated. Final models were selected based on the lowest AIC value of multiple fixed effect combinations and interactions. An estimated marginal means post hoc test was conducted when p-values of main effects were less than 0.05. Model residuals were visually analyzed for normality and heteroskedasticity.
For CH 4 , the independent variables of treatment, seasonal period, board period, and year were categorical variables assessed as fixed effects. Chamber location (1-8, with 1 being closest to the ditch and 8 being furthest) was initially tested as a continuous variable serving as a measure of distance from ditch edge but was not significant (p > 0.05), therefore distance was considered a random effect. Seasonal period was assigned a value based on inclusion within the data range of 1 April-30 June. These dates were selected because they corresponded to the beginning of the warmer spring season with wet soil conditions, prior to drier summer months. The final model was selected based on the lowest AIC value of multiple fixed effect combinations and interactions. An estimated marginal means post hoc test was conducted when p-values of main effects were less than 0.05. Model residuals were visually analyzed for normality and heteroskedasticity.

Regression Statistics
To examine how environmental variables controlled N 2 O emissions, a stepwise linear regression was created using soil temperature, soil water-filled pore space (WFPS), and nitrate concentration as the predictor variables for N 2 O emissions. The data were unaggregated, so each chamber flux measurement was an observation, for a total of 536 measurements with a complete suite of data for all dependent and independent variables. All data were tested for normality using the Shapiro-Wilk normality test and were log-transformed. For WFPS, a quadratic term was added to the model to account for an expected non-linear response [11]. The stepwise regression created multiple models with different combinations of predictor variables. The model with the lowest AIC was selected. Homoscedasticity and normality were assessed in the final model residuals.  Table S1, bounded by farm management events and observed seasonal impacts of patterns of measured fluxes. For periods in which no visually obvious N 2 O spikes occurred, an arithmetic mean was taken of all the flux measurements within that time period and multiplied by the number of days within that period. For periods in which a spike of N 2 O flux was observed followed by a sustained drop in flux values, usually within days or weeks after N fertilization in corn years or shortly after temperature and moisture increase in the soybean year, a decay curve and function were fitted to the data for the duration of those observations ( Figure S3). Daily fluxes were then summed to obtain a total period value. The emission values in each period (both arithmetically derived and decay curve derived) were summed to obtain an annual flux measurement. There were 3-4 periods per year with 2 to 5 flux measurements per period. This procedure was repeated to calculate an annual flux for each chamber. The chamber annual fluxes were averaged by treatment and a 95% confidence interval (CI) was calculated, representing spatial flux variability.
For CH 4 , each year was divided into 3 periods bounded by seasonal impacts on measured fluxes (Table S2). A mean flux was calculated for each period, multiplied by the number of days within that period, and period fluxes were summed to obtain an annual CH 4 flux. This procedure was repeated for each chamber and the chambers were averaged by treatment to obtain a 95% CI representing the spatial variation within the treatment.

Mixed-Effects Model Extrapolation
For an independent check on the interpolation method, annual CH 4 and N 2 O fluxes were also calculated by extrapolating the mean flux values reported in the mixed-effects model results for each year. This more simplistic approach did not account specifically for observed spikes in gas fluxes but represented the spatial and temporal variation accounted for in the mixed-effects model. An error term for each annual estimate was calculated based on the standard errors of the estimates of the mixed-effects model output.

Regression Model N 2 O Flux Interpolation
A regression-based model was used to obtain another independent annual N 2 O flux estimate using environmental variables of soil moisture, soil nitrate, and soil temperature. Linear interpolations were made between each of the measured data points for these environmental drivers. The interpolated daily variable values were applied to a stepwise linear regression equation. Predicted hourly N 2 O fluxes were summed to an annual flux. This was done for each chamber (n = 32) and an average and standard error were calculated by treatment.

Environmental Variables
The following section describes the climatic and soil variables measured to complement the flux measurements ( Figure 2).

Precipitation and Air Temperature
Daily precipitation and irrigation amounts were recorded over the course of this study (Figure 2A). The highest daily rainfall was 7.84 cm in February 2018. Seasonally averaged over three study years (2017-2019), March-May had the highest average rainfall with 31.8 ± 1 cm, followed by December-February with 26.4 ± 3.6 cm, June-August with 19.1 ± 8.6 cm, and September-November with 9.0 ± 1.9 cm. The daily average air temperature over three years was 14.1 • C and ranged from −13.0 to 30.1 • C (data not shown).

Soil Temperature and Moisture
Soil temperature and moisture at 5 cm depth were measured at the same frequencies as chamber fluxes. Soil temperature ranged from 2.1 to 29.1 • C and followed expected seasonality trends ( Figure 2B). The measured mean percent water-filled pore space (%WFPS) by treatment ranged from 13% to 67% over the study period ( Figure 2C). The highest recorded soil moisture contents occurred in the winter (December-February) and spring (March-May), and the lowest soil moisture contents occurred during the summer (June-August) and fall (September-November). There was a particularly dry period in July 2018, when the lowest soil moisture content of this study was recorded. Statistical analysis of WFPS can be found in results Section 3.4.

Soil Extractable Nitrogen
The soil extractable nitrate in the top 5 cm, measured at each chamber sampling, ranged from 0.38 to 34.6 mg N kg −1 ( Figure 2D). In April 2017 and 2019, 22 kg N ha −1 as UAN fertilizer were applied during planting and in June 2017 and 2019 200 kg N ha −1 as UAN were applied, as indicated by arrows in Figure 2. In general, there was a positive association between N fertilization and soil nitrogen concentration. Soil nitrate values increased soon after application of N fertilizer in April and June of 2017 and 2019. After this period, the soil nitrate concentrations decreased throughout the remainder of the year. In the 2018 soybean crop, the soil nitrate values rose in the spring but did not reach concentrations as high as the years with N fertilization. For soil extractable ammonium values, the concentrations remained below 10 mg N kg −1 during the study period, except one measurement in June 2018 (data not shown). Ammonium values were more variable than nitrate and higher during summer than during winter.

Soil Texture and Bulk Density
Soil texture varied spatially within and between the four fields. Augered soil cores in field C had higher clay and silt content compared to fields A, B, and D at depths ranging from 20 to 60 cm ( Figure S4 and Table S3). The 5 cm surface bulk density values, %C, and %N also varied by field, with higher %C and %N values associated with higher clay content (Table S4).

Piezometers
The daily mean piezometer water depth was aggregated by treatment and piezometer location during the inactive and active DWM DCS board periods (Figure 3). In the ditch piezometers, in the DCS inactive (unshaded) periods, the DWM depth values were generally deeper than the non-DWM plot values, indicating spatial heterogeneity across the field. During the DCS active (shaded) periods, however, the lines on Figure 3 tend to converge, indicating that the DCS active installation of boards in the DWM plots increased the ditch piezometer water levels. In the field piezometers, the effect of board placement is less clear. The statistical evaluation of these piezometer data were combined analysis of WFPS data and are found in Section 3.4.

Chamber Flux Summary
From 2017 to 2019, 959 chamber flux measurements were taken. Most fluxes were small and the data had right-skewed distribution. When calculating the fluxes using the change in slope over time, assessed by linear regression, 96% of the methane slopes were significantly different from zero (p < 0.05) and of those 21% were negative, indicating CH 4 consumption events. For N 2 O, 85% of the N 2 O slopes were significantly different from zero (p < 0.05) and 15% of those were negative. Because of the data skew, a log transformation was applied for parametric statistical analyses. An analysis of the minimum detectable flux (MDF) is provided in the supplementary materials (Table S5). Figure 2E shows the N 2 O-N flux data, averaged by treatment. The highest emissions occurred after N fertilization events (planting or broadcast) in the April-July period of 2017 and 2019, marked by arrows. In 2017, the larger N 2 O peak occurred during the second fertilization and in 2019, after the first fertilization. The large response to the first fertilization in April 2019 N application may have resulted from co-occurrence of a precipitation event the day following application. Soil N 2 O consumption events occurred at various times throughout this study at the individual chamber scale, but never as the average flux from all chambers at the field scale. This indicates small scale flux variability within a field, but not net negative fluxes at the field scale.
For CH 4 , the highest CH 4 emissions occurred in the spring period (April-June; Figure 2F). Several individual chambers showed net consumption of atmospheric CH 4 , but only rarely did this translate to net consumption as the mean at the field or treatment scale.

Non-Parametric DWM Model for Effects on Piezometers, WFPS, and Gas Fluxes
The non-parametric statistical analysis results for dependent variables of ditch piezometer, field piezometer, and surface WFPS showed differences in the interaction of treatment and board period (Figure 4). For the ditch piezometers, the mean water depth below ground was significantly deeper during the inactive DCS period than the active DCS period in the non-DWM sites. However, this trend was reversed in the DWM sites, with a higher mean water depth in the active DCS period compared to the inactive DCS period, indicating that the additional boards in the DWM sites raised the water table near the ditch, even during the dry summer. The field piezometer water levels showed the same trends as the ditch piezometers, except that the mean water level between active and inactive DCS periods was not significantly different in the DWM condition, indicating that the additional boards in the DWM sites had less of an effect on the water table in the mid-field, approximately 15 m from the ditch. The DWM only raised the water table in the mid-field piezometers to the same level as the non-DWM plots in the active DCS season, but did not exceed it, as was the case for the ditch piezometers ( Figure 4). The WFPS showed no significant differences by treatment or board period, indicating the DWM did not affect near surface moisture content. The same non-parametric analysis for gas fluxes showed no significant differences by treatment or board period, indicating the DWM did not impact the soil gas fluxes in the fields ( Figure 5).

Mixed-Effects Model Results for N 2 O and CH 4 Fluxes
The purpose of the mixed-effects model was to assess the impact of the DWM treatment and other dependent variables. For N 2 O fluxes, the effects of the post-N fertilization period, year, chamber location, and interactions were assessed using two sub-models. In the first sub-model, the N fertilization effect was evaluated using only 2017 and 2019 data (corn years) because 2018 data (soybean year) did not have a N fertilization period. In the second sub-model the effects of treatment and year were assessed using all years' data, but N fertilization variable was not included. The mixed model results (Table 1) show that DWM did not significantly affect N 2 O fluxes, which is consistent with the non-parametric model results (Figures 4 and 5).  Table 2). No interaction effects between independent variables were significant.  For CH 4 , DWM treatment, year and season were tested. CH 4 emissions did not have significant differences based on year or treatment ( Table 1). The spring seasonal period had higher CH 4 emissions compared to the non-spring periods ( Table 2).

N 2 O Multiple Regression
There was a weak positive linear relationship (r 2 = 0.12, p < 0.001, Figure S5) between soil nitrate and N 2 O fluxes. The largest N 2 O fluxes occurred when soil nitrate concentrations were high, but some large N 2 O fluxes also occurred at lower soil nitrate concentrations. For soil temperature, the higher N 2 O emissions occurred above 12 • C, although it should be noted that N fertilization events also only occurred during periods above this temperature. For WFPS, the highest N 2 O emission events occurred at intermediate WFPS values, but not all N 2 O flux events between those WPFS values resulted in high emissions. The highest N 2 O emissions typically occurred when the following combination of conditions existed: the soil temperature was above 12 • C, nitrate concentration was at least 10 mg N kg −1 , and WFPS was between 35% and 60% ( Figure 6). The stepwise linear regression yielded the lowest AIC when the model included all terms, including a quadratic term for WFPS. The most important predictor variable was soil nitrate followed by soil temperature, and WFPS. The full model (Equation (2)) and each predictor variable were significant (p < 0.001), and the full model yielded an r 2 value of 0.23 (Table S6). ln (N 2 O + 1) = (0.029 * NO 3 N) + (0.021 * soil temp) + (0.056 * WFPS) − (0.001 * WFPS 2 ) − 1.37 (2) where N 2 O is N 2 O flux (µg N 2 O-N m −2 min −1 ), NO 3 -N is concentration of nitrate (mg NO 3 -N kg −1 dry soil), WFPS is water-filled pore space (%), soil temp is soil temperature ( • C).

Annual Soil Emission Flux Estimates
Annual soil gas flux estimates were derived from the three methods described in Section 2.10 (Tables 3 and 4). All the methods indicated that the N 2 O-N emissions in the corn years of 2017 and 2019 were elevated in both the DWM and non-DWM conditions compared to the soybean year of 2018. While the annual estimates show higher emissions in 2019 compared to 2017 which is a result of a higher peak N 2 O emission and nitrate measurements after the first fertilization event in 2019 (Figure 2), the non-parametric and the mixed-effects model tests discussed above indicated that the difference between 2017 and 2019 fluxes was not significant. The seasonally-based interpolation method yielded higher annual flux estimates, but there was generally broad agreement between the three methods, given the uncertainty terms.  Table 4. Summary of comparison between the two annual CH 4 flux calculation methods. The SE associated with the seasonally-based interpolation represents the spatial heterogeneity among chambers.
The mixed-effects model mean extrapolation error was derived from the model output error term. The seasonally-based interpolation method permitted partitioning fluxes into periods. Individual seasonal period estimates (mean ± 95% CI) for N 2 O showed that brief (6-8 weeks) seasonal periods following the fertilization events, which were quantified by the decay function, contributed nearly half of the total annual N 2 O flux (Table S1).
There were small net positive annual CH 4 fluxes in each year and treatment (Tables 4 and S2). Consistent with the mixed-effects model, there were no clear differences between years and treatments, except that the seasonally-based interpolation estimates were higher for 2019. This may be an artifact caused by unusually high fluxes measured in a few chambers during a period with sparse temporal coverage.

Groundwater Level
The DWM treatment raised the near-ditch water table when the DCS were active (Figure 4). This effect was also observed in the field piezometers but was weaker, indicating that DWM had an impact on the groundwater storage, with the strongest effect close to the ditch. In contrast, the DWM did not affect the surface soil moisture within the agricultural fields. The observed effectiveness of DWM to reduce nitrate load in previous studies may be due to denitrification occurring in the ditches or in the elevated groundwater leading to the ditches [13,16,[44][45][46][47][48]. It is likely that the DWM reduced nitrate leaching into the ditches but did not significantly affect surface (5 cm) soil moisture and rates of denitrification within the surface soils of cropped fields between the ditches where the chamber flux measurements were conducted.

N 2 O and CH 4 Fluxes
Soil GHG emissions were not significantly impacted by DWM or the interaction of board period and treatment ( Figure 5 and Table 1). While DWM affected the depth of the water table, and preliminary results showed that hydrologic nitrate export was reduced [16], DWM did not impact the surface soil moisture or the GHG fluxes from the soil in the cropped areas ( Figure 4 and Table 1). Hence, the results do not support the hypothesis that DWM would increase CH 4 or N 2 O soil emissions. This null impact of DWM on GHG emissions confirms other literature findings, particularly those evaluating field scale DWM manipulations [30][31][32][33]. Other studies found some increase in CH 4 and N 2 O emissions in undrained fields, but were not evaluating the effect of DWM manipulations per se, but rather the effect of tile drainage in general, by comparing drained fields to undrained fields. In those instances, the undrained fields had higher N 2 O and CH 4 emissions, although the differences were seasonally variable and not consistent over multiple years [26,27].

Farm Management Impact on GHGs
For N 2 O fluxes, in both DWM and non-DWM fields, the mixed model results (Tables 1 and 2) demonstrated that N fertilization during corn years had the greatest impact on N 2 O emissions. Soon after each of the N fertilization events, we observed a pulsed increase in both soil nitrate concentrations and N 2 O emissions, which then declined over a period of 1-2 months as soil nitrate declined during crop development ( Figure 2D,E). This finding confirms existing literature that shows that N 2 O fluxes are highest in the days and weeks following N fertilization [42,[49][50][51].
For all N fertilization events in this study, fluxes were measured 1 day and again 3 or 4 days after. The largest fluxes were observed at 3-4 days after N fertilization. A study of corn/soybean in Brazil demonstrated a similar lag of a few days between fertilization and subsequent N 2 O flux spikes, which they attributed to the time it took for nitrification to convert fertilizer-N to nitrate [42].
The multiple regression results align with the published literature demonstrating that key variables of soil nitrate, temperature, and soil moisture play a role in controlling N 2 O emissions ( Figure 6; Equation (2)) [42,52]. The variation in these parameters fit within the context found for N-fertilized agricultural fields [53].
There was not a specific farm management action that most controlled CH 4 fluxes in this study. Seasonality explained most of the CH 4 flux trends. The highest CH 4 emissions occurred during the period from April to July in all three years of this study. This period corresponded with a time when the temperature had already increased above winter lows but WFPS still remained high (≥55%) due to spring rain and low rates of evapotranspiration before and for a few weeks following planting. The relatively warm and wet conditions were sufficiently conducive to methanogenesis within the soil to produce a period of somewhat elevated emissions relative to other seasons. No DWM effect occurred, and the spring seasonal effect occurred in both treatments.

Comparing Annual Soil Emission Estimates
For N 2 O, all three interpolation methods produced lower annual estimates in the 2018 soybean year compared to the corn years of 2017 and 2019. The higher emissions in the corn years were driven by N fertilization. There was an appearance of a possible difference in treatment means (albeit with overlapping standard errors) in the seasonally-based interpolation in 2017, but this was not observed in 2018 or 2019, nor in the other two annual estimate techniques ( Table 3). The non-parametric model and the mixed-effects model did not show significant differences in N 2 O gas fluxes between treatments in any year, which suggests that the interpolated treatment means in 2017 should not be considered different. The regression-based interpolation produced similar annual estimates to the mixed-effects model mean extrapolation, and both were somewhat lower than the seasonally-based interpolation. A lower annual estimate for the regression model is likely because non-linear responses to soil nitrate and lack of continuous nitrate data may cause underestimation of the N 2 O fluxes. Similarly, extrapolation of the means of the mixed-effects model gives more weight to mean responses rather than pulse (spike) events. Despite limitations for each of the three methods for estimating annual fluxes, they are in broad agreement and therefore represent the likely range of annual fluxes for each year and treatment. Because it captures the post-N fertilization N 2 O peaks, we used the seasonally-based interpolation estimate for comparison to the literature in Section 4.4.
For CH 4 , the seasonally-based annual interpolation had somewhat higher annual estimates compared to the mixed-effects model mean extrapolation, but all estimates indicated a very small net source, which is consistent with other studies of agricultural soils [54]. There were no differences in treatment in any year for either interpolation technique.

N 2 O Emissions Factor
The nitrous oxide emission factor (EF) was calculated by dividing the seasonally-based interpolated N 2 O annual flux, averaged over treatment and year, for the corn years reported in Table 3, by the average N fertilization rate of 213 kg N ha −1 . This resulted in an EF of 3.0%. This calculation method does not account for N 2 O emissions from a background, unfertilized agriculture site, likely leading to an EF overestimate. To partially account for this, the soybean N 2 O emissions in 2018 were used as an unfertilized control, and subtracted from the averaged N 2 O emissions in 2017 and 2019. This may overcompensate for the expected flux from an unfertilized soil control because the previous year's fertilization could have had a residual impact on soil N availability in 2018 and because biological N fixation likely added soil N during the soybean year. However, subtracting the entire annual flux from the soybean year provides a lower bracket estimate for the EF during the corn years. This correction yields an EF of 2.2%. An N 2 O EF in fertilized agriculture is influenced by a range of environmental variables such as soil type, N input type, and precipitation. The current IPCC guidelines use a default EF range of 1.3-1.9% for synthetic fertilizer agriculture in climates with greater than 1000 mm precipitation [55,56]. The EF estimate from this research is slightly higher than the upper limit of the uncertainty range of IPCC recommendation.

Nitrogen Balance
We calculated a nitrogen balance and nitrogen use efficiency (NUE) based upon average data from the two corn years. The average total amount of nitrogen applied during fertilized corn years was 213 kg N ha −1 . The average corn yield was 13,800 kg ha −1 (~220 bushels ac −1 ). The harvested corn had an average (across treatments) nitrogen content of 0.84 ± 0.18% (n = 8), resulting in removal of 91-141 kg N ha −1 , the range representing the variation in the measured nitrogen content. The nitrogen balance (also known as N surplus), defined as the amount of N applied minus the amount of N removed via harvest, was estimated to be 97 ± 20 kg N ha −1 . Using these estimates of N inputs and outputs, and assuming that annual N deposition rate at this site of 6 kg N ha −1 [57] had a negligible impact, and that N fixation was only applicable during soybean crop years, the average nitrogen use efficiency (NUE) for the corn years, defined as the N removed in the harvest divided by the N additions, was 54 ± 11%. This value is somewhat low for the U.S. but falls within the range of values calculated nationally and globally in irrigated maize systems [6,[58][59][60][61]. Ours is the first study known to us on N budget estimates from corn grown in the U.S. Mid-Atlantic region. Compared to the dataset compiled by [62] from the U.S. cornbelt, which supported a non-linear correlation between annual N 2 O emissions and N balance, the annual N 2 O emission estimates from the present study fall within the range reported by [62], but are higher than predicted based on their reported relationship with N balance ( Figure S6). This result is consistent with the lower than average NUE reported here.

Conclusions
The DWM condition did not increase soil N 2 O and CH 4 emissions. While there was evidence that DWM impacted the groundwater table near the drainage ditches, this impact was not evident in the surface soil moisture nor the surface soil GHG fluxes themselves. This finding spanned multiple years of corn and soybean crop rotation.
Three independent methods for deriving annual flux estimates revealed consistent estimates within their uncertainty ranges. Comparing the annual N 2 O flux estimates and related measures of N cycling to other studies on corn fields in the U.S., this farm had above average N surplus, below average NUE, and above average annual N 2 O fluxes, although these values fell within reported regional and global ranges.
Follow-up studies will address the possibility of N 2 O accumulation within the groundwater near the ditch areas and degassing of ditch water and will also report on the effect of DWM on export of dissolved N and P. For the cropland soil surface studied here, DWM did not have an unintended consequence of increased greenhouse gas emissions. This study provides evidence that DWM may be an effective BMP to reduce hydrologic nitrate export from crop fields without concerns of significant pollution swapping of increased greenhouse gas emissions from the soil surface.
Supplementary Materials: The following supporting information can be downloaded at: https:// www.mdpi.com/article/10.3390/nitrogen3010010/s1, Figure S1: Detail of DCS and DWM technique; Figure S2: Summary of diel variation of N 2 O fluxes; Figure S3: Visual of decay curve fitting for annual N 2 O estimate; Figure S4: Clay size fraction by depth and field; Figure S5: Summary of individual regression results for input variables into multiple linear regression; Figure S6: Plot of N balance versus N 2 O emissions; Table S1: N 2 O seasonally-based annual estimate details; Table S2: CH 4 seasonally-based annual estimate details; Table S3: USDA soil textural class description; Table S4: Soil characterization results; Table S5: Summary of minimal detectable flux calculations; Table S6: Summary results of multiple linear regression model. References [63][64][65][66]