Recently, significant die-back of nonnative common reed, Phragmites australis, has been reported in the Mississippi River Delta (MRD), Louisiana, USA. This dieback has been attributed to an invasive scale insect, Nipponaclerda biwakoensis. We test whether fungi are involved in the recent infestation by this insect and subsequent die-offs of Phragmites australis. Several haplotypes of P. australis occur in the MRD, and the European (M) and Delta (M1) haplotypes appear to experience differing levels of N. biwakoensis infestation. We tested whether these haplotypes differed in their fungal microbiomes in both their leaf and stem tissues, and whether differences in fungal community composition were linked to the level of infestation using a metabarcoding Internal Transcribed Spacer (ITS) amplicon sequencing approach. Our analyses showed differences in fungal community composition and diversity between haplotypes and tissue types, but none of these differences were directly correlated with N. biwakoensis infestation severity. However, we did find that the European haplotype hosted higher putative pathogen loads in stem tissues compared to the Delta haplotype, which may confer resistance to herbivory, though it is possible that differences in infestation between haplotypes are due to morphology.

Site and sample collection

On September 17, 2018, shortly after peak biomass, we traveled to the Passe-a-Loutre Wildlife Management Area of the Mississippi River, Louisiana, USA. Passe-a-Loutre is a coastal freshwater marsh with a mean salinity of less than 1 ppt. P. australis is the dominant vegetation in the region (Knight et al., 2020). 20 samples of each haplotype (European and Delta) were collected at two separate sites (Table 1) for a total of 80 individuals. The paired-haplotype plots at each site were identified by researchers at Louisiana State University on May 31, 2017 (Knight et al., 2020). Individuals of the European haplotype and Delta haplotype are easily differentiated in the MRD due to morphological and phenological differences, with the European haplotype exhibiting a much shorter, more gracile stature, earlier flowering time (all 40 European haplotype individuals collected were in flower at the time of sampling), and sparse or absent ligule hairs. The Delta haplotype is much taller, thicker-stemmed and possess dense ligule hairs (Hauber et al., 2011). No Delta haplotype individuals were in flower at the time of sampling. Culms of each haplotype were cut at the surface of the water if the plot was flooded (because scale insects do not occur under water) or at the soil in dry plots and placed individually in bags on wet ice for transport. Mean water depth at Site 1 was 19.8 cm and 2.6 cm at Site 2. Information on stand density and number of dead and living culms was determined for each individual using a circular quadrat which covers an area of ¼ m² placed with the selected culm in the center. On the day of collection, samples were placed in a 4°C refrigerator for holding and small subsamples of internodal (hereafter, “stem”) and leaf tissue from each of the 80 individual culms collected were placed in a -80°C freezer for metabarcoding analysis.

Using the remaining portions of the individuals, the number of N. biwakoensis insects per plant was determined by counting individuals every third internode and extrapolating to the entire culm. This was done by finding the average number of scales per node and multiplying by the number of internodes on the culm. Because the internode lengths differ within one individual, as well as between haplotypes, this value was standardized by the overall above-water height of the culm to determine scale density per centimeter of height. While the native scale insect Aclerda holci also occurs on P. australis in Louisiana, N. biwakoensis can be identified by its uniformly sclerotized and rounded abdomen (Knight et al. 2018) and was the only scale species present on the collected P. australis samples.

Surface sterilization

For metabarcoding, 10 culms from each plot were selected by arranging all collected culms within that plot in order of lowest scale density to highest. Starting with the lowest, every-other individual was chosen for analysis. This ensured we had a range of infestation represented for each plot and each haplotype. A 10cm section of healthy leaf tissue and a 4cm internodal section of stem from each sample were selected and placed in a tea strainer. Plant tissue was surface sterilized following protocol developed for Spartina patens (Lumibao et al., 2018). Tissues were sequentially soaked in 95% ethanol for 10 seconds, 0.5% sodium hypochlorite solution for 2 minutes, 70% ethanol for 2 minutes, and rinsed in sterilized deionized water for 2 minutes. Samples were then dried on UV sterilized KIMTECH Kimwipes before placing in gamma-sterilized cryovials for storage at -80°C.

DNA extraction

Following surface sterilization, samples were homogenized using liquid nitrogen and a mortar and pestle sterilized between samples. Genomic DNA extractions for the leaf samples were done following Qiagen DNeasy ^® PowerPlant^® Pro Kit protocol using 50mg of tissue, Phenolic Separation Solution and 250µl of Solution IR. Leaf samples were eluted using 100µl of Solution EB. Due to low yields for stem tissue, 75mg of stem tissue was used and elution was done using 50µl of Solution EB. Extraction concentrations were determined using a ThermoFisher Scientific Qubit fluorometer and all samples above a concentration of 10ng/µl were standardized to 10ng/µl by mixing extraction product with additional Solution EB.

PCR

To allow for maximum flexibility in choosing region-specific primers and dual-indexing barcode combinations, sequencing libraries were created in a two-step process following U’ren and Arnold (2017). PCR1, or amplification and primer ligation, was done in duplicate with an annealing temperature of 54.0° C. The primer sequences used were ITS1F 5’ – CTT GGT CAT TTA GAG GAA GTA – 3’ and ITS2R 5’ – GCT GCG TTC TTC ATC GAT GC – 3’ (Gardes & Bruns, 1993; White et al., 1990). Sample duplicates were then pooled and dual indexing primers were added so that no two samples contained the same combination of indexing barcodes. Five nanograms of DNA from each sample were pooled into a common library along with four negative controls which were included for identification of possible contaminant sequences. The pooled library was then purified and concentrated using Agencourt AMpure XP Beads with Dynamag-2 Magnet and sent to Duke University for Illumina MiSeq v3 sequencing (300 bp paired reads).

Sequence analysis

            Amplicon sequence variants (ASV) were identified following the DADA2 ITS pipeline (version 1.8) (Callahan et al., 2016) in R (R Core Team, 2019). ASVs are “phylotypes” that are a single DNA sequence (not clustered by a sequence similarity threshold like operational taxonomic units, OTUs). Primers and adapters were trimmed using the Biostrings package in R (Pagès et al., 2020) and data was denoised and paired reads joined using DADA2. No further sequence trimming was done, because DADA2 is robust to low-quality sequences through the incorporation of read quality information into its error model. Seven possible contaminant sequences were identified using the negative controls and removed using the decontam package in R (Davis et al., 2018). One sample was dropped due to very low reads. Taxonomy was assigned using UNITE 7.2 with a minimum bootstrap confidence of 70% (Nilsson et al., 2019) and functional guild assignments were obtained using FUNGuild (Nguyen et al., 2016).