Historical Differentiation and Recent Hybridization in Natural Populations of the Nematode-Trapping Fungus Arthrobotrys oligospora in China

Maintaining the effects of nematode-trapping fungi (NTF) agents in order to control plant-parasitic nematodes (PPNs) in different ecological environments has been a major challenge in biological control applications. To achieve such an objective, it is important to understand how populations of the biocontrol agent NTF are geographically and ecologically structured. A previous study reported evidence for ecological adaptation in the model NTF species Arthrobotrys oligospora. However, their large-scale geographic structure, patterns of gene flow, their potential phenotypic diversification, and host specialization remain largely unknown. In this study, we developed a new panel of 20 polymorphic short tandem repeat (STR) markers and analyzed 239 isolates of A. oligospora from 19 geographic populations in China. In addition, DNA sequences at six nuclear gene loci and strain mating types (MAT) were obtained for these strains. Our analyses suggest historical divergence within the A. oligospora population in China. The genetically differentiated populations also showed phenotypic differences that may be related to their ecological adaptations. Interestingly, our analyses identified evidence for recent dispersion and hybridization among the historically subdivided geographic populations in nature. Together, our results indicate a changing population structure of A. oligospora in China and that care must be taken in selecting the appropriate strains as biocontrol agents that can effectively reproduce in agriculture soil while maintaining their nematode-trapping ability.


