Management and climate variability effects on understory productivity of forest and savanna ecosystems in Oklahoma, USA

. The productivity of herbaceous and understory woody vegetation is critical for wildlife habitat, livestock forage, and biodiversity, and is in ﬂ uenced by both annual weather patterns and tree dominance. With the goals to inform management and understand climate change implications, we determined the effects of tree harvest, prescribed ﬁ re, and 31 yr of climate variability on understory aboveground net primary productivity (ANPP) for ecosystems ranging from mature forest to grassland in a long-term study in southeastern Oklahoma, USA. In 1984, starting with a mature forest dominated by Pinus echinata and Quercus stellata , replicated experimental units were created by various combinations of pine harvest, hardwood thinning, and subsequent ﬁ re return intervals (1 – 4 yr and none). Understory ANPP (forbs, grasses, and woody plants) was measured by clip plots at the end of each growing season. Stepwise regression models were developed between understory ANPP and tree basal area, litter accumulation, ﬁ re return interval, and monthly, seasonal, and annual weather variables. Understory ANPP was dominated by grasses and ranged from 27 g (cid:1) m (cid:3) 2 (cid:1) yr (cid:3) 1 for the mature forest to 374 g (cid:1) m (cid:3) 2 (cid:1) yr (cid:3) 1 for an annually burned grassland/savanna. In general, herbaceous ANPP was inversely related to tree dominance (basal area), litter accumulation, and early and late growing season temperatures, and positively related to June precipitation. Understory woody ANPP was in ﬂ uenced by tree basal area, positively in ﬂ uenced by precipitation, and negatively in ﬂ uenced by summer soil moisture de ﬁ cits. Our results indicate that prescribed ﬁ re, through its negative in ﬂ uence on tree basal area and litter accumulation, is critical to maintaining highly productive understories and that a three-year return interval is a threshold to stall redevelopment of forest. For herbaceous ANPP, timing of precipitation, especially mid-growing season, appears more important than total precipitation, and higher temperatures within the range our site experienced did not have a large negative effect. In contrast, understory woody ANPP was negatively in ﬂ uenced by drought indicating climate change may have variable effects on different functional groups.


