Mitogenomics, Phylogeny and Morphology Reveal Ophiocordyceps pingbianensis Sp. Nov., an Entomopathogenic Fungus from China

The new entomopathogenic fungus Ophiocordyceps pingbianensis, collected from Southeast China, was described by mitogenomic, morphological, and phylogenetic evidence. The systematic position of O. pingbianensis was determined by phylogenetic analyses based on six nuclear gene (ITS, tef1-α, nrSSU, nrLSU, rpb1 and rpb2) and 14 mitochondrial protein-coding gene (PCGs) (cox1, cox2, cox3, atp6, atp8, atp9, cob, nad1, nad2, nad3, nad4, nad5, nad6 and nad4L) data. Phylogenetic analyses reveal that O. pingbianensis was belonged to the Hirsutella nodulosa clade in the genus Ophiocordyceps of Ophiocordycipiaceae. This fungus exhibits distinctive characteristics which differed from other related Ophiocordyceps species with slender and geminate stromata, monophialidic conidiogenous cells with an inflated awl-shaped base, a twisty and warty phialide neck and a fusiform or oval conidia, as well as being found on a tiger beetle of Coleoptera buried in moss at the cave. The complete mitochondrial genome of O. pingbianensis was a circular DNA molecule 80,359 bp in length, containing 15 PCGs, 24 open reading frames genes (ORFs), 25 transfer RNA genes (tRNAs) and 27 introns. Ophiocordyceps pingbianensis, containing 27 introns, has the second largest mitogenome in Ophiocordycipiaceae and was next to O. sinensis. To our knowledge, this is the first report of the mitogenome from a new entomopathogenic fungus, and thus provides an important foundation for future studies on taxonomy, genetics and evolutionary biology of Ophiocordycipiaceae.

Asexual morphs of Ophiocordyceps consist of Hirsutella Pat., Hymenostilbe Petch, Paraisaria Samson and B.L. Brady, Sorosporella Sorokin, Stilbella Lindau, as well as Syngliocladium Life 2021, 11, 686 2 of 16 Petch [1,14]. The main anamorph of Ophiocordyceps refers to Hirsutella. The fungal genus Hirsutella was a pathogen of insects, mites and nematodes, which was erected by Patouillard in 1892, and Hirsutella was initially described by Petch in 1924. In 1998, Hodge introduced the monograph and molecular phylogeny of Hirsutella. To be specific, she reviewed the taxonomy of Hirsutella species and described 68 species to pertain to Hirsutella [16]. Hirsutella was considered to show relationships to the genus Ophiocordyceps typified by a sexual morph [1]. The morphological characteristics of this genus was considered to produce basally inflated phialides, tapering towards the apex, discontinuous hymenial layer and conidia embedded in a mucous sheath [15,17,18].
The mitochondrion acts as a master regulator of metabolism. As DNA sequencing technology leaps forward, the mitochondrial DNA (mtDNA) is recognized as an effective marker for phylogenetic analysis, with the advantage of exhibiting a faster mutation rate than nuclear DNA [19]. Moreover, mtDNA is also appropriate for phylogenetic analysis due to its high copy number and conservative gene functions [20].
Furthermore, there are differences in the size of fungal mitochondria in different species. From the current reports, the smallest mitogenome known in the fungi was 12.1 kb in Rozella allomycetis (Doweld) Letcher in Cryptomycota [21], while the largest was 272,238 bp in Morchella importuna M. Kuo, O'Donnell and T.J. Volk in Pezizomycetes [22]. The mitochondrial genome of the Ophiocordycipitaceae species is like an obvious closed circular structure, comprising naked double stranded DNA. It commonly contains 15 protein-coding genes (e.g., 3 subunits of cytochrome c oxidase (cox1-3), 3 subunits of ATP synthase (atp6, 8 and 9), cytochrome b gene (cob), 7 subunits of NADH dehydrogenase (nad1-6, nad4L) and 1 ribosomal protein S3 (rps3), as well as ribosomal RNA genes (rnl, rns) and transfer RNA genes (trn). Thus far, mitogenomes of only ten species in Ophiocordycipitaceae have been described, in which mitogenomes of Hirsutella species have been rarely reported. The existing sequence data of mitogenome from NCBI only include H. minnesotensis Sen Y. Chen, Xing Z. Liu [23], and the smallest was 52,245 bp in H. minnesotensis [24]. Different variations of mitogenome sizes of hypocrealean fungi were largely attributed to differences in the number/length of introns and the length of intergenic regions, and the mentioned introns most likely were obtained through horizontal transfer from other fungal species [25].
In the present study, an unknown species of Ophiocordyeps attacking a larva of tiger beetle was collected from Daweishan National Nature Reserve, Pingbian County, Yunnan Province, China. Mitogenomics, phylogeny and morphology of this fungus were determined, and its systematic position was established in Ophiocordycipitaceae, It was revealed a new species of Ophiocordyceps with a Hirsutella morph. The complete mitogenome of the new fungus was sequenced, assembled and annotated. Genome structures, gene contents, codon usage and gene arrangement were analyzed.

