Study Site
We conducted our experiment in Long Point Provincial Park, Long Point, Ontario, CA (42o35`01” N, 80 o22` 37” W; Appendix A). Long Point is a freshwater coastal marsh located on the north shore of Lake Erie and is a designated World Biosphere Reserve, Ramsar Wetland, and a globally significant Important Bird Area. Wetland vegetation in these marshes assemble along a water depth gradient from deep standing water (up to 70 cm) to shallow standing water or saturated soils.
This range of water depths is subdivided into three distinctive vegetation communities: emergent marsh (deepest), sedge meadow (intermediate) and meadow marsh (shallowest). Each of these vegetation communities comprises multiple plant species, but evenness is typically low (Robichaud and Rooney 2021), with each community dominated by a characteristic rhizomatous perennial that reproduces sexually and asexually. Coincidentally, these are traits that are shared by P. australis. For emergent marsh, the dominant species is most likely Typha x glauca, but identification based on morphology is difficult in the field due to extensive hybridization (Travis et al. 2010; Bansal et al. 2019) so we instead refer to “Typha spp.” throughout. Sedge meadow was most commonly dominated by Carex aquatilis (water sedge [Wahlenb.]) and meadow marsh by Calamagrostis canadensis (Canada bluejoint grass). Hereafter, Typha spp., C. aquatilis, and C. canadensis are referred to as “resident species.” Like P. australis, Typha x glauca generates dense monocultures (e.g., Galatowitsch et al. 1999), and C. canadensis is taxonomically closely related to P. australis as they both belong to the Poaceae family.
Experimental set up
For our experiment we selected replicate phytometers (an individual stem/ramet used to measure plant responses to experimental manipulations) of P. australis and the three resident species representative of each remnant vegetation community. All phytometers were situated along the leading edge of a P. australis stand (Appendix A). This best approximates realistic competition in invaded areas, as P. australis relies mostly on clonal expansion once it has established itself through seeds or clonal propagules (Kettenring et al. 2016). To minimize intraspecific competition, we selected resident phytometers that were growing surrounded primarily by P. australis within the stand and P. australis phytometers that were growing surrounded by each resident species. All of the phytometers were established in the same area, with a maximum distance of 150 m between plots (Appendix A).
We selected 96 phytometers on 23-May-2016 (Appendix B) and 27-May-2017 (Appendix C), for a total of 192 phytometers. Each year, these included 24 phytometers of Carex aquatilis and 24 of C. canadensis. We included only 12 phytometers of Typha spp. each year as these were more robust and less prone to damage during the field season. We included 36 phytometers of P. australis (12 each to compare with the three resident species). Phytometer pairs of equivalent height were identified for each species and one member of each pair was assigned at random to the ‘with interspecific aboveground competition’ treatment while the other was assigned to the ‘without interspecific aboveground competition’ treatment. This ensured that both treatments included phytometers spanning the full range of early-growing season ramet heights (Appendix D).
To create the ‘without interspecific aboveground competition’ treatment, we clipped all plants growing in the 1 m2 area surrounding the phytometer to within 2 cm of the soil. For the ‘with interspecific aboveground competition’ treatments, we did not alter the surrounding above-ground biomass within 1 m2 around the phytometer. However, to limit potential clonal subsidy and standardize belowground interactions, as all the phytometer species are rhizomatous clonal species, we severed roots and rhizomes in both the ‘with competition’ and ‘without competition’ treatments by sawing the perimeter of the 1 m2 plots to a depth of 50 cm with a hand saw. We performed this once, as pilot work determined that severing the below-ground material repeatedly caused physical disturbance detrimental to the phytometers.
Approximately every 10 days, we measured the height of phytometers and re-clipped the surrounding vegetation in the ‘without competition’ treatment. Over the course of the experiment all phytometers were subject to natural herbivory and physical stresses. Phytometers that were consumed or died are reported in Appendix B & C.
Once phytometers had matured (in July, but dates varied between years with interannual differences in weather), we measured the carbon assimilation rate (µmol CO2 s-1 m-2) (A) and photosynthetic water use efficiency (CO2 mmol s-1 m-2 H20) (WUE) of each phytometer using a CIRAS-3 true differential gas analyzer with a PLC3 Universal LED Light Unit (RGBW) and PLC3 narrow leaf cuvette (PP Systems, Amesbury, MA, USA). We selected a fresh, entire (e.g., no damage) leaf growing with maximum sun exposure from each phytometer and then measured a photosynthesis-irradiance (PI) curve in the field. The PI curve began by exposing the leaf to 1500 µmol s-1 m-2 of photosynthetically active radiation, equivalent to an average full-sunlight day during the growing season, and slowly reduced PAR to 1000, 500, 200, 100, 50, and 0 µmol s-1 m-2 while simultaneously measuring carbon assimilation and photosynthetic water use efficiency. Measurements at each PAR level were taken until carbon assimilation rates plateaued, which typically occurred within two to three minutes. We took these measurements from July 26th to August 2nd in 2016 and from July 4th to July 14th in 2017, with phytometers of the same species measured on the same day to reduce potential temporal differences in performance between treatments. We also measured the amount of PAR reaching the top of each phytometer relative to the incident PAR above the canopy using a LI-1500 Light Sensor Logger coupled with two LI-190R quantum sensors (Li-Cor Biosciences, Lincoln, NE, USA). These sensors were deployed to take simultaneous readings from above the canopy and at the top of the phytometer to most accurately calculate the percent of incoming PAR intercepted by the canopy. PAR measurements were taken on cloudless days, between 09:00 and 15:00 h.
To compare a proportion of the aboveground biomass produced by each species over the growing season, we we clipped each phytometer at the base of the stem right before peak aboveground biomass occurred in early August of 2016 and mid-July of 2017. Research in these marshes determined that peak biomass for emergent marsh, meadow marsh, and P. australis vegetation communities occurred in mid-August 2016 and late July 2017 (Yuckin and Rooney 2019). Due to clonal origin of the phytometers and the extent of interweaving of roots and rhizomes within the upper 40 cm of sediment (Lei et al. 2019), below-ground biomass could not be accurately determined.
Resource Measurements
To characterize niche overlap among species, we collected environmental variables from sites dominated by each of the phytometer species: P. australis (n = 15), Typha spp. (n = 15), C. aquatilis (n = 15) and C. canadensis (n = 15) for a total of 60 sites. For resident species, we selected patches of remnant vegetation at equivalent water depths but not experiencing direct interactions with P. australis to best approximate their realized niche within the marsh. Sites for each species were a minimum of 10 meters from one another, and all of the sites were situated within 1000 m of the phytometers (Appendix A). At each site we collected a 10 cm deep soil core to measure soil nutrients, measured soil moisture using a WET sensor kit and HH2 moisture meter (Delta-T Devices, Burlington, ON), and the measured the percent of incident PAR intercepted by the canopy by deploying the Li-Cor sensors to take simultaneous readings from the top of the canopy and the soil or water surface on cloudless days between 09:00 and 15:00 h.
Laboratory Analyses
Each phytometer was dried at 80oC for 24 hours and then weighed on an analytical balance to the nearest 0.0001 g to determine the aboveground biomass. For a subset of 48 phytometers from the 2016 season: ten individuals of each of the resident species (5 ‘with competition,’ and 5 ‘without competition’) (10 x 3 = 30 samples), and six P. australis phytometers (3 ‘with competition,’ 3 ‘without competition’) for each of the neighbouring species (6 x 3 = 18 samples), we also measured the carbon (% dry weight) and nitrogen (% dry weight) content from selected leaves. To relate these nutrient content measures to photosynthetic performance (e.g., Hirose and Werger 1987; Hirtreiter and Potts 2012), we also measured the δ13C isotopic composition of each selected leaf. The plants in our study are all C3 photosynthesizers which have δ13C values that range between -20 to -37‰ (Kohn 2010). The C3 photosynthetic pathway discriminates against the heavier 13C isotope during stomatal diffusion and carboxylation by Rubisco (Fry 1992) - plants that discriminate less between C isotopes typically photosynthesize more efficiently and have a higher (less negative) δ13C value (Farquhar et al. 1989; McAlpine et al. 2008).
Selected leaves were ground into a homogenous powder and a subsample of 1 mg collected for analysis of C, N and δ13C by the Environmental Isotope Laboratory at the University of Waterloo. First, samples underwent combustion conversion to gas through a 1108 Elemental Analyzer (Fisons Instruments) coupled to a Delta Plus XL (Thermo-Finnigan, Germany) continuous flow isotope ratio mass spectrometer. The %N and %C element content is a bulk measurement based on the sample weight against known certified elemental standards. The δ13C values were corrected to the primary reference scale of Vienna Pee Dee Belemnite. Every fifth sample was duplicated for precision quality control/quality assurance. Three of the duplicate samples (two C. canadensis, one P. australis) were outside of the calibration range and were removed from the dataset, leaving seven duplicates for precision analysis. Analytical precision was measured using relative percent difference between duplicates (EPA 2014), and the average precision was 0.002% (± 0.113 st. error) for δ13C, 3.25% (± 3.12 st. error) for %C, 4.89% (± 3.83 st. error) for %N.
Soil samples collected for resource measurements were dried at 80oC, ground into a homogenous mixture, then analyzed for soil pH, phosphorus (mg/Kg), carbon (% dry weight), total nitrogen (% dry weight), calcium (mg/Kg), potassium (mg/Kg), magnesium (mg/Kg), sodium (mg/Kg), copper (mg/Kg), iron (mg/Kg), zinc (mg/Kg dry), manganese (mg/Kg), and sulfur (µg/g). Plant available phosphorus was measured using sodium bicarbonate-extractable phosphorus following Reid (1998). Total nitrogen (TN) and carbon were measured using thermal conductivity detection (Reid 1998). The K, Mg, Ca, and Na samples were extracted using 1.0 N Ammonium Acetate solution, following Simard (1993). Copper, iron, and zinc samples were extracted using a 0.005M DTPA solution and the filtrate was analyzed by ICP-OES following Liang and Karamanos (1993). Manganese was measured using 0.1 N phosphoric acid as the extracting solution, following Reid (1998). For sulfur, homogenized samples were closed-vessel microwave digested with nitric acid and hydrochloric acid, then the microwave digested sample was brought to volume with Nanopure water and quantification was performed using ICP-OES (AOAC 2011.14). Nitrogen and carbon analyses were done at the Biogeochemical Analytical Service Laboratory at the University of Alberta while the other nutrient analyses were conducted by the Agriculture and Food Laboratory at the University of Guelph.
Statistical Analyses
We determined that year did not influence plant biomass (general linear model F1,164 = 2.40, p = 0.123) or carbon assimilation rates at full sunlight (1500 µmol s-1 m-2) (general linear model F1, 164 = 0.004, p = 0.948), so we pooled the data from both years for these variables. We used two-way ANOVAs to compare carbon assimilation rates and water use efficiency at 1500 µmol s-1 m-2 among species and between the treatments. We ran four models in total with either carbon assimilation rate or water use efficiency as the response variable, with an interaction between phytometer species (for resident phytometers) or phytometer neighbours (for P. australis phytometers) and treatment (with or without competition). We used the same model design to assess differences in δ13C. Duplicates in the δ13C data were averaged to one value for analyses. We used Type III sums of squares, unless an interaction was not significant, then we report Type II sums of squares. If a fixed factor was significant, without a significant interaction term, we used Tukey’s HSD post-hoc test to assess differences among levels of the factor. Analyses were performed using the car package (Fox and Weisberg 2019) and agricolae (de Mendiburu 2020).
To evaluate the allocation of resources to above-ground biomass between treatments, we compared the yield of phytometers growing with competition to those growing without competition. Using the phytometers that were paired by height at the beginning of the experiment, we calculated differences in above-ground yield using the relative competition index (RCI) approach (Grace 1995; Goldberg et al. 1999; Vilà and Weiner 2004):
(Ywithout competition – Ywith competition)/Ywithout competition
Yield (Y) represents the above-ground biomass of each phytometer. As the weights are standardized, values greater than 0, with a maximum of 1, indicate that the above-ground biomass of the plant growing with competition was lower than its counterpart growing without competition. Values < 0 indicate the above-ground biomass of the plant growing with competition was higher than its counterpart growing without competition. This allows us to compare the differences in above-ground yield among the species while accounting for variation in size among species.
The soil data, except pH and light, collected from the unmanipulated sites were converted to ppm (i.e., mg/Kg) and log transformed to improve normality. To control for collinearity among environmental variables, we summarized the underlying correlation structure using principal components analysis (PCA). We created a matrix of soil nutrients, pH, soil moisture, and proportion of incident PAR reaching the ground and conducted the PCA, with a correlation matrix, using the `rda` function in vegan (Oksanen et al. 2020). The PCA scores were then multiplied by the proportion of variance explained by each axis, to give them appropriate weight, and were used as an indicator of ecological niche to quantify trophic niche region and overlap among the plant species using nicheROVER (Lysy et al. 2017). To estimate pairwise niche overlap, nicheROVER employs a Bayesian framework to calculate the probability that an individual from species A is found in the niche region (a 95% probability region in multivariate space) of species B (Swanson et al. 2015). All analyses were performed using R v. 4.0.3 (R Core Team 2020).