Miscanthus sinensis Andress. is a seed-fertile crop introduced from Japan for horticultural use in the 1800s. It became naturalized in the eastern U.S. by the mid-20th century, and is locally invasive (Quinn et al., 2010). M. sinensis is both a candidate biomass crop as well as a parent species to other candidate varieties (e.g., Arnoult and Brancourt-Hulmel, 2015). We collected M. sinensis root-stock and seeds from roadside and forest opening patches in Daniel Boone National Forest, Powell County, KY, USA, in September 2009. Miscanthus × giganteus ‘Illinois’ clone (hereafter, M. giganteus) J.M. Greef & Deuter ex Hodkinson & Renvoize is a seed infertile hybrid of M. sinensis and M. sacchariflorus (Christian and Haase, 2001), and is one of the most widely planted cellulosic biofeedstock in the U.S. (Anderson et al., 2015). We obtained M. giganteus root stock for experimental plantings from the Chicago Botanic Garden. Utilizing plugs rather than seeds allowed us to evaluate plant survival and provide an estimate of long-term persistence once plants were introduced, while controlling accidental seed introductions into our study communities. We also estimated seed-based recruitment and overwintering persistence within our plots. Because the ‘Illinois’ clone is seed-infertile, we obtained seeds harvested from a pilot plantation of a pre-release, seed-fertile, tetraploid M. giganteus cultivar (‘PowerCane,’ Mendel Biotechnology, Hayward, CA, USA; see references in Bonin et al., 2017) to evaluate the potential for seed-based recruitment into study areas. Multiple studies have examined the seed and seedling viability and persistence of both ‘PowerCane’ and M. sinensis under various conditions and in various habitats (e.g., Smith and Barney, 2014; West et al., 2014a; Hager et al., 2015a; Smith et al., 2015; Bonin et al., 2017); these studies provide an understanding of the demographic contributions of seed to Miscanthus invasion potential. For our purposes, seed-based recruitment provides an additional measure of habitat suitability for escapes.

We established Miscanthus in three old field and three floodplain forest sites. Our old field habitats, Phillips and Trelease Prairies and the Vermillion River Observatory, and two floodplain forest sites, Nanney and Richter Tracts, are owned and managed by the University of Illinois Urbana–Champaign. The last floodplain forest habitat, Homer Lake, is part of the Champaign County, IL Park District. Consistent with management practices in our region, old field sites were mowed annually in the spring to a height of 7.5–10 cm to inhibit woody encroachment. Floodplain forests were unmanaged and subject to frequent and occasionally prolonged flooding.

Both Miscanthus species were propagated in the greenhouse prior to planting. We divided potted Miscanthus into approximately 10-cm diameter plugs with 10 to 15-cm long shoots. We hardened them off for a week and then transported them to the field. We introduced plugs into eight 10-m × 10-m single-species plots per site (four plots per species) in a split-split plot design (Supplement Figure 1A). Main single-species plots were divided into four 5-m × 5-m subplots that were each randomly assigned one of four density treatments: high (n = 16 plants with 1-m spacing, 1 plant/m²); medium (n = 9 plants with 1.25-m spacing, 0.56 plant/m²); low (n = 4 plants with 1.67-m spacing, 0.25 plant/m²); or control (n = 0 plants/m²). Plantings were positioned a minimum of 1-m from the subplot edge, and planting layouts were centered within the subplots (Supplement Figure 1B).

Introductions were initiated in April 2010. Plugs were planted into 25 cm deep holes and covered with a 25-cm × 25-cm heavy-gauge plastic mesh secured with sod staples to prevent movement due to flooding or animals. To maintain density treatments, plants that did not resprout by the spring census in 2011 and 2012 were removed and replaced. Nearly 75% of plugs (1037 out of 1392) were replanted in 2011 due to mortality; only 27% (379 out of 1392) were replanted in 2012. To control for effects of planting disturbance on community comparisons, we did sham plantings, which consisted of digging similar sized holes within the control subplots and then replacing the soil. The number of sham plantings was set as the median number of plantings (or replantings) out of all densities within that plot. We minimized soil disturbance during introduction, and any excess soil remaining after planting was transported out of the plot. Because Miscanthus establishment is sensitive to water limitation (Zub and Brancourt-Hulmel, 2010; Anderson et al., 2015), plugs were watered at the time of planting (or replanting), and periodically for the following month. Therefore, our data on establishment reflect a best-case scenario in terms of moisture conditions.