Fungal Materials and Isolation
In this study, live specimens were collected on a tiger beetle of Coleoptera buried in the moss in a cave in Daweishan National Nature Reserve, Pingbian County, Yunnan Province, China. Before examination, the specimen was transferred to a laboratory and then stored at 4 • C. To detect and further study the culture, stromata were photographed and then measured with an Olympus SZ61 stereomicroscope. A stroma of the fungus growing from the beetle larva fell to small segments nearly 5 mm long. The respective segment was surface-sterilized in 30% H 2 O 2 for 5 min, then soaked in 70% alcohol for 2 min and subsequently rinsed sufficiently in sterilized water, and lastly dried on the sterilized filter paper. Afterwards, the segments were placed on Potato Dextrose Agar (PDA) plates. The specimen was deposited in Yunnan Herbal Herbarium (YHH), Yunnan University. The cultures of this fungus were deposited in Yunnan Fungal Culture Collection (YFCC), Yunnan University.

Morphological Observations
The fresh specimen, including stromata and host, were photographed with Olympus CX40 and BX53 microscopes. The colonies on PDA plates were cultured at 4 • C for 4 weeks, and the colonies characteristics (size, texture and color) were photographed with the Cannon 700D camera to characterize the morphology of colonies. For to asexual morphological descriptions, microscope slide cultures were prepared by placing a small amount of mycelia on 5-mm diameter PDA medium blocks that were overlaid by a cover slip. The cultures on the slants were transferred to PDA plates and then cultured with an incubator for 25 d at 4 • C. Next, the colonies were photographed and then measured every four days. The micro-morphological observations and measurements were performed under Olympus CX40 and BX53 microscopes, as well as a FEI QUANTA200 scanning electron microscope [26].

DNA Extraction, PCR and Sequencing of Nuclear Genes
Axenic living cultures of the new species were collected from a PDA plate to prepare for DNA extraction. The DNA extraction was performed with the CTAB method in Liu et al. [27]. The nuclear ribosomal small subunit (nrSSU) was amplified with the primer pair nrSSU-CoF and nrSSU-CoR [28]. The nuclear ribosomal large subunit (nrLSU) and translation elongation factor 1α (tef-1α) were amplified with the primers LR5 and LR0R [29,30] and EF1α-EF and EF1α-ER [1,30], respectively. To amplify the largest subunits of RNA polymerase II(rpb1) (rpb2), the primer pair RPB1-5 F and RPB1-5 R, as well as the primer pair RPB2-5 F and RPB2-5 R, were applied respectively [31]. The nuclear ribosomal internal transcribed spacer region (ITS) was amplified using the primer pair ITS4 and ITS5 [32]. The polymerase chain reaction (PCR) assay was performed as described by Wang et al. [33]. PCR products were sequenced with an automatic sequence analyzer (BGI). When PCR products could not be sequenced directly, cloning was performed with the TaKaRa PMD™18-T vector system (TaKaRa Biotechnology, Dalian, China).