Introduction
Plant-parasitic nematodes (PPNs), especially the root-knot nematodes in the genus Meloidogyne, are widespread pests that cause crop yield losses worth more than US $157 billion worldwide each year [1]. For decades, the control of Meloidogyne spp. has heavily relied on chemical nematicides. However, resistance to chemical nematicides has emerged, and the environmental and human health impacts of nematicide residues are becoming increasingly recognized [2]. Therefore, currently available chemical nematicides are being phased out, and an increasing number of biocontrol agents are being introduced to help control PPNs [3].
At present, about 700 fungal species are known to be capable of attacking living nematodes (juveniles, adults, and eggs). These fungi are taxonomically diverse but traditionally divided into five groups based on the mechanisms by which they attack nematodes: Microorganisms 2021, 9,1919 2 of 23 (i) nematode-trapping fungi (NTF), (ii) endoparasitic fungi, (iii) egg-and cyst-parasitic fungi, (iv) toxin-producing fungi, and (v) fungi with special nematode-attacking devices [4]. Among these groups, the broad adaptability and flexible lifestyles of NTFs make them ideal agents for controlling PPNs [5][6][7]. NTFs are abundantly distributed in a broad range of habitats, especially in temperate agricultural pastures, coniferous leaf litters, and coastal vegetations [8]. Among the NTFs, A. oligospora has generally been considered to be the most common nematode predator in temperate soils [9]. In soil environments, A. oligospora naturally encounters a broad range of nematodes and behaves as a generalist predator with the characteristic ability of forming adhesive trapping nets when its mycelia are in contact with nematodes. Unlike endoparasitic fungi, NTFs are able to grow saprophytically. However, how they maintain and balance their saprophytic and predatory lifestyles are not known. In addition, how geographic populations interact with each other and how ecological factors impact population dynamics remain poorly understood. Such knowledge is important not only for understanding the diversity and evolution of these fungi but also for identifying the most appropriate isolate(s) for commercial biological-control applications [3,6].
Most ecological studies of nematophagous fungi have been restricted to surveys on their geographical and seasonal distributions, including examining the effects of abiotic (e.g., soil conditions) and biotic (mainly nematode density) factors on their distributions [8,[10][11][12][13][14][15]. Though molecular techniques have increasingly been used to examine phylogenetic relationships among nematode-trapping Orbiliales [16][17][18][19][20] and to investigate the molecular genetics of fungus-nematode interactions [21][22][23][24][25][26][27][28][29], they have scarcely been used to study the patterns of genetic variation present in NTF populations, including the spatial and temporal distributions of genotypes in specific species. The authors of an earlier study [30] analyzed 97 A. oligospora strains from several geographic locations and ecological niches in China using DNA sequences at three gene fragments: its (internal transcribed spacer region of the ribosomal RNA), tub (β-tubulin), and rpb2 (RNA polymerase II subunit). That study identified a large number of unique alleles and genotypes, as well as their limited geographic distributions, consistent with a large effective population size of A. oligospora and its potential for genetic differentiation among geographic populations, likely driven by ecological adaptation [30]. Further analyses of more samples using restriction fragment length polymorphism (RFLP) genotyping revealed a similar pattern [31]. Recently, the heritability of A. oligospora phenotypic variation in the response to nematodes was assessed by analyzing the genetic variation of genomic SNPs using 18 wild isolates [21]. However, due to limited sampling at most geographic locations and/or the small number of molecular loci used, A. oligospora's overall patterns of genetic variations, the relationships between genetic variations and phenotypic variations, and its mode of reproduction remain to be fully described.
For this study, we conducted a broad range of sampling on the model NTF species A. oligospora in China. These samples were analyzed with a panel of 20 short tandem repeat (STR or microsatellite) markers newly developed for this study. Combined with the DNA sequences of six housekeeping genes (its, tub, tef-1, rpb2, mitogen-activated protein kinase (mapk), and subtilisin-like serine protease (sp)) and mating type (MAT) genes (mat 1-1 and mat 1-2), the population structure, phylogenetic relationships, and mode of reproduction of the nematode-trapping hyphomycete A. oligospora in nature were investigated. In addition, we analyzed the patterns of dispersal and potential ecological niche adaptations of strains in China. Furthermore, representative strains from each of the robust phylogenetic lineages [30] were used to explore the intraspecific differentiation by studying the size and shape of the conidia, the growth rate of the colony under the experimental conditions, and their nematicidal activity. Our analyses revealed evidence for historical differentiation, as well as recent dispersion and hybridization. The implications of our results regarding the application of A. oligospora as a biocontrol agent against PPNs in agricultural fields are discussed.

Strains
This study analyzed a total of 239 strains of A. oligospora, including 104 previously analyzed strains [30,31] and 135 strains newly isolated from seven provinces. These isolates were obtained from agricultural and forest soils from 19 geographic populations in 14 provinces, covering seven large regions in China; the detailed geographic information of sampling sites and population size are given in Figure 1 and Table 1. The isolation and morphological identification of A. oligospora were based on the standard method described by Zhang et al. [31]. Briefly, for each sample, 1-2 g of soil/water were sprinkled onto sterile corn meal agar (CMA; 20 g of cornmeal, 18 g of agar, and enough water to reach a final volume of 1000 mL) and then free-living Panagrellus redivius nematodes were added onto the culture plate and inoculated at room temperature. Each sample was screened three independent times, resulting in 1500 agar plates for a total of 500 soil/water samples. After one month, nematodes were examined using a dissecting microscope. Single spores were directly isolated from the mycelia that trapped nematodes and then cultivated on a CMA Petri dish at 25 • C. The taxonomic characters were examined from cultures on CMA at 7-14 days after inoculation. The initial identification of A. oligospora was based on species-specific microscopic and macroscopic colony morphological features. The sizes of conidia and conidiophores were observed using an Olympus BX51 microscope (Olympus Corporation, Tokyo, Japan).
results regarding the application of A. oligospora as a biocontrol agent against PPNs in agricultural fields are discussed.

Strains
This study analyzed a total of 239 strains of A. oligospora, including 104 previously analyzed strains [30,31] and 135 strains newly isolated from seven provinces. These iso lates were obtained from agricultural and forest soils from 19 geographic populations in 14 provinces, covering seven large regions in China; the detailed geographic information of sampling sites and population size are given in Figure 1 and Table 1. The isolation and morphological identification of A. oligospora were based on the standard method de scribed by Zhang et al. [31]. Briefly, for each sample, 1-2 g of soil/water were sprinkled onto sterile corn meal agar (CMA; 20 g of cornmeal, 18 g of agar, and enough water to reach a final volume of 1000 mL) and then free-living Panagrellus redivius nematodes wer added onto the culture plate and inoculated at room temperature. Each sample wa screened three independent times, resulting in 1500 agar plates for a total of 500 soil/wate samples. After one month, nematodes were examined using a dissecting microscope. Sin gle spores were directly isolated from the mycelia that trapped nematodes and then cul tivated on a CMA Petri dish at 25 °C. The taxonomic characters were examined from cul tures on CMA at 7-14 days after inoculation. The initial identification of A. oligospora wa based on species-specific microscopic and macroscopic colony morphological features The sizes of conidia and conidiophores were observed using an Olympus BX51 micro scope (Olympus Corporation, Tokyo, Japan).   Table 1.

Identification of Tandem Repeat Loci and Development of Amplification Primers
STRs (or microsatellites) have been extensively used for molecular ecology and population genetic studies due to their high-level polymorphisms and high reproducibility in eukaryotic genomes [32][33][34][35][36]. In this study, the STR markers were developed based on the genome sequence of A. oligospora strain ATCC24927 (GenBank accession ADOT00000000), using the MISA software (Available online: https://webblast.ipk-gatersleben.de/misa/ (accessed on 12 April 2021)). From the candidate STR markers in the genome, we selected three di-, tri-, and tetranucleotide repeats based on loci with the highest repeat numbers and an interruption (max difference: 2 SSRs) of 100 bp within the boundaries of potential PCR primer regions [37]. Primer pairs were designed from flanking sequences of each STR locus (Available online: http://pgrc.ipk-gatersleben.de/misa/primer3.html (accessed on 13 April 2021)), and primers with lengths of 20-23 bp, a T m value about 60 • C, and product sizes in the 150-300 bp range were selected for further screening. For selected primer pairs, fluorescent tags were added to the forward primers ( Table 2). The fluorescent PCR products were examined via capillary electrophoresis assay on an Applied Biosystems (ABI 3730xl) sequencer (Applied Biosystems, Foster City, CA, USA). The GeneMapper 4.1 software (Available online: http://www.appliedbiosys-tems.com/Genemapper (accessed on 15 June 2021)) was employed to determine the polymorphism information content (PIC) for each STR marker.
PCR reactions were performed using 1 µL of genomic DNA template in a 50 µL PCR mixture (1.5 units of Taq DNA polymerase, 10 × PCR buffer (100 mM KCl; 200 mM Tris-HCl, pH 8.8; 100 mM (NH4) 2 SO 4 ; 20 mM MgSO 4 ; 1% Triton X-100; 1 and mg/mL BSA), 1.5 mM MgCl 2 , 2.5 mM dNTPs, and 1.5 mM primers) under the following amplifying conditions: 1 min of initial denaturation at 94 • C, followed by 30 cycles of 1 min of denaturation at 94 • C, 1 min of primer annealing at 50 • C, 1.5 min of extension at 72 • C, and a final extension period of 10 min at 72 • C. The PCR products were visualized under UV light exposure, and the purified fragments were directly sequenced on both strands with the same primers that were used for amplification.

Data Analysis
The genotyping of A. oligospora isolates was performed with both selected STR loci and single nucleotide polymorphic (SNP) loci from multilocus sequence typing (MLST). Alleles at the each STR locus were combined with SNPs to generate the multilocus genotype for each strain. All genotypes were entered into Microsoft Excel to conduct population genetic analyses by GenAlEx version 6.1 (The Australian National University, Canberra, Australia) [42]. For the LOCPRIOR model with admixture and correlated allele frequencies, the program STRUCTURE version 2.3.3 (Stanford University, Stanford, CA, USA) [43] was used to explore the number of genetic clusters (K) occurring in the sample. A total of 10 replicates were performed of each simulation for K = 1-12, with a burnin of 10,000 and an MCMC of 100,000 iterations for the best fixed value of K. Furthermore, a minimum spanning tree was constructed with default settings based on STR genotypes (BIONUMERICS 8.0, Applied Maths, Sint-Martens-Latem, Belgium) to investigate the genetic relationship between isolates from China. It has been reported that the corresponding teleomorph (sexual stage) of A. oligospora is Orbilia auricolor [44]. Interestingly, three other morphologically distinct anamorphs have been linked to the same teleomorphic species [45]. Theoretically, there are 15 different possible links between sexual and asexual morphs [46]. Thus far, the sexual stage has been rarely reported in laboratory cultures or nature for A. oligospora. Thus, its importance in nature remains largely unknown. In a previous study, a population in a stressful environment was shown to have significant evidence of recombination, likely due to sexual reproduction [31]. Here, to examine the reproductive mode in natural populations of this species, the distributions of mating type idiomorphs (mat1-1 and mat1-2) were analyzed from our A. oligospora strains and populations, and the following two complementary tests were conducted. First, we calculated the index of association (I A ) and rBarD. The null hypothesis for I A is that there is random association (recombination, linkage equilibrium) among alleles at different loci [47], and I A was standardized by the number of loci, with the rBarD algorithm adjusting for the numbers of loci. In the second test, we calculated the proportion of pairwise loci that were phylogenetically compatible (PrC). PrC is an indicator of recombination at the population level. In contrast, a lack of phylogenetic incompatibility implies asexual reproduction. Both tests were conducted using MultiLocus version 1.3 (Department of Biology, Imperial College, Ascot, UK) [47].

Comparison of Mycelial Growth and Conidial Shape
In our previous study, three divergent lineages (cryptic species) were suggested for samples from China, and these were also included in this study [30]. The comparisons of our intraspecific phenotypic variation were all based on divergence of the three lineages. To compare the growth rate and conidial yield of samples from the three different lineages (A-C) under different nutritional and temperature conditions, 7 mm diameter hyphal disks punched from the edges of plate colonies were inoculated onto three agar media (PDA; CMA; and tryptone, yeast extract and glucose agar (TYGA)), with each treatment incubated at three temperatures (25, 28, and 30 • C, respectively) for 3-7 days. The mycelial growth rate and colony morphology were examined and quantified at specific time intervals [48]. The aforementioned hyphal disks were also incubated on corn meal yeast extract agar (CMY) plates for 14 days at 28 • C. The conidial yield of each strain was then measured as previously described [49]. The morphological descriptions and photographs of hyphal and conidial morphology were observed and measured using fresh specimens growing on CMA. Both the lengths and widths of each conidium were observed and measured with an Olympus BX51 microscope with differential interference contrast. Measurement data were based on 50 random conidia of each strain, and the length-width ratios were used for comparative analyses among samples.

Trap Formation and Bioassay
The Caenorhabditis elegans employed in this study were cultured on nematode growth media (NGM) with Escherichia coli (OP50) as food source and maintained following standard procedures [50].
To investigate their nematode-trapping abilities, freshly harvested conidia (2 × 10 4 ) of strains of A. oligospora from the three different lineages were evenly spread on water agar (WA) and incubated at 28 • C for 3-4 days. About 200 nematodes were added to each WA plate to induce trap formation. Traps and captured nematodes were examined and quantified (per plate) under a microscope at 12 h intervals according to our previous observation of trap formation: the appearance of immature traps and mature ones could be found at every 12 h after nematode induction. Nematodes crawling on the wall of the plate were excluded because they were not in contact with the fungal hyphae. The experiments were performed in triplicate.

Nematicidal Activity of Fermentation Broth
To determine the potential nematicidal activity of fermented broth of each isolate, 1 cm diameter hyphal disks punched from the edges of plate colonies were inoculated into potato-dextrose broth medium (PDB) for 15 d at 28 • C and 180 rpm. The fermentation broth was then centrifuged at 11,000 rpm for 30 min. The supernatants were collected, and their nematicidal effects were measured in 24-well microtiter plates. In each well, 1 mL of a supernatant was mixed with 5 µL of live nematodes (at 50-100 L2 nematodes per µL). At room temperature, the number of dead nematodes and the total amount were counted at 2, 4, 8, 12, 24, 48, and 72 h, respectively. A nematode was considered dead if it was stiff or bent irregularly for a long time. A PDB medium and an M9 buffer were used as controls. The mortality of nematodes was defined as the ratio of dead nematodes over the tested nematodes, and the corrected mortality rate was calculated by the ratio of tested mortality minus those of controls over "1-mortality of controls".

Statistical Analysis
To compare the differences among strains and populations, we calculated the means and standard deviations for each treatment. A one-way analysis of variance followed by Tukey's multiple comparison test (when necessary) were used to analyze data, and p values < 0.05 were considered statistically significant. All statistical analyses were conducted using IBM SPSS statistical software V22.0 (IBM, Armonk, NY, USA).

Genetic Diversity of the Chinese Samples of A. oligospora Detected by STRs and MLST
Based on 200 initially designed STR primer pairs from the genome of A. oligospora (ATCC 24927), we selected 20 to genotype strains of A. oligospora. Information about these 20 STR markers in our sample is shown in Table 2. A previous study classified the polymorphism level of STR markers based on PIC values into low (PIC < 0.25), medium (0.5 > PIC > 0. 25), and high (PIC > 0.5) [51]. Based on this grouping, the 20 STR markers we developed for A. oligospora were found to have a moderate to high level of polymorphism, with 13 having high polymorphism, six having medium polymorphism, and only one having a low PIC value ( Table 2). The highest discriminatory power value for a single locus was obtained with A191, with a gene diversity value of 0.87 and 21 different alleles in our sample (Tables 1 and 2).
Among the 239 A. oligospora isolates from 19 different geographical populations, a total of 188 alleles and 178 multilocus genotypes were found based on the 20 STR loci. The number of alleles range from 3 to 21, with an average of 9.4 alleles per locus. Of the 188 alleles, 140 were shared between at least two of the 19 geographical populations in China, the remain 48 alleles were found only in one geographical population each ( Table 1). The 19 geographical populations also differed in their total numbers of alleles, which ranged from 21 (Kanas Lake, Xinjiang) to 88 (Heijing, Yunnan). Except for their absence in three geographic populations (Inner Mongolia, Qinghai, and Kanas Lake, and Xinjiang populations), private alleles were found in all other 16 populations (Table 1). Of the total 178 STR multilocus genotypes, only four were shared by two or more geographical populations, and the other 174 were only found in one geographical population each. The most shared genotype was found among two populations from Xinjiang and three populations from Dianchi, Yunnan, Sichuan, and Hubei. To better visualize the relationships among strains from the 19 geographic populations, a minimum spanning tree (MST) was constructed based on their STR genotypes. Figure 2a shows information on the geographic origins of the strains superimposed on the STR genotype relationships. Interestingly, strains from many geographic locations were often intermixed with each other ( Figure S1). However, samples from Hainan and Southwestern China (Yunnan, Sichuan, and Tibet provinces) formed several tight clusters on the minimum spanning tree, consistent with both localized genetic differentiation and long-distance dispersal among geographic populations in China. Notes: Each circle corresponds to a unique genotype, and the size of the circle proportionally represents the number of isolates with that genotype. Connecting lines correspond to the number of differences between the genotypes. Short bold line, 1 difference; black line, 2 differences; long grey line, 3 differences; dotted line, 4 or more differences. Notes: Each circle corresponds to a unique genotype, and the size of the circle proportionally represents the number of isolates with that genotype. Connecting lines correspond to the number of differences between the genotypes. Short bold line, 1 difference; black line, 2 differences; long grey line, 3 differences; dotted line, 4 or more differences.
The length of the combined sequences of the six gene was 2106 bp, with individual lengths ranging from 184 to 465 bp. Among the 2106 aligned nucleotides, 125 were polymorphic among the 239 strains of A. oligospora, and the percentages of nucleotides being polymorphic at each locus ranged from 0.0243 (tef-1) to 0.0995 (mapk), with an average 0.0594 (125/2106) (Table S1). Here, we found no clear correlation between the length of the sequenced gene fragment and the number of polymorphic nucleotide sites (R = 0.720; p = 0.106) among these six sequenced fragments. Based on MLST, a total of 59 multilocus genotypes were found among the 239 isolates. Among these 59 multilocus genotypes, 47 were only found in one geographical population each, while the remaining 12 were shared among geographic populations. Of the 59 MLST genotypes, Genotype #23 was the most frequently found at 29.7% (71/239), and it was found in 13 of the 19 geographic populations. The second and third most frequent MLST genotypes were Genotype #7 (9.2%; 22/239) and Genotype #1 (7.1%; 17/239), which were found in five and seven geographic populations, respectively. Except for populations in Zhejiang, Jilin, Kanas Lake and Urumqi Botanical Garden from Xinjiang, Guangxi, Gejiu in Yunnan, and Sichuan, private genotypes were found in 12 geographic populations (Table 3). Population from Yimen in Yunnan had the most private genotypes (n = 10), followed by Hainan and Heijing (Yunnan) populations, both containing seven private genotypes. Among the 19 geographic populations, the diversity of multilocus genotypes ranged from 0 to 0.981. Overall, the Heijing, Yunnan population had the highest genotypic diversity (0.981), followed by Qinghai (0.956).

Genetic Differentiation and Population Structure
To assess the relationships among the 19 A. oligospora geographic populations from China, we analyzed the genetic differentiation of these isolates by using the STR and MLST types of molecular markers (SNPs at all loci). After clonal correction, the sample sizes in the STR and MLST datasets were 186 and 100, respectively. The populations of Jilin and Xinjiang had only one genotype in the MLST dataset and thus could not be used for further population genetic analysis based on DNA sequences. However, based on the combined STR alleles and SNPs at all six MLST loci, 188 genotypes were found, and the combined genetic information was used for clonal correction for the following population genetic analyses.
Overall, in the MLST dataset, 19% of the total genetic variation was attributed to geographic separation among populations, and 81% was found within populations (PhiPT = 0.192; p = 0.006). The pairwise comparisons of the same dataset showed that among the 120 local population pairs, 29 pairs showed statistically significant differentiation (p < 0.05). The biggest differentiation was found between the Tibet and Hubei populations (PhiPT = 0.736; p = 0.004), followed by that between the Tibet and Guangxi populations (PhiPT = 0.647; p = 0.016) ( Table S2). The Mantel test further showed that there was significant correlation between geographical distance and Nei genetic distance (p = 0.01) ( Figure S2a). However, those from the STR dataset were quite different from the above-mentioned results revealed by MLST. AMOVA showed that more genetic variation (31% vs. 19% for the MLST dataset) was distributed among populations, and the remaining 69% genetic variations were found within populations in the STR dataset (PhiPT = 0.308; p = 0.001), indicating a higher level of genetic differentiation among populations by the new STR markers that that by sequence variation at the house-keeping gene based MLST. This pattern was similarly observed in the pairwise comparisons of genetic variance between populations, where the STR dataset showed that 95 out of 171 local population pairs were significantly differentiated (p < 0.05), which was more than those indicated by the MLST dataset. Based on STR markers, the Sichuan population was significantly different from all other geographic populations (p < 0.05) (Table S3). Moreover, the specific local population pairs with high genetic differentiation values were very different between the two datasets. The highest differentiation based on STR markers was between Inner Mongolia and Kanas Lake, Xinjiang populations (PhiPT = 0.994; p = 0.101, not significant), followed by that between Hubei and Sichuan populations (PhiPT = 0.869; p = 0.001) (Table S3). However, the correlation between geographical distance and Nei genetic distance indicated by the STR dataset was not significant (p = 0.24) (Figure S2b), which was not consistent with that of the MLST dataset.
AMOVA based on the combined dataset revealed that 72% of genetic variation was found within local populations, 28% of genetic differentiation was found among populations, and 132 out of 171 local population pairs showed statistically significant differentiation (p < 0.05); more of these pairs were found from populations from Southwest China (Sichuan, Guizhou, Tibet, and Yunnan including Dianchi, Heijing, Yimen, and Gejiu) and Hainan Island, as compared to those from other parts of China (Table S4). However, though there was a positive correlation between geographical distance and Nei genetic distance, such a correlation was statistically insignificant (p = 0.16) ( Figure S2c). Though significant genetic differentiations were found among many geographic populations, our analyses also showed evidence of gene flow among certain geographic populations. The existence of gene flow and genetic differentiation was further supported by results from the STRUCTURE analysis that showed that there were two genetic clusters (K = 2) in the Chinese population of A. oligospora, and most geographic populations contained alleles and genotypes of both genetic clusters. The population structure indicated by STRUCTURE software and PCoA was very similar, and both the MLST and STR datasets identified two distinct genetic clusters in all 239 A. oligospora isolates (Figure 3). Specifically, most samples from Southwestern China (Sichuan, Guizhou, Tibet, and Yunnan including Dianchi, Heijing, Yimen, and Gejiu) and Hainan Island formed one cluster, while those from other parts of China formed another cluster (Figure 3a-c). Populations in the former clade were found to have more frequent genetic exchanges than those in the second cluster (Figure 3d,e). However, the composition of two different genetic clusters indicated by the two kinds of molecular markers were quite different in some geographic populations. For example, STRs detected the existence of genetic elements of the second cluster in geographic populations from Hubei, Inner Mongolia, Jilin, Xinjiang, and Guangxi, while they were absent in the MLST dataset (Figure 3e). Interestingly, most genotypes containing elements of both genetic clusters (blue and red alleles in Figure 3d,e) were distributed in Southwestern China and Hainan. These genotypes likely represent hybrids of those two genetic clusters.

Phylogenetic Divergence
Phylogenetic analyses based on SNPs and MLST identified two distinct clades (I and II) within A. oligospora from 19 geographic populations in China (Figure 4, Figures S1 and S4; Table S5), consistent with the two genetic clusters revealed by the population structure analyses. Clades I and II of the MLST phylogeny included 38 and 21 multilocus sequence types, respectively. Most A. oligospora samples from different geographic populations (except for Guizhou) were widely distributed across Clade II. Interestingly, different geographical populations had different distribution patterns between the two clades: five geographic populations (Hubei, Jilin, Guangxi, and two populations from Xinjiang) were only found in Clade II, and samples from Tibet were only found in Clade I. Two populations (Hainan and Heijing from Yunnan) were widely distributed in both clades.

Phylogenetic Divergence
Phylogenetic analyses based on SNPs and MLST identified two distinct clades (Ⅰ and Ⅱ) within A. oligospora from 19 geographic populations in China (Figures 4, S1, and S4; Table S5), consistent with the two genetic clusters revealed by the population structure analyses. Clades Ⅰ and Ⅱ of the MLST phylogeny included 38 and 21 multilocus sequence types, respectively. Most A. oligospora samples from different geographic populations (except for Guizhou) were widely distributed across Clade ⅠI. Interestingly, different geographical populations had different distribution patterns between the two clades: five geographic populations (Hubei, Jilin, Guangxi, and two populations from Xinjiang) were only found in Clade ⅠI, and samples from Tibet were only found in Clade I. Two populations (Hainan and Heijing from Yunnan) were widely distributed in both clades. Significant differences in the relationships among strains among the analyzed gene fragments were identified based on the pairwise comparisons of topologies among all single gene phylogenetic trees ( Figure 5 and Figure S3). The results indicated evidence of potential allele-sharing and/or hybridization between the genetic lineages of this important NTF. For example, strains containing both blue and red alleles in the STRUCTURE analysis (colored in green on trees) formed tight and independent clusters on trees of mapk, rpb2 and sp fragments, suggesting they belonged to populations that are diverging from the other two clades. However, the same strains had mixed relationships with others based on sequences at its, tub, and tef-1 loci. Together, these differences in phylogenetic relationships among strains by different gene fragments are consistent with recombination in nature.  Table S5. Significant differences in the relationships among strains among the analyzed gene fragments were identified based on the pairwise comparisons of topologies among all single gene phylogenetic trees (Figures 5 and S3). The results indicated evidence of potential allele-sharing and/or hybridization between the genetic lineages of this important NTF. For example, strains containing both blue and red alleles in the STRUCTURE analysis (colored in green on trees) formed tight and independent clusters on trees of mapk, rpb2 and sp fragments, suggesting they belonged to populations that are diverging from the other two clades. However, the same strains had mixed relationships with others based on sequences at its, tub, and tef-1 loci. Together, these differences in phylogenetic relationships among strains by different gene fragments are consistent with recombination in nature.  Table S5.

Clonality and Recombination
A total of 197 A. oligospora strains were randomly selected, and their mating type (MAT) genes were successfully amplified. Among the 197 strains, 84 and 113 belonged to mat1-1 and mat1-2, respectively, in a ratio of~1:1.3, consistent with A. oligospora being a heterothallic fungus. Furthermore, the minimum spanning tree (MST) with the 197 strains based on the STR dataset showed that the A. oligospora strains containing mat1-1 and mat1-2 were both overall broadly distributed on the MST (Figure 2b). Indeed, some strains of different mating types shared the same STR genotype.
The rBarD and phylogenetic incompatibility tests for detecting recombination were conducted for (i) the total sample, (ii) samples from each geographical sample groups by both the STR and MLST datasets, and (iii) samples representing MLST genotypes in each of the phylogenetic clades (Table 4); the allelic sequences at each of the six MLST loci were treated as alleles in these tests. Though predominantly clonal population structures were detected for these sample sets, we found unambiguous evidences for recombination among the 20 STR loci (PrC = 0.01, p < 0.001; rBarD = 0.355, p < 0.001) and six MLST loci (PrC = 0, p = 1; rBarD = 0.6398, p < 0.001) in the total sample, including all 239 A. oligospora isolates. As indicated by the phylogenetic incompatibility test, variable levels of recombination were found within different sample groups. Interestingly, population groups from Southwestern China had the highest levels of phylogenetic incompatibility. It is worth noting that phylogenetic incompatibility was also identified in all samples where random recombination was rejected, indicating that low levels of recombination are common in samples of A. oligospora in China.

Conidial Morphology
For conidial morphology analyses, 38 strains of A. oligospora-including 30, 6, and 2 strains from the divergent Clades A, B, and C identified in our previous study [30], respectively-were selected to observe. Clades A and B respectively correspond to our Clades I and II here, and Clade C falls into the basal branch of Clade II. After conidia sporulation on CMA, the length and width of 50 mature conidia of each strain were observed and measured, and the length-width (L/W) ratio was selected as the quantitative measure for the following statistical analyses.
First, we tested variance homogeneity among the data from three different phylogenetic clusters. A Levene's p value of greater than 0.05 was obtained, suggesting there was no significant difference in the total variance among the three clusters. Second, an analysis of variance was conducted. The L/W ratio of Clade C strains was found to be significantly higher than those of Clades A and B by Duncan's multiple comparison test (Figure 6c). The L/W ratios of strains from Clades A and B were not significantly different from each other. Together, our results indicated some morphological differentiation among the clades within A. oligospora strains from China.

Mycelial Growth
A total of 12 strains representing three clades (six from Clade A, four from Clade B, and two from Clade C) from different geographic populations were selected to investigate their potential differences in mycelial growth rate and colony morphology on three different media at 25, 28, and 30 • C. The fungal colonies on the CMA medium were generally very loose and their aerial hyphae were sparse. By contrast, those on TYGA medium produced extremely dense mycelia and their aerial mycelia grew most robustly. The growth status of hyphae on the PDA medium was intermediate between CMA and TYGA (Figure 6b).
Overall, the optimal growth temperature of different A. oligospora strains was 25 • C. At 25 • C, all strains could grow normally and showed the fastest growth rate among the three tested temperatures. However, the growth rates of isolates from Clades B and C (except YMF1.02797 from Clade C) were slower than those of Clade A at 28 • C, and the aerial mycelia deviated from radial spread. In fact, isolates from Clades B and C were unable to grow at 30 • C (Figure 6a,b). Together, these results suggest clade-specific differences in colony morphology and mycelial growth, regardless of geographic origins. no significant difference in the total variance among the three clusters. Second, an analysis of variance was conducted. The L/W ratio of Clade C strains was found to be significantly higher than those of Clades A and B by Duncan's multiple comparison test (Figure 6c). The L/W ratios of strains from Clades A and B were not significantly different from each other. Together, our results indicated some morphological differentiation among the clades within A. oligospora strains from China.

Mycelial Growth
A total of 12 strains representing three clades (six from Clade A, four from Clade B, and two from Clade C) from different geographic populations were selected to investigate their potential differences in mycelial growth rate and colony morphology on three different media at 25, 28, and 30 °C . The fungal colonies on the CMA medium were generally very loose and their aerial hyphae were sparse. By contrast, those on TYGA medium produced extremely dense mycelia and their aerial mycelia grew most robustly. The growth status of hyphae on the PDA medium was intermediate between CMA and TYGA ( Figure  6b).
Overall, the optimal growth temperature of different A. oligospora strains was 25 °C . At 25 °C , all strains could grow normally and showed the fastest growth rate among the

Conidial Yield, Trap Formation and Nematode-Trapping Bioassay
A. oligospora, which can capture nematodes by forming three-dimensional networks from specialized hyphae, is considered the most abundant nematode predator in temperate soils [11]. In addition, genetic differences have been reported to have effects on the response of A. oligospora to its prey, though robust correlation has not been observed [21]. Therefore, the conidial yield, trap formation, and nematode-trapping bioassays were further investigated to examine whether there were notable differences among nematicidal activities among the three phylogenetic clades identified in a previous study.
Overall, we found significant differences in conidial yield among strains, but the differences among phylogenetic clades were not statistically significant. Strains from each of the phylogenetic clades have variable abilities to produce conidia. For example, three strains from the three clades (with one representing each clade)-YMF1.03101 (Clade A), YMF1.03135 (Clade B) and YMF1.03254 (Clade C)-were found to have significantly higher conidial yields than all other tested strains from either the same or different clades.
The number of traps produced by different strains of A. oligospora was most distinguishing when exposed to C. elegans for 12 h. The highest number was produced by strain YMF1.03063 of Clade B at >47% more than others, while the numbers produced by YMF1.03056 of Clade A and YMF1.02797 of Clade C were the least, only accounting for about 20% of that produced by YMF1.03063. The number of trappers produced by other tested strains were similar (Figure 6d). Overall, there was no significant correlation between the various abilities of the trap formation and phylogenetic division of A. oligospora strains from China. Correspondingly, limited difference regarding nematode-trapping abilities was found among the three phylogenetic clusters (Figure 6e).
We further detected the pathogenicity of fermentation broth from different representative strains in the three phylogenetic clades. Forty-eight hours after the nematodes were added into the fermentation supernatants, the lethal effects of different A. oligospora strains on nematodes was significantly different. Generally, strains YMF1.02856 (Clade A), YMF1.03095 (Clade B), and YMF1.03254 (Clade C) had obvious nematicidal effect, and the corrected mortality caused by strain YMF1.03095 was as high as 71.5%. However, strains YMF1.03056 and YMF1.02797, which showed the distinguishing ability to form traps and capture nematodes on the solid media, produced fermentation supernatants with weak pathogenicity (Figure 6f). Overall, we found no significant difference among clades in terms of the nematicidal activity of their fermentation broth.

Development of Novel STRs for NTF
For decades, molecular tools have acted as complementary techniques to differentiate species and examine phylogenetic relationships and systematics for NTF with traditional morphological methods [17][18][19][20]38]. Recently, MLST and RFLP have been used for ecological studies of nematophagous fungi [30,31], and Orbiliales-specific PCR primers were developed to directly detect NTF in environmental samples [52]. However, limited information is available on the patterns of genetic variations within and among geographic populations. In this study, we described a new panel of 20 STRs for the exact and highresolution genotyping strains of A. oligospora from China. About 2/3 of these STR markers were identified as highly polymorphic, and abundant allelic and genotypic diversities were detected among isolates of A. oligospora from China, thus indicating the high discriminatory power of our newly developed STRs.
Traditionally, the sprinkle plates method is used to isolate and identify NTFs based on their colony morphology and type of nematode traps produced. However, it is often difficult to analyze such data to infer the relationships among strains and populations. For example, the sprinkle plate method may be biased toward detecting species of nematophagous fungi that are adapted for growth in culture and have good predacious ability [53]. Here, we showed that the STR markers were highly discriminatory for identifying alleles and genotypes, and they enabled the generation of reproducible research results that could be easily shared among investigators when analyzing NTFs. Furthermore, the interactions among fungi, nematodes, and other soil organisms in different soils and the population dynamics of NTFs could be better uncovered with these markers.
In this study, the overall genetic clustering and phylogenetic separation indicated by the two kinds of molecular marks (STR vs. MLST) were similar. However, differences concerning genetic differentiation and the level of genetic exchange (hybridization) were also found. Specifically, more genetic variation among populations and frequent gene exchanges were revealed by STRs, which was likely due to the faster evolutionary rate of STRs than MLST. Indeed, the STR markers allowed for the detection of hybridization in most geographic populations, suggesting the hybridization events likely occurred very recently. Similarly, STRs failed to identify a statistically significant positive correlation between genetic and geographic distances or to discriminate the ancient population divergence that the MLST revealed. In addition, though the STR markers showed overall higher population differentiation than the MLST markers, they also showed greater gene flow among many of the 19 geographical populations. Together, due to the high mutation rate, the STR markers allowed us to capture the recent, fine-scale population genetic events that the MLST markers failed to reveal. Future studies should focus on analyzing the potential factors contributing to the differences in population structures, as revealed by MLST and STR markers, including factors such as ecological niche adaptation and host nematode distributions and preferences.

Historical Population Differentiation
Our study identified two distinctly differentiated genetic clusters among A. oligospora strains in China, indicating historical population divergence in this species. One of the major genetic cluster was composed of samples from Southwestern China and Hainan Island. Geographic populations from this region had a large number of private genotypes and unique alleles. The widespread genetic diversity and uniqueness in these areas suggested historical geographic isolation and, likely, the accumulation of locally adapted genotypes within and among local populations. Factors such as landscape features (e.g., mountain ranges), climate, and habitats in different geographic populations could all impact local ecology, genetic drift, selection, and adaptation. Specifically, Southwest China has a large number of mountain ranges and river systems that create a diversity of ecological niches, making it one of the most important biodiversity hot spots on earth [54]. Indeed, the high species diversity and endemism of fungal species and genotypes have been frequently revealed in this area [55]. For example, our recent population genetic investigations revealed multiple independent origins of the Aspergillus fumigatus pathogenic hyphomycetes in Yunnan, and they were significantly different from those in other parts of the world [56,57]. Furthermore, Hainan is the only island population with a tropical monsoon oceanic climate, and the geographic isolation between Hainan and the rest of China may result in the accumulation of unique genetic diversity. Together, the geographic and climatic factors likely contribute to the observed unique genetic diversities within individual local populations.

Recent Hybridization and Recombination
Our analyses revealed evidence for genetic exchanges among the strains and populations in Southwestern China and Hainan Island. As indicated by the population structure and gene genealogical analyses, there is extensive evidence for allele-sharing among lineages of A. oligospora, most likely due to hybridization. In addition to inter-lineage hybridization, we also found unambiguous evidence for recombination within individual lineages, though most of the samples showed evidence of non-random recombination. The frequent sharing of alleles but limited or no overlap in multilocus genotypes were further consistent with frequent recombination in nature.
The similar frequencies of mat1-1 and mat1-2 idiomorphs in our samples were also consistent with sexual reproduction and heterothallic lifecycle of A. oligospora in nature. Indeed, the widespread distributions of both mat1-1 and mat1-2 idiomorphs across geographic populations made mating between strains with different genetic elements possible. The finding that a couple of strains with different mating types shared the same STR genotype is also a strong piece of evidence for their recombinant origins. Our results here complement those in an earlier study that showed A. oligospora populations from stressful environments to be in greater linkage equilibrium than those from other environments in the same geographic areas [31].
In addition to extensive allelic and multilocus genotypic sharing (hybrids) in our samples, our other analyses such as gene genealogical comparison and STRUCTURE also suggested frequent hybridizations, with different strains showing variable proportions of genetic elements from the two genetic clusters. Similar population structures have been found in several other fungal populations in Southern China, such as those of the opportunistic human fungal pathogens Aspergillus fumigatus in Yunnan [56] and Candida tropicalis in Hainan Island [58], the plant pathogen Fusarium oxysporum in South-Central China [59], and the wild edible mushroom Russula virescens Ally from Southwestern China [60]. Our results were consistent with the idea of gene flow being a major force shaping the A. oligospora population from China. As indicated in our previous studies, the long-distance dispersal of asexual spores by natural forces, anthropogenic activities, the transport of agricultural and forestry products, and animal (soil nematode and insects) vectors all could have contributed to the spread of genetic elements of A. oligospora to different ecological and geographic niches [30]. With increasing anthropogenic influences such as human travel and the transport of good and services among regions, gene flow among geographic popu-lations of A. oligospora will likely increase, which could further facilitate hybridization and recombination, potentially accelerating the distribution of adaptive alleles originating in one area to broader geographic reaches [61].

Intraspecies Diversification
The Fungal kingdom is among the most diverse and specious groups of eukaryotes. At present, we are still far from knowing the full extent of fungal diversity. One major finding based on DNA sequences over the last two decades has been the existence of multiple genetically divergent cryptic species within many previously defined species [61][62][63][64][65]. Within NTF, morphology-based classifications of strains and species have also likely underestimated the true magnitude of species diversity in this group of fungi. Though the globally isolated 22 isolates of a promising biocontrol agent Duddingtonia flagrans have a single recent common ancestor with limited genetic differentiation among geographic populations [66], several species within the NTF genera Arthrobotrys and Monacrosporium have been found to contain significant genetic variation and cryptic species among morphologically similar isolates [67].
Over the last two decades, the phylogenetic species concept (PSC) has been frequently used for identifying fungal species [64]. According to the PSC, species are diagnosed as a cluster of individuals that are sufficiently differentiated from other clusters, as revealed by DNA sequences. Given the presence of genetically differentiated clusters of strains in the phylogenetic tree, we sought to further determine whether the observed genetic divergence within A. oligospora was associated with phenotypic characters such as morphology and nematode-trapping abilities. As indicated by shape of conidia (L/W), the morphological differentiation of A. oligospora strains was found to be associated with phylogenetic lineage differentiations. However, we also observed large variations among strains within each of the lineages, a result consistent with significant observed genetic variations within individual lineages of this species. The observed variations suggest in field applications as biocontrol agents against plant-parasitic nematodes, a mixture of genotypes from the locally sourced strains might be more effective than a single strain.
The genetic clusters and phylogenetic lineages identified here represent speciating populations within A. oligospora. At the same time, the identification of hybrids in most geographic populations suggests that strains in different clusters/lineages are sexually compatible. Indeed, the comparable distributions of mat1-1 and mat1-2 genes suggest that sexual mating is likely common in natural populations of this species. Analyses of laboratory crosses between strains from the same and different lineages could help reveal the underlying genetic mechanisms responsible for the observed phenotypic variations. Indeed, we are in the process of identifying strains of opposite mating types and conducting laboratory crosses for strains between lineages to examine their cross-fertility and the extent of potential reproductive isolation. Species and cryptic species have been defined based on many different criteria [61], including geographical origin [68], ecological adaptation [69,70], and host specialization [71]. Host preference has been found in pairs of tightly interacting hosts and pathogens due to shared common evolutionary histories. However, in facultative pathogens, the selective pressure exerted by both hosts and habitats could be equally important [72,73]. At present, the extent of nematode-host specificity and how such specificity may be related to phenotypic diversification and ecological niche adaptation remain largely unknown.
In selecting strains for field applications as biocontrol agents against PPNs, a common approach is to select strains with the strongest nematode-trapping ability, as demonstrated in laboratory settings. However, we believe many other factors, including the specificity of NTF-PPN interactions and the growth and reproductive abilities of the selected strains in field conditions, should be considered. The large number of strains and genotypes identified here for A. oligospora from China will allow us to test the importance of some of these factors for future biocontrol applications.

Supplementary Materials:
The following are available online at https://www.mdpi.com/article/10 .3390/microorganisms9091919/s1. Figure S1: Genetic clustering of 239 A. oligospora isolates from China based on 20 STR loci by Nei's distance. Notes: Scale bar indicates percentage identity. Figure S2: Mantel tests of the relationships among genetic differentiation (Nei values) and geographical distance. (a) MLST dataset (p = 0.01), (b) STR dataset (p = 0.24), and (c) combined dataset (p = 0.4). Figure S3: 14 Pairwise tanglegrams between single gene phylogenies. Notes: Hybrids indicated by structure analysis are colored in green on the tanglegram. Figure S4: Bayesian phylogeny for the MLST SNPs. Table S1: Length and polymorphism information of six sequenced DNA fragments.  Figure 4.