We measured a combination of soil fertility and soil water conditions to represent abiotic habitat differences. To quantify soil fertility, we collected soil samples to a depth of 10 cm from the center of each subplot with a 10 cm diameter soil corer in summer 2011. These samples were dried and sieved to remove non-soil particles, and analyzed by A&L Great Lakes Laboratories for plant macro- and micronutrients (P, K, Mg, Ca, S, Zn Mn, Fe, Cu, NO₃^-, NH₄⁺, and B), pH, soil organic matter, and cation exchange capacity (CEC). We measured soil redox potential as an integrated measure of saturated soil water conditions over the growing season using the Indicator of Reduction in Soil (IRIS) procedure (Castenson and Rabenhorst, 2006; Jenkinson and Franzmeier, 2006). We left IRIS tubes in the control subplot of each plot from April to August 2014, and recorded the amount of ferrihydrite paint lost from the tube surface at the end of the season. To quantify surface area exposed, we took digital images of each plane of the tube surface, and combined the multiple views into a single image. We then adjusted combined images to a consistent size (2745 × 675 pixels), and counted the number of pixels that lacked paint using MATLAB (Mathworks, Natick, MA, USA). We used these counts to calculate the percent paint surface area lost to redox reactions. Additionally, we did a pulse measurement of soil water content by sampling two 10-cm soil cores per sub-subplot 24–36 h after a rain of more than 2.5 cm in July 2013. We weighed wet samples, dried them for 72 h, took the dry weight, and subtracted the difference to estimate gravimetric soil moisture. These two measurements allowed us to compare relative differences in soil water status among plots.

To simplify the inclusion of soil nutrient conditions in the examination of habitat differences, and to account for strong covariance among the soil nutrients, we created a composite soil variable. We identified the optimal group of uncorrelated soil variables necessary to adequately distinguish habitats using a linear discriminant analysis (subselect R package: ldaHmat and eleaps functions). Improvement in correlation with a first canonical axis peaked with four soil factors (Mg, CEC, Fe, and B), and these variables successfully predicted habitat membership with less than a 5% error. Therefore, we combined them into a composite soil fertility variable quantified as the first axis scores from a principal components analysis of the four identified variables.

To represent relative differences associated with light availability within our two habitats, we evaluated light conditions both above and below understory vegetation (vegetation cover below any existing tree canopy). We measured photosynthetically active radiation above the understory canopy (approximately 1.5 m from the ground, in μmol m^-2 s^-1, PARA) and light transmittance (% PARA). These factors were quantified at the subplot level, and then averaged for analysis at the plot level. We measured both above and below understory canopy PAR with a linear ceptometer (LP-80 Accu-PAR, Decagon Devices, Court Pullman, WA, USA) as the average PAR at four points around each subplot. We quantified transmittance (the amount of above understory PAR penetrating to ground level) as below understory PAR/above understory PAR.

To test habitat type effects on seedling emergence, we established caged plots (to deter herbivory and seed predation) in one randomly selected corner of each plot in late fall 2011 (Supplement Figure 1B). Seeds were planted in 10 cm × 10 cm seed trays filled to a depth of 5 cm with soil from the receiving site. Each tray was placed in a cylindrical 1-mm mesh cage 40 cm in diameter and 30 cm high, which was additionally filled with site-collected soil to allow the seed tray to lie flush with the surrounding soil surface. The base of each cage was buried approximately 10 cm and secured in place. Because of site-use and material transfer agreement restrictions, we were unable to plant the seed-fertile M. giganteus in the field. However, previous work on the regeneration niche of this pre-release cultivar indicates fertile M. giganteus seed performs similarly to M. sinensis (Smith and Barney, 2014; West et al., 2014a). We monitored seedling emergence monthly March–November 2012, and seedlings were removed after each count to avoid confounding measures of emergence and survival.