INTRODUCTION
Forest-savanna-grassland transition zones occur throughout the world where precipitation is marginal for tree growth. At this interface, shifts may occur between tree-dominated forest and herbaceous-dominated grassland based on disturbance, primarily fire and drought (Scholes and Archer 1997, Loveland et al. 2000, Oliveras and Malhi 2016. At intermediate precipitation levels and without strongly differentiated seasonal distribution, tree cover may be either high or low with fire frequency being the determinant between savanna and forest (Staver et al. 2011). The forest-grassland transition zone in the south-central United States historically encompassed 407,000 km 2 and includes north-central Texas, central Oklahoma, eastern Kansas, northern Missouri, and much of Illinois (https://www. worldwildlife.org/ecoregions/na0804). The region is characterized by its mixture of savanna, prairie, and woodlands.
Forest-grassland transition zones comprise significant variability in energy flow, soil moisture, hydrology, and biogeochemical cycling (Breshears 2006). The relative abundance of trees vs. herbaceous vegetation strongly influences the microclimate, ecological processes, fire regime, and fauna. Tree canopy cover and basal area are the most often used structural characteristics for defining grassland, savanna, woodlands, and forests, yet definitions vary somewhat in the scientific literature. For consistency, we define grassland as <10% tree canopy cover, savannas 10-30% tree canopy cover, woodlands 30-80% tree canopy cover, and forests >80% tree canopy cover (Dey et al. 2017).
The interplay between fire regimes and climate variability causes considerable temporal and spatial variation in species composition, tree abundance, biomass, and net primary productivity within forest-grassland transition zones (Lehmann et al. 2014, Scheiter et al. 2015. This is particularly true for understory vegetation, that is, herbaceous vegetation and woody vegetation shorter than 1 m, which is regulated in large part by tree density (Masters et al. 1996, Peterson et al. 2007, Zou et al. 2007, Feltrin et al. 2016) in addition to direct effects of climate variability and fire (Masters et al. 1993, Briggs and Knapp 1995, Heisler et al. 2003. Understanding processes related to understory vegetation is critical because understory vegetation is an important component of ecosystem productivity, vital for energy flow and nutrient cycling (Ojima et al. 1994, Wise and Schaefer 1994, Lapointe 2001, and essential for the maintenance of biodiversity, wildlife habitat, forage for domestic animals, and habitat for pollinators and other insects (Masters et al. 1996, Oliveras andMalhi 2016).
Within forest-grassland transition zones, tree mortality associated with drought and other disturbance causes rapid shifts in ecosystem composition, structure, and functioning and increases the prevalence of savanna and grasslands (Albertson and Weaver 1945, Rice and Penfound 1959, Balch et al. 2015. Climate change may reduce productivity and cause an eastward shift in the forest-grassland transition zone of the southcentral United States (Notaro 2008, Seager et al. 2018. The average temperature throughout the south-central United States is expected to increase by 2.5-4.0°C by the latter half of this century (Collins et al. 2013). Larger and more intense precipitation events along with increased duration and intensity of drought are projected, particularly in Oklahoma and north-central Texas (Abatzoglou and Brown 2012, Collins et al. 2013, Qiao et al. 2017, Wuebbles et al. 2017. Given that precipitation is a strong determinant of plant growth, projected drought frequency and intensity as well as large rainfall events triggered by climate change will likely impact productivity within forestgrassland transition zone (Knapp and Smith 2001, Fay et al. 2008, Wu et al. 2018. In addition to climate effects, frequent burning (<3-yr intervals) shifts closed-canopy forest ecosystems toward open woodlands and savannas and increases understory aboveground net primary productivity (ANPP) in part by killing or top-killing fire-intolerant trees, reducing tree growth rates, and slowing recruitment of new trees Reich 2001, Moore et al. 2016). Fire also has a direct effect on understory ANPP by removing litter, releasing and mobilizing nutrients, and encouraging growth of nitrogenfixing legumes Seastedt 1986, Ojima et al. 1994). In contrast, reduction in fire frequency results in greater tree canopy cover (Bond et al. 2005), which limits the abundance and productivity of understory herbaceous vegetation and shifts fuels from grass-dominated to tree leaf litter, but can increase diversity and productivity of understory woody plants Archer 1989, Feltrin et al. 2016).
Management to reduce tree abundance, that is, mechanical, herbicide, or both, along with prescribed fire can be used to create and maintain woodland and savanna ecosystems (Masters and Waymire 2012). Within forest-grassland transition zones, climate change and increased water stress might increase tree mortality and increase the abundance of more open ecosystems (Albertson andWeaver 1945, Balch et al. 2015). A critical unknown is the interaction between management and climate variation and how it might affect the future productivity of the herbaceous and woody understory. While modeling efforts have examined the potential effects of climate change drivers on ecosystems and shifts among ecosystems in South Africa and elsewhere (Woodward et al. 2004, Bond et al. 2005, Bond 2008, Scheiter et al. 2015, the majority of empirical studies considering global change impacts within forest-grassland transition zones have focused on individual ecosystems (Smoliak 1986, Paruelo et al. 1999, Oesterheld et al. 2001, Nippert et al. 2006, Patton et al. 2007, Peterson et al. 2007, Hsu and Adler 2014 and studies examining potential shifts among ecosystems due to the interactions of climate and management are scarce (Sala and Maestre 2014).
The overall objective of our study was to identify the major factors that influence the interannual variation in understory ANPP of different ecosystem structures within a forest-grassland transition zone in southeastern Oklahoma, USA. Different ecosystem structures were created in 1984 by thinning (a combination of commercial pine tree harvest and thinning of hardwoods by herbicide injection), and subsequent introduction of fire at varying frequency (annual to four-year fire return intervals and fire exclusion). We explored the effects of interannual and intraannual variation in weather, time since fire, and ecological factors related to tree density. Understanding how past climate variability influenced understory productivity will help inform the future response of climate change by comparing periods of drought and above-average temperature to productivity during less stressful periods. Knowledge of how fire frequency and tree density affect understory productivity is critical to establish management regimes and overstory target densities to achieve objectives related to ecological restoration, forage production, and wildlife habitat as well as to better understand ecosystem carbon dynamics.
Several expectations guided our work: (1) Across a forest-grassland continuum, climate variability, fire, and tree canopy cover combine to affect woody and herbaceous ANPP, (2) herbaceous ANPP increases and understory woody ANPP decreases with increased fire frequency because more frequent fire increases canopy openness and reduces litter, (3) herbaceous ANPP increases in years immediately after fire and declines as time since fire increases and leaf litter accumulates, and (4) herbaceous ANPP is correlated with timing of precipitation, particularly early growing season. In contrast, multimonth summer drought has a strong negative effect on ANPP of understory woody vegetation.

Study area
The study was conducted on a 53-ha area (34°31 0 40″ N, 95°21 0 10″ W) located within the 7690-ha Pushmataha Wildlife Management Area (WMA), owned and administered by the Oklahoma Department of Wildlife Conservation (Fig. 1). The Pushmataha WMA is located in the Kiamichi Mountains along the western edge of the Ouachita Mountains in southeastern Oklahoma. Soils are derived from sandstones and shales and belong to the Carnasaw (fine, mixed, semi-active, thermic Typic Hapludults) and Stapp (fine, mixed, active, thermic Aquic Hapludults) series (NRCS 2019).
The climate of Pushmataha WMA is characterized as semi-humid with hot summers and moderate winters. The 31-yr average  annual minimum and maximum temperatures of the study area were 7.5-22.2°C (Thornton et al. 2018). During this same time period, the study area received average annual rainfall of 1436 mm (Thornton et al. 2018).
The study area was initially dominated by post oak (Quercus stellata), shortleaf pine (Pinus echinata), blackjack oak (Q. marilandica), and hickory (Carya spp.). Pretreatment understory condition was primarily leaf litter with occasional herbaceous plants. However, the understory woody vegetation was dominated by sparkleberry (Vaccinium arboreum), winged sumac (Rhus copallinum), poison ivy (Toxicodendron radicans), and greenbriers (Smilax spp.) as well as seedlings of the overstory tree species. Posttreatment groundlevel vegetation was mainly composed of C 4 grasses, little bluestem (Schizachyrium scoparium), big bluestem (Andropogon gerardii), and Indian grass (Sorghastrum nutans) and to a lesser extent v www.esajournals.org panicum grasses (Panicum spp., Dichanthelium spp.) and sedge (Carex spp.). Posttreatment woody vegetation most frequently encountered was post oak, blackjack oak, and shortleaf pine resprouts or saplings, winged sumac, and greenbriers. Before treatment, the study area had been logged, grazed, and frequently burned until acquisition in the mid-1940s, and later transitioned to a closed-canopy forest. Based on ages of the shortleaf pine determined from cores, the last major disturbance was logging in the 1930s.

Treatments
In summer of 1983, 28 management units ranging between 0.8 and 1.6 ha were established with a randomized experimental design across homogenous closed-canopy forest ( Fig. 1; Masters and Waymire 2012). During the summer of 1984, treatments were applied to each unit so that 23 units represented eight cultural treatments with three replications of each except the treatment (HT3), which had only two replicates ( Fig. 1). The treatments consisted of combinations of harvesting shortleaf pine trees with diameter at breast height (dbh) 11.4 cm or larger (H), thinning of hardwoods to approximately 9 m 2 /ha basal area using single-stem injection of herbicide (T), and fire return intervals (1-4 yr and fire exclusion; Table 1). Prescribed fires were conducted from January to early April with most occurring in March, which is considered late dormant season. The specified fire frequency treatments were maintained throughout the experiment. Prescribed fires ranged in fire line intensity from 44 to 5150 kW/m 2 . By 1987, when we began to collect understory biomass data, the experimental units that had been thinned and burned were grassland and savanna (Table 1).

Sampling
Within each treatment unit, ten permanent 4 9 4 m plots were established at 20-m intervals along two randomly located parallel transects perpendicular to the contour. To measure ANPP  Table 1 for the description of treatments (Source: Google Earth Pro 2019).
v www.esajournals.org of understory vegetation, we sampled 0.25-m 2 subplots adjacent to each of the 10 permanent plots in each of 17 units. A different area was randomly sampled every year outside the 4 9 4 m plots, but was consistent among all plots. From 1987 to 2017, sampling was conducted in September or October before first frost at the time of peak standing biomass. We quantified the understory ANPP by clipping aboveground vegetation and separating it into different functional groups (understory woody, grass, nonlegume forbs [hereafter forb], sedge, and legume). For woody understory vegetation, we collected annual growth (leaves plus currentyear shoots) for stems shorter than 1.4 m in height. In addition, we collected litter (dead understory plant materials, bark, leaves, and branches <2.5 cm diameter) from the same subplots. The samples were kept separate and dried at 60°C until constant weight and weighed. During the 31-yr interval, we did not sample in 1991, 1992, 1996, 2002, 2004, and 2009. We also measured basal area and canopy cover of overstory trees annually between September and October. Basal area was measured with a 10factor wedge prism for trees ≥5 cm dbh with the plot centered in each permanent plot. Overstory canopy cover was estimated initially using a 5point grid and later switched to a 9-point grid in a sighting tube with vertical and horizontal levels. Estimates were made at plot center and cardinal points at 2 and 4 m from each plot location, for nine sample points (Mueller-Dombois and Ellenberg 1974).
Using heteroscedasticity-adjusted ANOVA with Tukey's multiple comparisons, we tested for differences between various treatments for the following variables: (1) ANPP of understory herbaceous vegetation, which included forbs, grass, legumes, and sedge, (2) ANPP of understory woody vegetation, (3) ANPP of all understory vegetation (sum of above-mentioned two categories), (4) litter, (5) total basal area, and (6) canopy cover. As appropriate, we tested for the effect of number of growing seasons since most recent prescribed fire. We used the term growing years since fire (GYSF), where 1 represents data collected at the end of the first growing season after burning, 2 at the end of the second growing season, 3 at the end of the third growing season, and 4 at the end of the fourth growing season. We used the average for each treatment unit, as the experimental unit in these analyses, that is, treated the 10 plots per unit as subsamples. Notes: Hardwood thinning and shortleaf pine (Pinus echinata) harvesting were completed during the establishment of the experimental site at the Pushmataha Wildlife Management Area, Oklahoma, USA, during 1984. Grassland <10% canopy cover; savanna 10-30% canopy cover; woodland 30-80% canopy cover; forest >80% canopy cover.

Relationship between understory ANPP and weather and ecological variables
We related understory ANPP (herbaceous, woody, herbaceous + woody) to (1) ecological variables, which included litter, overstory basal area (total, pine only, hardwood only), and canopy cover, (2) weather variables, which included various measures of temperature, precipitation, and drought, and (3) fire management, that is, GYSF. Specifically, we included monthly, seasonal, and annual intervals for average, average maximum, and average minimum temperature and monthly, seasonal, and annual precipitation for time periods within 1987-2017. The monthly values for climate data were computed from daily data obtained from Daymet (Thornton et al. 2018) for all months. We supplemented the temperature and precipitation variables with the Keetch-Byram Drought Index (KBDI). The KBDI estimates the dryness of the soil for a wide range of climatic and rainfall conditions (Keetch and Byram 1968). The KBDI value ranges from 0 to 800; a value of 0 represents a complete saturation of the soil, while a value of 800 represents an absolutely dry soil. Keetch-Byram Drought Index was estimated using daily precipitation and temperature following an equation described by Dolling et al. (2005).
The aforementioned variables include a total of 89 potential predictors of understory ANPP (Appendix S1). To minimize the burden of an over-parameterized model for predicting each of the three understory ANPP components (herbaceous, woody, and total), we first performed a simple linear regression between each potential predictor and the selected response variable of understory ANPP. This was performed separately within each treatment and by aggregating data across treatments. In addition to testing relationships between current-year variables and understory ANPP, we also tested for relationships with previous year's variables. Based on this initial analysis, we dropped the predictor variables that were not significantly correlated (P > 0.05) with the dependent variable of interest. For these analyses, we used the mean annual response for ANPP of each treatment, that is, averaged ANPP for the two or three replicates of each treatment for each year.
Starting with the significant predictors based on bivariate regression, we built a series of multiple regression models with bidirectional elimination of parameters (stepwise regression). The competing models were evaluated based on the Akaike information criterion adjusted for small sample size (Buckland et al. 1997, Calcagno andde Mazancourt 2010). At each step of bidirectional elimination, a new predictor variable was included in the model only if the partial R 2 associated with the variable was statistically significant and provided unique explanatory power to explain the response beyond its correlation with other predictors. The models for the various treatments occasionally included different but closely correlated variables, for example, July maximum temperature vs. July average temperature. To simplify the presentation of results, we eliminated one of these variables and reran the analysis, keeping the new result only if the total R 2 was not reduced upon recalculation (Data S1).

Climate
Between 1987 and 2017, annual precipitation (based on the October-September water-year) varied from 971 to 1935 mm with a mean of 1436 mm and below average precipitation for 17 yr (Fig. 2a). The study site experienced two periods of extended drought during the course of this study, from 2003 to 2006 and from 2009 to 2013. The annual pattern based on monthly averages was for greater precipitation in the early growing season followed by hot and dry summers resulting in late growing season water deficits (Fig. 2b). As a result, the site typically experienced seasonally dry soils from July to October (maximum average KBDI in August and September). Given the cooler temperatures and increased precipitation in winter and spring, soil moisture recharged over the winter (minimum average KBDI in February and March). Regression analysis showed no statistical trends of change in annual precipitation or temperature over time (P < 0.05).

Ecosystem structure
Application of thinning, harvesting, and fire altered the basal area and canopy cover of v www.esajournals.org overstory trees across the experimental site. Before treatment, the entire site was similar to the control, which was initially classified as woodland or forest (Table 1). In 2017, the control treatment consisted of closed-canopy, mature forest dominated by larger diameter pine trees (23 cm diameter at breast height) and mix of larger and medium-sized hardwood trees (Table 1). The oldest shortleaf pines measured were approximately 94 yr old as estimated from tree cores. When measured at the end of the 1985 growing season, one year after harvesting, and the year of the first burn, canopy cover of individual units ranged from 0.2% to 28.2% for the harvested, thinned, and burned treatments and units were classified as grassland or savanna ( Table 1). The HNT1 treatment tended to have higher canopy cover due to the greater presence of residual hardwoods and was classified as savanna or woodland. By 2017, burning every four years resulted in the transition to an uneven-aged woodland due to growth of residual trees and periodic recruitment of new trees (canopy cover increased from 7.3% to 52.4%; Table 1). In contrast, burning every 1-3 yr resulted in smaller increases in canopy cover, mostly due to residual tree crown expansion, resulting in grassland, savanna, or woodland (a) Average annual precipitation and average temperature over 31 yr of the study period and (b) average monthly temperature, precipitation, and Keetch-Byram Drought Index (KBDI). The dotted black line represents average precipitation during study period. Temperature and precipitation were computed from daily data obtained from Daymet (Thornton et al. 2018), and KBDI was calculated using temperature and precipitation (Keetch and Byram 1968

Understory production
When comparing the ANPP of different functional groups (i.e., sedge, legume, forbs, grass, and woody) after 31 yr of treatment in 2017 (Fig. 3), the proportion of herbaceous biomass was greatest in HNT1 (99%) and the least in control (47%). Total understory ANPP ranged from 104 gÁm À2 Áyr À1 in the control to 747 gÁm À2 Áyr À1 in the HT1 treatments. Woody understory ANPP was significantly lower in the annually burned treatments compared with the other treatments. Mean grass ANPP was significantly lower in control than all other treatments, lower in the HT4 treatment than the other burned treatments, and greater in the HT1 treatment among all treatments.
Mean understory ANPP in control treatment was consistently low throughout the years (Fig. 4, Table 2). The maximum herbaceous ANPP for the HT1 treatment was 556 g/m 2 in 1999, and minimum was 182 g/m 2 in 2005. The other treatments exhibited similar large variation in productivity over time as herbaceous ANPP responded to differences in weather and timing of fires. One notable trend was the decrease in herbaceous ANPP for the HT4 treatment that occurred as the uneven-aged woodland developed. The trends for woody and herbaceous ANPP over time were not similar, and total understory ANPP was driven mostly by changes in the larger herbaceous component (Fig. 4). When averaged across the 31-yr study period for the entire study area, herbaceous plants contributed 76% to total understory ANPP. The contribution of herbaceous plants to total understory ANPP increased with increasing fire frequency, that is, control (54%), HT4 (70%), HT3 (73%), HT2 (76%), HT1 (93%), and HNT1 (92%). Among the functional groups in the herbaceous category, grasses made the largest contribution to herbaceous ANPP (89% on average).
Averaged throughout the study period, mean understory ANPP was 255.6 gÁm À2 Áyr À1 across the study site (Table 2) and varied between 27.4 gÁm À2 Áyr À1 for the closed-canopy mature forest (control) and 374.3 gÁm À2 Áyr À1 for annually burned savanna (HT1; Table 2). The inclusion of  Table 1 for the description of treatments.
v www.esajournals.org Fig. 4. Interannual trends in mean understory, herbaceous, and woody aboveground net primary productivity (ANPP) in different treatments over time . See Table 1 for the description of treatments.
fire significantly increased (P < 0.05) total ANPP, and the HT1 was significantly greater than the other treatments including fire (Table 2). Mean herbaceous ANPP of study site was 206.8 gÁm À2 Áyr À1 over the experiment (Table 2). Differences in herbaceous ANPP were similar as for total ANPP with the exception that ANPP of the HT4 treatment was lower than the other burned treatments. Mean woody ANPP across the study site averaged 48.8 gÁm À2 Áyr À1 (Table 2). Mean woody ANPP was significantly lower in the control, HT1, and HNT1 treatments compared with the HT2, HT3, and HT4 treatments.

Climate and ecological controls of understory primary production
Understory ANPP was correlated with many of the ecological (Table 3) and weather variables (Data S1). Because some of the ecological variables (Table 3) and weather variables (Appendix S2) were correlated with each other, different but often related variables were included in multiple regressions depending on the model (Tables 4-6). Across all treatments, herbaceous ANPP was the dominant component (76% on average) with grass ANPP composing 89% of the herbaceous category. Therefore, when total understory ANPP (herbaceous + woody) was considered, model results largely reflected the herbaceous component, but generally had poorer fit because the woody and herbaceous components responded differently (Tables 4-6). Therefore, we focus on the results of the herbaceous and woody vegetation separately.
In general, herbaceous ANPP was inversely related to tree dominance (basal area; Fig. 5a) Table 2. Mean value AE standard error of ANPP for total understory, herbaceous, woody, and grass (gÁm À2 Áyr À1 ), total basal area (m 2 /ha), canopy cover (%), and litter biomass (g/m 2 ) across study area and for different treatments averaged over time   v www.esajournals.org and litter accumulation (Fig. 5b) and positively related to growing season precipitation and early and late season temperatures (Data S1). Total R 2 for the multiple regression including all treatments was 0.63 and ranged between 0.41 and 0.81 for individual treatments (Table 4). Considering all treatments together, which span conditions from grassland to closed-canopy forest, herbaceous ANPP was negatively related to litter accumulation (controlled by both tree inputs and fire frequency) and total basal area. The herbaceous ANPP of the control was negatively related to KBDI in September (late season drought; Fig. 6) and total basal area as well as positively related to October temperature (Table 5). In contrast, the herbaceous ANPP of the HT4 treatment was primarily influenced by litter accumulation (inverse relationship) and also related to total basal area, precipitation in February, and September temperature ( Table 5). The only factor related to herbaceous ANPP for the HT3 was precipitation in September (also treatment with the lowest R 2 ) (0.15, P = 0.052). The herbaceous ANPP of the HT2 treatment was negatively influenced by total basal area and GYSF because herbaceous ANPP was greater the year immediately after fire as compared to the second year after fire. The HT1 and HNT1 treatments were burned annually, which eliminated variation for factors such as GYSF and litter. In these two treatments, herbaceous ANPP was most influenced by a positive relationship with June precipitation (Figs. 7a, b) and February temperature (Table 5). While woody understory ANPP was only 24% of the total, it was more important in the forest treatments (46% of control and 30% of HT4). Total R 2 of the model of different treatments ranged 0.41-0.72 (Table 6). In general, understory woody ANPP was influenced by tree basal area, positively influenced by precipitation, and negatively influenced by KBDI in July and November temperature. When considering all treatments together, the model for understory woody ANPP  v www.esajournals.org was relatively poor (R 2 = 0.30) because individual treatments responded to different factors. Understory woody ANPP of the control treatment was strongly influenced by pine basal area and August precipitation, and negatively related to July KDBI. Woody understory ANPP of the HT4 treatment was most strongly influenced by June temperature and June precipitation. Understory woody understory ANPP of the HT3 treatment was positively influenced by hardwood basal area and negatively influenced by July KDBI. The woody understory ANPP of the HT2 treatment was positively influenced by litter and June/ August precipitation. The woody understory ANPP of the HT1 treatment was negatively influenced by hardwood basal area and November temperature. The woody understory ANPP of the HNT1 treatment was negatively influenced by June and November temperature (Table 6).

Relationships with the most influential variables
Basal area and litter.-Considering all treatments together, the only individual variables correlated (negative relationships) with herbaceous ANPP were variables associated with basal area, canopy cover, and litter accumulation (Data S1). Both basal area and canopy cover are measures of tree dominance and were highly correlated with one another (R = 0.98). Total basal area was correlated with herbaceous ANPP for all individual treatments (Data S1) and was included as a variable for the majority of stepwise regression models (Table 4). The basal area of the control was consistent throughout the experiment, averaging 27.1 m 2 /ha (Fig. 5c, Table 2). The basal area of the HT4 treatment increased from 3.1 to 20.7 m 2 /ha during the experiment as an unevenaged woodland developed. In contrast, basal areas of the other treatments were consistently low (<10 m 2 /ha). Comparing just the effects of basal area on herbaceous ANPP resulted in a linear relationship (herbaceous ANPP [g/ m 2 ] = 306.75 À 10.83 9 basal area [m 2 /ha] R 2 = 0.58; Fig. 5a). The result is that for every 10 m 2 /ha increase in basal area, herbaceous ANPP decreased by 108.3 gÁm À2 Áyr À1 . Woody understory ANPP of the various treatments was not as strongly correlated with total basal areas (R 2 < 0.20; Data S1).
Litter accumulation was primarily a result of woody canopy cover, residual dead grass from previous years, and longer fire return intervals. Litter was consistently high in the control treatment, consistently low in the annually burned treatments, and more variable in the periodically burned treatments (Fig. 5d). When all treatments were considered together, herbaceous ANPP was strongly and nonlinearly correlated with litter such that herbaceous ANPP showed a sharp decrease as litter biomass increased between 0 and 500 g/m 2 , but much less of a decrease at higher litter accumulations, approaching 0 gÁm À2 Áyr À1 herbaceous ANPP at approximately 2000 g/m 2 litter (herbaceous ANPP [g/ m 2 ] = 402.70 9 exp [0.01 9 litter (g/m 2 ) À 49.2, R 2 = 0.62]; Fig. 5b). In contrast to the herbaceous ANPP, woody understory ANPP was not correlated with litter (Data S1).
Precipitation.-Some variables related to precipitation and drought were included for most individual treatment models. However, the various other influencing factors related to management make it hard to isolate these effects. The HT1 and HNT1 treatments provide the best opportunity to isolate the effects of precipitation as litter accumulation and periodicity of fire are constant year to year. For these two treatments, herbaceous ANPP was best related to June precipitation (Figs. 7a, b). The slopes were 0.78 and 0.56 gÁm À2 herbaceous ANPPÁmm À1 precipitation for the HT1 and HNT1 treatments, respectively. The HT1 treatment produced more herbaceous ANPP than the HNT1 treatment because of the presence of large post oak trees in the HNT1. The net effect was that for every 100mm increase in June precipitation, herbaceous ANPP increased by approximately 78 and 56 gÁm À2 Áyr À1 for HT1 and HNT1, respectively. Woody understory ANPP also was positively correlated with June precipitation for several treatments (Data S1), but June precipitation was not included in the multiple regression because other variables related to a wet spring explained a greater amount of variation.
Drought.-September KBDI had a strong influence on herbaceous ANPP in the control treatment (Fig. 6). The slope of À0.04 indicates that for every 100 unit change in KBDI, herbaceous ANPP decreased by 4 gÁm À2 Áyr À1 . The relationship was significant with 95% confidence intervals for the slope of the regression line ranged from À0.06 to À0.03 (P < 0.0003). In contrast to the control treatment, KBDI was not as strongly correlated with herbaceous ANPP for the other treatments with maximum R 2 all less than 0.23 for simple regressions (Data S1). For woody ANPP, July KBDI was included in several of the stepwise regression models (Table 6). For other months, KBDI was occasionally correlated with woody understory ANPP, but the R 2 values were generally low and the correlations lacked a clear pattern (Data S1).
Fire. -Herbaceous ANPP decreased with time since fire for the HT2 and the HT4 treatments, but not for the HT3 treatment (Fig. 8). Herbaceous ANPP decreased from 273.5 to 192.5 gÁm À2 Áyr À1 in the HT2 treatment from the first year to the second year after fire. For the HT4 treatment, herbaceous ANPP progressively decreased from 235.3 to 167.6, 143.0, and 118.2 gÁm À2 Áyr À1 after one, two, three, and four years of fire, respectively (Fig. 8). Woody understory ANPP did not significantly vary with GYSF.

DISCUSSION
When investigated over 31 yr, fire, climate variability, and management all were important drivers of ecosystem structure, composition, and understory productivity within the forest-grassland transition zone of southeastern Oklahoma. While both herbaceous and woody ANPP were affected by ecological factors and weather, their responses differed. As expected, we found that herbaceous ANPP increased with shorter fire return intervals because fire maintained low tree basal area and canopy cover and consumed litter. Also as expected, herbaceous ANPP was greater the year immediately after fire for treatments burned every two or four years, likely related to litter consumption and possibly nutrient release. Herbaceous ANPP was most sensitive to midgrowing season precipitation, that is, June. In contrast, annual burning reduced understory woody ANPP (<1.4 m tall), but longer fire return intervals and time since fire did not have an effect. Unlike herbaceous ANPP, understory woody ANPP was not consistently related to overstory basal area and was positively related to litter accumulation. Understory woody ANPP was more influenced by longer-term patterns of precipitation and soil moisture rather than the timing of precipitation.
At this study site, relatively fast recovery of perennial grasses following treatment in the mid-1980s was in large part due to existing rhizomes, which were a legacy of the savanna and woodland ecosystems that prevailed before fire exclusion and subsequent forest development (as observed by Ronald Masters). Burning every three years was an ecological threshold that forestalled redevelopment of forest conditions. In contrast, the four-year fire return interval was long enough that resprouting trees developed sufficiently thick basal bark or a protective whorl of sprouts to insulate the dominant stem from top kill. The longer fire return interval also caused reduction in fire intensity. Burning at intervals of three years or less maintained sufficient C 4 grass fuels to increase fire line intensity v www.esajournals.org Fig. 5. (a) The relationship between basal area and average herbaceous aboveground net primary productivity above 475 kW/m 2 thus top-killing small woody stems and saplings >1 m tall (Sparks et al. 1999).
The effects of fire and tree dominance, that is, basal area and canopy cover, at our site were inter-related, and both had direct and indirect effects on understory ANPP. Increasing tree cover reduces light availability for the understory (Feltrin et al. 2016), increases competition for water and nutrients (Eastham et al. 1990, Kumar et al. 2010, increases litter production (Facelli andPickett 1991, Feltrin et al. 2016), and changes fire behavior (Baker and Bunyavejchewin 2009). Increased fire frequency reduces tree dominance, consumes litter, and increases nutrient availability (Moore et al. 2016, Ficken andWright 2017).
When considering all treatments together, about two-thirds of the variation in herbaceous ANPP was explained by litter accumulation and tree basal area. While increasing shade cast by the overstory undoubtedly was a factor, the effects of fire on litter reduction were likely as important. At this same research site, Feltrin et al. (2016) found a downward shift in the relationship between available light and understory productivity for forest compared with savanna ecosystems, which in part was likely due to the greater litter accumulation and reduced frequency of fire in the forest ecosystems. Litter suppresses herbaceous ANPP by preventing germination, especially for species that require mineral soil for germination, casting shade at the soil surface, changing microclimate, and adding allelopathic chemicals (Sydes and Grime 1981, Facelli and Pickett 1991, Xiong and Nilsson 1999, Hiers et al. 2007. A study of frequently burned Pinus palustris ecosystems found that tree litter had a stronger influence on herbaceous ANPP than did the shade cast by mid-story oak trees (Hiers et al. 2014). From a management standpoint, this emphasizes the need for periodic fire to reduce litter accumulation regardless of overstory density. (ANPP) over time , (b) the relationship between litter accumulation and herbaceous ANPP over the study period, (c) trends in total basal area change, and (d) trends in litter accumulation, across different treatments during study period. (Fig. 5. Continued) v www.esajournals.org Considering individual treatments, herbaceous ANPP was negatively influenced by basal area or litter except the annually burned treatments and the HT3 treatment (no significant relationships). Lack of basal area and litter influence in the HT1 and HNT1 was expected since these treatments have low tree basal area and fire removes litter on an annual basis. The lack of statistical significance for the relationship between herbaceous ANPP and predictor variables in the HT3 treatment could be due to the HT3 treatment only having two replications. The interrelationship between fire and litter discussed above was apparent when considering the HT2 and HT4 treatments. Although ANPP of both treatments was correlated with GYSF (Data S1), the stepwise regression models included GYSF for HT2 and litter for HT4 indicating that litter and GYSF explained much of the same variation in understory ANPP response.
Herbaceous ANPP of the annually burned treatments was most influenced by June precipitation. In our study area, a disproportionate amount of annual precipitation occurs in spring (33%) and hot, dry summers are the norm (23% of annual precipitation with 33°C average temperature) such that most growing seasons begin with ample soil moisture and end up with v www.esajournals.org summer water deficits (Fig 2). As a consequence, August on average experiences maximum soil water deficit (KBDI = 507) and the largest increase in KBDI happens between June and July. Therefore, June is typically the transition between high soil moisture and low soil moisture such that above average June rainfall likely extends the period of ample soil moisture for grass growth, specifically C 4 grasses, as many of the tallgrass prairie species bloom and produce seed later in the summer. While precipitation positively influenced herbaceous ANPP in savanna ecosystems of our study site, similar to previous studies (Oesterheld et al. 2001, Nippert et al. 2006, Kanniah et al. 2013, late season water deficit, that is, September KBDI, had a strong negative influence on herbaceous ANPP of the forested control treatment (see 11th paragraph of Discussion). The positive relationships between herbaceous ANPP and early or late season temperatures likely were related to extension of the growing season (Reyes-Fox et al. 2016. Understory ANPP can be influenced by legacy of the previous year's (or longer) environmental conditions (Lauenroth and Sala 1992, Wiegand et al. 2004, Monger et al. 2015. While we tested for the effects of carryover from previous years, we did not find any that were statistically significant. This is different from other studies where previous-year conditions accounted for part of the variation in understory ANPP such that dry legacy years reduced ANPP and wet legacy years increased ANPP (Smoliak 1986, Jobb agy and Sala 2000, Oesterheld et al. 2001, Sala et al. 2012, Sternberg et al. 2017. Since the magnitude of legacy is related to differences between previousyear precipitation and current-year precipitation (Reichmann et al. 2013), it is possible that the relative annual differences in precipitation at our study sites were not enough to produce significant lag effects. Our study site receives rainfall between 1000 and 1900 mm (annual averagẽ 1400 mm during the experiment). In contrast, studies showing a strong legacy on understory ANPP typically had lower annual precipitation (Lauenroth and Sala 1992).
As opposed to herbaceous ANPP, which was negatively influenced by litter and tree basal area, understory woody ANPP was positively related to litter (overall and HT2 treatment) and inconsistently influenced by basal area. While all understory woody plants in the frequently burned treatments appeared to be resprouts following top kill, they differed in time of establishment and spatial pattern. Most of the hardwood understory tree stems, such as post oak, likely have been resprouting since initiation of treatments, while shrub species that were abundant in the HT2 and HT3 treatments, such as sumac, and shortleaf pine saplings in the HT4 treatment, likely established later. Initial spatial pattern of hardwoods combined with some infilling of shrubs and saplings may have caused some of the differences in woody ANPP response to overstory basal area among treatments.
In the control treatment, understory woody ANPP was positively related to pine basal area. In this closed-canopy forest, there may be more light availability under the large pine trees because their canopies are generally higher and because foliar arrangement of pines allows greater passage of light (McGuire et al. 2001, Battaglia et al. 2003. The positive relationship between hardwood basal area and woody ANPP for the HT3 treatment may be related to proximity of seed source, that is, greater regeneration near overstory trees, as a three-year burning regime may allow occasional seedling and sapling recruitment over time. Alternatively, a change in fuels and fire behavior under hardwood trees in the HT3 treatment, as opposed to grass litter, may reduce fire intensity and facilitate seedling and sapling establishment. Likewise, the positive relationship between woody ANPP and litter for the HT2 treatment and for all treatments combined may be due to greater tree litter accumulation and change in fuel and fire behavior directly under trees (correlation between hardwood basal area and litter R = 0.73).
Like the herbaceous ANPP, woody understory ANPP was positively related to precipitation, but the month included in the models varied by treatments, likely indicating a general benefit of greater precipitation. While precipitation and drought are obviously related, drought had a stronger influence on woody than herbaceous ANPP as woody understory ANPP decreased in relation to July soil moisture deficit (KBDI) for all treatments combined and for the control and HT3 treatments. Also, the HT4 treatment was correlated with July KDBI (Data S1). Deeper rooting of woody vegetation might buffer fluctuations in soil moisture near the soil surface such that more extreme reductions in soil moisture throughout the soil profile are needed to influence woody understory ANPP.
One of our goals was to use the past responses of understory ANPP to higher than average temperature and drought to inform the potential impacts of climate change. A few patterns regarding the effects of temperature and precipitation on understory ANPP stand out. When removing the confounding influences of time since fire and mostly eliminating the influences of tree overstory, June precipitation had the strongest influence on herbaceous ANPP, which was mainly composed of grasses. This indicates that timing of rainfall is perhaps more important than total rainfall. Recent trends in the southcentral United States indicate greater spring rainfall (Oklahoma Mesonet 2019), which might increase productivity of grass-dominated systems even if total rainfall does not change or decreases. However, in the forested control treatment, herbaceous ANPP might be most affected by late season soil moisture given the fivefold decrease in herbaceous ANPP we measured with increasing September KBDI. The differences between the forest (control) and savanna/woodland ecosystems (HT1 and HNT1) could be related to different understory species. For instance, the grass component of the control treatment was 20% rosette panicum grasses (Panicum spp., Dichanthelium spp.), while the other treatments had less than 5%. Another possibility is that the trees exerted greater competition for water (Ludwig et al. 2004) or the tree canopies and deep litter layers intercepted more precipitation (Zou et al. 2015), which worsened the effects of late season soil moisture deficit on herbaceous ANPP in the control treatment.
While higher temperatures contribute to drought due in large part to greater vapor pressure deficits (Breshears et al. 2013), we found very few negative correlations between temperature and herbaceous ANPP (Data S1), indicating the general resistance of the herbaceous ANPP to a potential upward temperature shift at our study site. The positive correlations with February and October temperatures may extend the growing season, perhaps such that herbaceous ANPP might increase with higher temperatures in some instances. In the previous studies, negative impacts of summer season temperature on herbaceous ANPP were reported (Smoliak 1986, Patton et al. 2007, Ren et al. 2012. Woody understory ANPP had a greater number of negative v www.esajournals.org correlations with temperature (Data S1), and these appear related to drought as the stepwise regression models more often selected KBDI as opposed to temperatures (drought and higher temperatures were often correlated). We did not measure tree seedling or sapling survival so we cannot assess the potential impacts of hotter, drier conditions on regeneration or persistence.
During the course of our experiment, atmospheric CO 2 concentration increased from~350 to~405 ppm. This had the potential to benefit the woody and forb understory components (C 3 pathway) more so than the grass component (C 4 pathway; Tissue et al. 1995). Increasing atmospheric CO 2 favors the tree component in savannas and can facilitate expansion of woody plants Bond 2015, Scheiter et al. 2015). Our study focused on understory productivity and accounted for changing overstory conditions by incorporating basal area. However, elevated CO 2 can increase growth and water use efficiency of C 3 forbs and woody plants (Wayne Polley et al. 2002, Hatfield andDold 2019), which may increase understory ANPP or help ameliorate the negative effects of water stress. We did not examine the potential impact of CO 2 because the dominant understory vegetation across our site was C 4 grasses, the woody understory component was mainly composed of resprouts that were periodically top-killed and which were temporally and spatially variable, and the increase in CO 2 was confounded with time since treatment.

CONCLUSIONS
Forest-grassland transition zones are widespread globally, provide critical ecosystem services, and are sensitive to both climate variability and management. We found that the indirect effects of fire that stalled successional processes (maintain low tree cover) and removed litter were the most important factors influencing understory productivity. At our study site, we found that fire on a three-year or shorter return interval was needed to maintain savanna. From our results, tree basal area less than 10 m 2 /ha and litter accumulations less than 600 g/m 2 would be appropriate management targets to maintain productive understories and provide abundant forage and browse for livestock and wildlife.
The grass-dominated herbaceous ANPP of savanna systems was influenced by spring precipitation while understory ANPP in forest ecosystems by late summer soil water deficits. The effects of temperature were both positive and negative. Overall, in relatively moist southeastern Oklahoma, management decisions regarding the use of prescribed fire and overstory reduction will likely be more important than the direct effects of climate change. However, climate change and increased drought may lead to overstory mortality (Allen et al. 2010, Anderegg et al. 2013) and lower tree seedling survival and growth, which would have a positive feedback on understory productivity (Bond 2008).

ACKNOWLEDGMENTS
We thank the Oklahoma Department of Wildlife Conservation and Jack Waymire. We thank the numerous Oklahoma State University faculty, students, and staff and the staff of the Kiamichi Forestry Research Station for assistance with data collection and fieldwork. Comments from two anonymous reviewers were highly appreciated. The study was sponsored by USDA-National Institute of Food and Agriculture (Grant # 2017-06177). Additional funding was provided by Oklahoma Agricultural Experiment Station, McIntire-Stennis Project Number OKL0 2796, and Oklahoma Forestry Services.