Sequencing, Assembly and Annotation of Mitogenome
After the DNA extraction and the quality detection, the extracted genomic DNA was transported on dry ice for high-throughput sequencing at Novogene Co., Ltd. (Beijing, China) by applying dry ice for high-throughput sequencing. Following standard procedures, the extracted genomic DNA was sequenced on the Illumina HiSeq 3000 platform using a 350 bp paired-end reads, which resulted in 3 GB of data and 100× sequencing depth. Using the high-throughput sequencing data, mitogenome of the new species was assembled with the software SPAdes v. 3.11.0 [34], which involved three processes, i.e., (1) initial splicing, (2) contig screening and extension and (3) hole filling and re-splicing. The complete mitogenome of the new species was annotated by combining the results of a MFannot tool [35] and an ARWEN web server, with artificial correction. The content of the annotated sample mainly included the protein coding genes (PCGs), the rRNA genes, the tRNA genes and the open reading frames (ORFs). Next, the graphical map of the complete mitogenome of the new species was generated with the software Organellar Genome DRAW tool [36].

Phylogenetic Analyses
Raw trace files were inspected and edited using MEGA7. Two species of Tolypocladium inflatum W. Gams and T. ophioglossoides J.F. Gmel were designated as the outgroup to root the trees. Five gene (nrSSU, nrLSU, tef-1α, rpb1 and rpb2) sequences of Hirsutella, Ophiocordyceps and the outgroup, were retrieved from GenBank. The majority of ITS sequences of Hirsutella were also retrieved from GenBank. Some manual adjustments were made as the nrSSU gene contained an intron (Table S1). The sequences were aligned with the programmes Clustal X2.0 and MEGA7 [37,38]. Phylogenetic analyses were performed by combining the concatenated 14 PCGs sequences generated here with other mitogenomic sequences of species in Hypocreales downloaded from NCBI. Neurospora crassa Shear and B.O. Dodge and Podospora anserine Rabenh were designated as the outgroup. The 14 mitochondrial PCGs were aligned using MUSCLE. [39]. Maximum likelihood (ML) and Bayesian Inference (BI) phylogenetic trees were estimated using RaxML 7.0.3 and MrBayes v.3.1.2. [26,40,41].

Phylogenetic Analyses Based on Nuclear Genes
For the ML and BI trees estimated for the five nuclear genes, the new species was closely related to H. liboensis X. Zou

Phylogenetic Analysis Based on Mitochondrial Genes
The ML tree and BI tree were estimated for the mitochondrial PCG dataset of 52 species in Ascomycota from GenBank (Table S2). The six well-supported clades of Hypocreales were recognized by ML and BI analyses, including species of the families Bionectriaceae, Clavicipitaceae, Cordycipitaceae, Hypocreaceae, Nectriaceae and Ophiocordycipita-

Phylogenetic Analysis Based on Mitochondrial Genes
The ML tree and BI tree were estimated for the mitochondrial PCG dataset of 52 species in Ascomycota from GenBank (Table S2). The six well-supported clades of Hypocreales were recognized by ML and BI analyses, including species of the families Bionectriaceae, Clavicipitaceae, Cordycipitaceae, Hypocreaceae, Nectriaceae and Ophiocordycipitaceae ( Figure 3). As indicated from the phylogenetic analysis, the new species was a member of the family Ophiocordycipitaceae. The new species was clustered together with H. minnesotensis and H. vermicola in the genus Ophiocordyceps of Ophiocordycipitaceae.