To examine habitat type effects on seed viability after overwintering, additional seeds were cold-stratified in situ within stainless steel mesh packets buried in the field next to seedling plots in November 2012 (Supplement Figure 1B). Because it did not require field germination, we were able to use both species of Miscanthus for this test. We staked each 20 cm × 20 cm 0.5-mm wire mesh bag containing 100 seeds of the appropriate species in each corresponding subplot in November 2012. Bags were placed on bare ground, and any detritus moved to accommodate the bag was replaced to emulate site conditions. Bags remained in place until late April 2013, when they were collected and germinable seed fraction determined by counting the number of overwintered seeds that germinated under greenhouse conditions.

We recorded plant survival and tiller number twice annually from spring 2010 to 2014 to quantify plant performance and the integrity of density treatments over time. Census timing during the year was variable due to phenological fluctuation; thus, spring measurements occurred in April–May, and fall in October–November. Per year plant measurements, such as growth and survival, involved the period from spring to spring each year. In 2014, the second census was conducted in late July to optimize eradication efforts.

In 2014, we also quantified Miscanthus biomass within plots. We clipped and weighed all aboveground Miscanthus biomass per individual plant in the field. A subset of these plants were taken back to the lab, dried for 48 h at 45 degrees Celsius, and weighed to determine the relationship between field and dry weights (see footnote to Table 1). Additionally, we measured the area covered by each plant by measuring the widest axis of tiller extent and the one perpendicular to it, and calculating the area as an ellipse. Any flowering panicles produced in 2013 were collected, and the number of caryopses produced by habitat and by plant were quantified by weight relative to a 100 caryopsis weight standard. Viable seed production by M. giganteus ‘Illinois’ clone is inhibited by incomplete gametophyte development that results in sterility (Słomka et al., 2012). Therefore, these estimates merely represent a quantification of potential seed production within these habitats given the possibility of fertile genotypes being introduced for agronomic purposes (i.e., Bonin et al., 2017). Because our management agreements required panicle collection before dehiscence (to avoid unintentional spread and naturalization at the study sites), we were unable to reliably quantify viable seed produced for M. sinensis, and can only present relative differences in reproductive effort.

Plant community data were obtained by randomly sampling a total of 2 m² within each density subplot in late June-early July 2011–2014. Each year, we quantified the total number and percent cover of species in eight 25 cm × 25 cm quadrats along four randomly placed transects within each density subplot, and combined these 32 small-scale estimates of plant cover to represent community metrics at the density subplot level. Species richness (S) was the total number of unique species recorded within the density subplot. We calculated Pielou’s species evenness (J), where P_i is the relative contribution of ith species to total cover, and S is species richness [J = -Σ (P_i ^∗ ln(P_i)/ln (S)], i.e., Maron et al., 2014). Abundance observations were combined for the density subplot by converting cover estimates to area approximations (e.g., 80% of a 25 cm × 25 cm quadrat = 5 cm²) and summing them for each species.

We used non-metric multidimensional scaling (NMDS) to display the range of site variation in environmental variables (function metaMDS, vegan package, with Gower metric dissimilarities to account for variables measured at different scales), and its relation to plant survival per single species plot (overlaid using ordisurf function, see below for calculation). We evaluated the correlation strength of five environmental parameters (PARA, gravimetric soil moisture, soil redox potential, transmittance, and composite soil fertility) with the NMDS axis scores (envfit function), and plotted the significant habitat vectors with an r² greater than 0.5. Because these data were analyzed at the plot level, density was not included as a factor.

We quantified potential recruitment differences between habitats as: (1) overwintering seed viability (M. sinensis and M. giganteus); and (2) field emergence (M. sinensis). We analyzed mean recruitment differences between habitats separately for each species using linear mixed effects models (function lme in the nlme package) with plot nested within site as a random effect.

We quantified overall survival as the proportion of total individuals planted per plot that were alive in 2014. This total number of individuals planted per plot includes the initial 29 plants established in 2010, as well as any replants introduced in 2011 or 2012. Because these data were analyzed at the plot level, density was not included as a factor. We evaluated the influence of plant density on the resultant area occupied by Miscanthus in subplots by the end of the experiment (2014) using a linear model with pairwise Tukey tests to evaluate treatment differences.

To evaluate individual survival and growth, we analyzed differences in: (1) age-related survival probability [i.e., g(x) in Gotelli, 2008] and growth; and (2) plant performance between habitats. Because density had no significant effects on growth and survival parameters (see Results), it was excluded as a factor from these analyses.

Because plant survival was extremely low in floodplain forest habitats (see Results below), we limited our analysis of community impacts to the old field sites. We used a linear mixed effects model (lme in R nlme package) to examine relative impact of Miscanthus presence and density on community evenness and richness after 5 years (that is, in 2014). Miscanthus density × Miscanthus species were entered as fixed effects, and plot was nested within site as a random effect in the model.

We calculated the relative impact of Miscanthus on evenness and richness (two components of diversity) within each subplot using the relative impact (RI) equation from Vilà et al. (2006) (adapted from Armas et al., 2004).

A negative RI value indicates an increase of the dependent variable (V, e.g., richness or evenness) associated with invader presence (positive impact of invasion). Conversely, a positive value means that invader presence decreases V (negative impact of invasion). A zero value indicates that the invader presence does not have a significant effect on the parameter (Vilà et al., 2006). The RI for each Miscanthus introduction subplot was computed using values from the control subplot within the same plot. To further determine whether Miscanthus presence significantly affected overall evenness or richness of the plant community, we performed a single sample t-test of the relative impact on each of these metrics within Miscanthus subplots with a null hypothesis of μ = 0 (i.e., RI = 0). In other words, if the mean RI value is significantly different than 0, there is a significant impact of the invader, relative to the control, for the variable measured. Further, if the RI is significant, positive values indicate greater declines in the measured variable (negative impacts) associated with invader presence. Negative RI values indicate potential increases in the variable associated with invader presence.

We examined the influence of Miscanthus introduction and density on community dynamics to evaluate potential implications for old field community composition and structure over time. We determined annual species turnover (shifts in community constituents) and change in species rank abundances (shifts in species hierarchies) from 2011 to 2014. Two species (Festuca arundinacea and Solidago canadensis) present in all old field plots formed the dominant matrix vegetation (>50% cover in most old field plots). We excluded these two from community dynamics analyses to maximize our ability to detect initial short term impacts on species composition.

Species annual turnover was quantified as both the proportion of species dropping out of plot censuses between 1 year and the next (‘- species’) and species added between 1 year and the next (‘+ species’) (turnover in R codyn package). Mean change in species rank was quantified as the average difference in species rank abundance between consecutive years among species that were present across the entire measurement period, and represents community shifts in relative abundance over time (Hallett et al., 2016, rank_shift in R package codyn). To account for the possibility of changes in functional group abundance over time with Miscanthus introduction, we examined the effect of Miscanthus density on the relative prevalence of four functional groups (grasses, forbs, legumes, shrub/woody species) over time (adonis in R vegan package, strata = site). We also evaluated each variable for overall impacts (RI) from Miscanthus presence (single sample t-test with null hypothesis of RI = 0).

PARA, soil moisture, and soil fertility were strongly associated with the NMDS axis scores [r² = 0.92 (PARA); 0.66 (soil moisture); 0.84 (soil fertility), p < 0.01 for all]. These variables also provide good separation between the two habitat types (Figure 1). Although redox potential and transmittance also significantly separated between NMDS axes, their association was weaker [r² = 0.22 (redox); 0.32 (transmittance), p < 0.01 for both]. After the third year, we did not detect any Miscanthus in one floodplain site (Richter), and only two individuals of M. giganteus in another (Nanney). Several plants of both species did persist in the third floodplain site (Homer Lake). In contrast, all old field plots had many surviving individuals. Any small scale differences in how habitat characters within sites might have affected survival were obscured by the overwhelming survival difference between habitats (see below). Therefore, habitat characteristics associated with plant survival could not be statistically compared.