Mitogenomic Characteristics of Ophiocordyceps Pingbianensis
The complete mitogenome of O. pingbianensis was an obvious circular molecule with a length of 80,359 bp ( Figure 5). It had an overall GC content of 29.89%. The mitogenome of O. pingbianensis was a rather high capacity and compacted genome with genic regions (including intronic regions) accounting for 71.54% and intergenic regions taking up 27%. The protein-coding regions took up 92.19% of the whole mitogenome of O. pingbianensis, containing 2 rRNA genes (i.e., rnl and rns) genes, 15 standard protein-coding genes (PCGs), 24 ORFs, 25 transfer genes (tRNAs) ( Table S3). The sizes of rns and rnl were 1652 bp and 6444 bp, respectively. The total length of the 15 PCGs was 14682 bp, including 3 subunits of cytochrome c oxidase (cox1, cox2 and cox3), 3 subunits of ATP synthase (atp6, atp8 and atp9), one cytochrome b gene (cob) and seven subunits of NADH dehydrogenase complex (nad1, nad2, nad3, nad4, nad5, nad6 and nad4L) and 1 ribosomal protein S3 (rps3). The lengths of 25 tRNAs ranged from 69 to 84 bp and 24 ORFs (orf146, orf132, orf293, orf276, Sexual stage: The fruiting body was not observed. Notes: Ophiocordyceps pingbianensis is characterized by geminate stromata, slender, unbranched and a surface color from light deep yellow to dark brown; monophialidic conidiogenous cells with an inflated awl-shaped base, twisty and warty phialide neck, fusiform or oval conidia. Morphological comparison showed that O. pingbianensis was distinct from other closely related species (Table 1), which was concordant with molecular phylogenetic results.

Transfer RNAs
In the mitogenome of O. pingbianensis, the 25 tRNAs coded for all 20 amino acids (Table S6). As indicated from the results, all tRNAs adopt a typical cloverleaf structure except for trnS-GCT, trnY-GTA, trnL-TAA, trnL-TAG and trnS-TGA ( Figure 6). As shown in Table S6, some of the 25 tRNAs existed in multiple copies, the trnM-CAT gene for methionine appeared in three copies. Two different tRNAs for arginine (trnR-ACG and trnR-TCT), Serine (trnS-GCT and trnS-TGA) and Leucine (trnL-TAA and trnL-TAG) were found.

Discussion
By morphological and molecular phylogenetic researches, O. pingbianensis was found as a novel species, being belonged to the H. nodulosa clade. Morphological differences still exist between O. pingbianensis and its related species, such as the conidiogenous cells of H. liboensis with being smooth and polyphialidic [43] [47]. The mitogenome size of O. pingbianensis was the second largest because it contained 27 introns. The mitochondrial genome of H. minnesotensis exhibited significant intron degeneration for the nad4L and cox1-i1 genes, which might be caused by unexpected stop codons or frame shifting [24]. In H. rhossiliensis, there were 10 group I introns and one unclassified intron in six genes (i.e., rnl, cob, cox1, cox3, nad1 and nad5) [48]. This indicated an intron presence/absence dynamics cause mitogenome size variations in the Ophiocordyceps species. On the whole, 27 ORFs were identified in O. pingbianensis, and except for orf144, orf131, all ORFs were transcribed at the identical orientation. In the mitogenomes of five Ophiocordyceps species reported, the number of ORFs was also different, e.g., the H. thompsonii mitogenome containing 3 ORFs [25], H. rhossiliensis mitogenome containing 5 ORFs [48], H. minnesotensis mitogenome containing 4 ORFs [24], as well as the H. vermicola mitogenome containing 3 ORFs [47].

Discussion
By morphological and molecular phylogenetic researches, O. pingbianensis was found as a novel species, being belonged to the H. nodulosa clade. Morphological differences still exist between O. pingbianensis and its related species, such as the conidiogenous cells of H. liboensis with being smooth and polyphialidic [43]. O. pingbianensis has a warty phialide neck of conidiogenous cells, whereas a larger periclinal protuberance was produced near the apex of conidiogenous cells of O. unituberculata H. Yu, Y.B. Wang and Y.D. Dai [18]. Conidia of O. pingbianensis were smooth, fusiform or oval, inconsistent with O. unituberculata with larger lanceolate to fusiform conidia and H. nodulosa with non-mucoid conidia. The host of O. pingbianensis was a tiger beetle of Coleoptera; however, these hosts of its closely related species H. satumaensis, H. nodulosa, H. liboensis and O. cochlidiicola are Lepidoptera.
The mitogenome of O. pingbianensis was a circular DNA molecule with a length of 80,359 bp, and the mitogenome size was the second largest within the mitogenomes of the Ophiocordycipitaceae species reported. The mitogenome sizes were different among species in Ophiocordyceps attributed to differences in the number of introns. For instance, the H. thompsonii mitogenome exhibited 15 introns [25], the H. rhossiliensis mitogenome had 13 introns [48], the H. minnesotensis mitogenome had 13 introns [24] and the H. vermicola mitogenome exhibited 7 introns [47]. The mitogenome size of O. pingbianensis was the second largest because it contained 27 introns. The mitochondrial genome of H. minnesotensis exhibited significant intron degeneration for the nad4L and cox1-i1 genes, which might be caused by unexpected stop codons or frame shifting [24]. In H. rhossiliensis, there were 10 group I introns and one unclassified intron in six genes (i.e., rnl, cob, cox1, cox3, nad1 and nad5) [48]. This indicated an intron presence/absence dynamics cause mitogenome size variations in the Ophiocordyceps species. On the whole, 27 ORFs were identified in O. pingbianensis, and except for orf 144, orf 131, all ORFs were transcribed at the identical orientation. In the mitogenomes of five Ophiocordyceps species reported, the number of ORFs was also different, e.g., the H. thompsonii mitogenome containing 3 ORFs [25], H. rhossiliensis mitogenome containing 5 ORFs [48], H. minnesotensis mitogenome containing 4 ORFs [24], as well as the H. vermicola mitogenome containing 3 ORFs [47].
To gain more insights into the variability of the mitochondrial gene among Ophiocordyceps species, the genomic arrangements of O. pingbianensis with other five mitogenomes in Ophiocordyceps were analyzed comparatively. The lative positions of most PCGs and most tRNAs of the six Ophiocordyceps species were relatively conservative, except for several tRNAs and PCGs distributed in different positions, which resulted from the variability of the relative gene order of these mitogenomes.
As revealed from genomic synteny, the six Ophiocordycep mitogenomes falled to 4-5 homologous regions. Homologous region B only existed in the mitogenome of O. pingbianensis and H. minnesotensis, however, it was lacked in H. thompsonii, H. vermicola, H. rhossiliensis and O. sinensis mitogenomes. The relevant orders of four homologous regions (i.e., A, C, D and E)were conserved across the six Ophiocordyceps mitogenomes and homologous regions B were distributed in different positions among the mitogenomes; as a result, the diversity of gene order was generated among six mitogenomes. Gaining insights into gene arrangements, component and genomic synteny of Ophiocordyceps species might help assemble mitogenomes and trace the evolutionary history of other Ophiocordycipitaceae species in future research.

Supplementary Materials:
The following are available online at https://www.mdpi.com/article/ 10.3390/life11070686/s1, Table S1: Specimen information and GenBank accession numbers for sequences used in this study. Table S2: Specimen information and GenBank accession numbers for mitogenomes used in this study. Table S3: General features in the mitogenome of Ophiocordyceps pingbianensis. Table S4: Gene component of the mitogenome from Ophiocordyceps pingbianensis. Table S5: Codon usage of protein-coding genes in the mitogenome from Ophiocordyceps pingbianensis. Table S6: tRNAs in the mitogenome from Ophiocordyceps pingbianensis.