We tested the effect of overwintering on seed germination for both Miscanthus species. Seed germination after overwintering differed between Miscanthus species; however, there was no significant effect of habitat type (t = -2.82; df = 4, 40; p = 0.007, Table 1A). Although, we could not test M. giganteus for habitat differences in seedling emergence, M. sinensis did not differ in seedling emergence between habitats (t = 2.12; df = 4, 42; p = 0.1, Table 1A).

Miscanthus sinensis had both higher survival in old fields, and lower survival in floodplains, compared to M. giganteus (p_interaction < 0.01, df.residual = 42, z = -3.5, Table 1). M. giganteus survival was reduced nearly 85%, and M. sinensis survival over 90%, in floodplain forests compared to old fields (Table 1B). Overall, survival differed substantially between habitats for both species (Figure 1; df.residual = 20, M. giganteus: p < 0.01, z = -3.1; M. sinensis: p < 0.01, z = -3.8).

Increased planting density did increase the area of the plot occupied by Miscanthus species by 2014 (Figure 2). There was a significant difference in Miscanthus area with density treatment (p < 0.001, F = 5.44 on 5 and 64 df); however, the low and medium treatments were not significantly different from each other in final Miscanthus area (p > 0.01 for high versus low and medium; p = 0.18 for low versus medium) by 2014. Additionally, there was no difference between Miscanthus species (p = 0.28). There was also no significant effect of plant population density on either growth or survival (e.g., Figure 3; p > 0.1 for interactions of density with habitat, survival probability, and number of tillers). Density was therefore excluded from subsequent growth and survival analyses, and we present growth and survival at the whole plot level (out of 29 individuals, plus replants, planted in each single species plot).

Miscanthus presence was associated with marginal increases in richness relative to control plots (p < 0.001, t = -3.90, 71 df, Table 2A). Excluding Miscanthus additions, introduction plots had 1.8 (±0.4 SE) species more than control plots. Relative richness impacts (RI) were unrelated to Miscanthus species identity (p = 0.46, t = -0.74, 46 df) or average density within the subplots (p = 0.57, t = -0.56, 46 df) compared to controls. Similarly, Miscanthus presence did not affect evenness (p = 0.98, t = -0.03, 71 df), regardless of Miscanthus species (p = 0.83, t = 0.21, 46 df, Table 2B) or density (p = 0.06, t = -1.9, 46 df). Overall, we did not detect strong impacts of Miscanthus introduction on community structure metrics.

Species turnover did vary by year (- species: p < 0.001, t = -3.8; + species: p = 0.002, t = 3.1; 189 df). Additionally, turnover was unaffected by Miscanthus species identity (- species: p = 0.9, t = -0.2; + species: p = 0.2, t = -0.1; 92 df) or density (- species: p = 0.4, t = 0.8; + species: p = 0.2, t = -1.2; 189 df). Miscanthus introduction in general did not appear to exclude species (RI - species: p = 0.12, t = 1.55, 215 df), or increase new species occurrences (RI + species: p = 0.4, t = 0.8, 215 df; Figure 4) relative to control plots.

Miscanthus introduction did increase species rank shifts relative to control plots overall (RI rank shift: p = 0.03, t = -2.1, 215 df). However, mean species rank shift did not differ between years, species or with Miscanthus densities (linear model: Year: p = 0.6, t = 0.5, 189 df; Species: p = 0.2, t = 1.4, 92 df; Density: p = 0.3, t = 1.0, 189 df). The abundance of functional groups within plots varied by year and between Miscanthus species (Year: p = 0.006, F = 2.5, 1 df; Species: p = 0.004, F = 4.3, 1 df; Figure 5), but did not vary with Miscanthus density (p = 0.4, F = 4.3, 1 df). Overall, Miscanthus had some association with relative abundance of species and functional groups in the community. In addition to Miscanthus, introduction subplots had around two species more relative to control subplots, and RI varied by year and species. Miscanthus did not appear to exclude species from the community.