Evidence for the Role of CYP51A and Xenobiotic Detoxification in Differential Sensitivity to Azole Fungicides in Boxwood Blight Pathogens

Boxwood blight, a fungal disease of ornamental plants (Buxus spp.), is caused by two sister species, Calonectria pseudonaviculata (Cps) and C. henricotiae (Che). Compared to Cps, Che is documented to display reduced sensitivity to fungicides, including the azole class of antifungals, which block synthesis of a key fungal membrane component, ergosterol. A previous study reported an ergosterol biosynthesis gene in Cps, CYP51A, to be a pseudogene, and RNA-Seq data confirm that a functional CYP51A is expressed only in Che. The lack of additional ergosterol biosynthesis genes showing significant differential expression suggests that the functional CYP51A in Che could contribute to reduced azole sensitivity when compared to Cps. RNA-Seq and bioinformatic analyses found that following azole treatment, 55 genes in Cps, belonging to diverse pathways, displayed a significant decrease in expression. Putative xenobiotic detoxification genes overexpressed in tetraconazole-treated Che encoded predicted monooxygenase and oxidoreductase enzymes. In summary, expression of a functional CYP51A gene and overexpression of predicted xenobiotic detoxification genes appear likely to contribute to differential fungicide sensitivity in these two sister taxa.


Introduction
Boxwood blight is a disease of native and commercially grown boxwood (Buxus spp.) caused by two emerging fungal plant pathogens: Calonectria pseudonaviculata (Cps) and C. henricotiae (Che) [1,2]. During infection, plants undergo stem and foliar lesioning, rapid defoliation, and die back, and preventative fungicides are required for disease management. Cps infects Buxus plants globally, while Che is restricted to Europe [1][2][3][4]. Both species were previously classified as distinct genotypes of Cps, but multi-locus phylogenetic analysis, comparative genomics, and phenotypic differences supported their placement into separate sister species, with the erection of Che [3,5].
Various classes of fungicides and antifungal drugs are used to control plant and animal diseases caused by fungal pathogens, including boxwood blight [6][7][8]. One group of antifungals, the azole class (or demethylation inhibitors/DMIs), impede production of ergosterol, a vital component of fungal cell membranes [9,10]. This class of antifungal chemicals has broad application in both agriculture and medicine and is used to treat and prevent a broad range of plant and animal diseases. Azoles inhibit fungal growth by affecting membrane integrity and are targeted to the active site of sterol 14α-demethylase, a cytochrome P450 protein encoded by the CYP51 gene, which catalyzes demethylation of lanosterol in the ergosterol biosynthesis pathway. Differential sensitivity towards select azoles, such as the triazoles tetraconazole and propiconazole, has been demonstrated by in vitro comparisons between Cps and Che, where the latter species tolerated significantly higher exposure to the chemicals in a poison plate assay [3].
As with known antimicrobial compounds, repeated exposure to azoles places selective pressure on fungi, by which the CYP51 target-site may gain resistance-conferring mutations that inhibit azole binding [11][12][13]. In addition to mutations in the active site, overexpression of the CYP51 gene is also known to confer varying levels of resistance to azoles in diverse fungi [13,14]. Select lineages, including Cps and Che, are reported to possess three CYP51 paralogs, denoted as CYP51A, CYP51B, and CYP51C [9,13,15,16]. From a comparative analysis of a worldwide collection of Cps and Che genomes, a CYP51A pseudogene was found in all examined Cps isolates from populations spanning four continents, but not in Che [16]. A premature stop codon occurring upstream of the translated active site, presumably disabling CYP51A function, was hypothesized to play a role in differential azole sensitivity seen in Cps relative to Che [16]. However, ergosterol biosynthesis is a multistaged process involving over 20 genes, overexpression of which may serve as putative determinants of decreased azole sensitivity [17,18].
Transcriptomics studies have reported that when exposed to antifungals, fungi with reduced azole sensitivity display increased expression of genes from the xenome, a threephase detoxification pathway found in fungi and other eukaryote lineages [19][20][21][22]. During the first phase, proteins oxidize harmful substances entering the cell [20]. In this context, cytochrome P450 monooxygenase enzymes aid metabolism of xenobiotics, such as antifungal drugs and plant secondary metabolites [20,21]. Although specific enzymes in this family can contribute to routine cellular activities/functions, such as lipid processing and synthesis of secondary metabolites, their role to aiding xenobiotic detoxification can be vital to non-target site mechanisms of tolerance to antifungals [21,22].
The second phase of detoxification involves deployment of proteins that perform chemical conjugation (modification) of xenobiotic substances [20]. Through conjugation the oxidized drug/toxin can be rendered more soluble, facilitating removal from the cell. Alternatively, the compound could be primed for breakdown or compartmentalized within the fungal cell. Among the proteins involved in this stage of detoxification are transferases, largely glutathione transferases [20,23,24]. The latter behaves by conjugating glutathione onto toxin molecules, facilitating efflux and thus reducing toxicity [23,24].
The final phase, known as efflux, is the removal of toxins, minimizing the opportunity for toxins to bind and inhibit cellular functions [13,20]. The ATP-binding cassette (ABC) transporters, ATP-dependent pumps that transport substances across membranes, are known to export toxins, including azoles, from the cell, and are often overexpressed in azole-resistant fungi [22,25,26]. Additionally, members of the major facilitator superfamily (MFS) of active transport proteins efflux azoles in parallel with ion movement along an electrochemical gradient [27,28]. Evidence in Alternaria has demonstrated that the loss of functional MFS transporters was linked with increased sensitivity towards fungicides, demonstrating the importance of efflux transporters in antifungal tolerance [28].
In light of previous findings on differential fungicide sensitivity in these two closely related sister species causing boxwood blight disease, we set out in the present study to examine gene expression patterns that may underpin their differential azole fungicide sensitivity [3,16]. We analyzed RNA-Seq expression data of control and tetraconazole-treated Cps and Che to examine CYP51 paralog expression patterns, and to identify differentially overand underexpressed genes in Cps and Che likely associated with xenobiotic detoxification or ergosterol biosynthesis [12,20,21]. Such genetic determinants could serve as potential targets for future functional studies and aid in teasing apart mechanisms of antifungal resistance to inform development and screening of novel antifungal compounds.

Sensitivity of Cps and Che to Select Azole Fungicides
Che isolate CB45 had insignificantly higher EC 50 and EC 85 values when treated with either myclobutanil or tebuconazole (Table 1). Cps CT1 had an insignificantly higher EC 50 value and an insignificantly lower EC 85 value following treatment with triadimefon when compared to Che CB45 (0.86 times larger and~1.41 times lower, respectively). However, Che isolate CB45 demonstrated substantially larger EC 50 values when treated with tetraconazole (five times higher) or propiconazole (three times higher) compared with Cps CT1. CB45 also had a significantly larger EC 85 value following tetraconazole treatment (over 16 times larger).

Predicted Detoxification Genes in CT1 and CB45
CT1 and CB45 had 504 and 544 genes, respectively, possessing at least one of the four InterPro IDs associated with phases one and three of fungal xenobiotic detoxification (Tables S2 and S3) [20]. Using the broad GO name search, we identified 831 and 940 genes from CT1 and CB45, respectively, predicted to encode proteins with transferase activity (acetyl-CoA, glutathione, acetyl, etc.) (Tables S4 and S5).

Differentially Expressed Genes in Tetraconazole-Treated Cps
For all expression analyses, log 2 FC (fold change) values greater than 1.5 were considered overexpressed, as previously specified by [21], while log 2 FC values below −1.5 were considered underexpressed. Sixteen genes were statistically, significantly overexpressed in tetraconazole-treated CT1 compared to untreated CT1 (Table 2), and 55 genes whose expression significantly declined in CT1 following exposure to tetraconazole. Using our multi-blast analyses, we saw that 64 of these 71 total genes in Cps CT1 had a corresponding hit in Che CB45. None of the genes significantly overexpressed in CT1 were substantially underexpressed in CB45 (Table 2; Figure 1A). Rather, seven genes underexpressed in CT1 matched with overexpressed genes in CB45, though only one of which (a predicted 2,4-dienoyl-CoA reductase precursor, g12476) mapped to an overexpressed CB45 gene (log 2 FC > 1.5). Though generally not reaching our minimum threshold for a meaningful change in gene expression, genes in Cps CT1 following tetraconazole treatment were more likely to decline in expression.   We observed that only six of the genes had log 2 FC values greater than 1.5 (Table 2; Figure 2A). An amine oxidase gene was the most overexpressed among treated isolates (log 2 FC ≈ 2.90), though two of the genes were identified as encoding cytochrome P450 monooxygenases (one of which with a log 2 FC > 1.50), which can be putatively linked to chemical detoxification [20,21]. At least three genes were predicted to produce dehydrogenases (none with log 2 FC > 1.50), and one of the genes was predicted to yield a polyamine transporter protein (log 2 FC ≈ 1.59), which may furthermore be associated with detoxification [20]. There were 18 significantly underexpressed genes in CT1 whose log 2 FC values fell below -1.5 (Table 2; Figure 2B). Various proteins predicted to possess diverse functions declined in expression following azole treatment, such as a prion associated with epigenetic activity (g4209). However, three protein transporters exhibited large reductions in gene expression: a predicted efflux pump (g12477; log 2 FC ≈ −1.61), a putative transporter (g11475, log 2 FC ≈ −1.77), and a potential MFS transporter (g11706; log 2 FC ≈ −2.77).  . Hypothetical protein was abbreviated to Hyp. protein for simplicity. Note that the Cps isolates had over double the number of differentially underexpressed genes than overexpressed genes, some of which were associated with efflux and detoxification activity. The genes are arranged according to increasing p-value. The data was plotted using GraphPad Prism version 9.1.2 for Windows (GraphPad Software; San Diego, CA, USA, www.graphpad.com).
The overexpressed genes in tetraconazole-treated Cps CT1 possessed a variety of different GO terms, though only the GO annotations associated with the membrane/integral component of the membrane appeared relatively common (occurring in eight of the 16 genes) ( Figure 3A). Tetraconazole treatment appears to impact a large array of proteins in Cps CT1, such as proteins putatively associated with the nucleus (GO label occurring in five of 55 genes) or carrying out metal/zinc binding (occurring three and five times, respectively) and redox reactions (seen six times among 55 genes) ( Figure 3B). Similar to the overexpressed genes, tetraconazole appears to significantly impact genes annotated to synthesize membrane proteins (observed in 19 of the 55 underexpressed genes).

Differentially Expressed Genes in Tetraconazole-Treated Che
Seventeen genes were statistically, differentially overexpressed in tetraconazole-treated CB45 compared to untreated CB45 (Table 3), and another 17 genes were significantly, differentially underexpressed in the tetraconazole-treated Che CB45. We observed that 22 of these 34 differentially over-and underexpressed genes successfully matched to a gene in Cps CT1 in our multi-blast analysis ( Figure 1B). No gene underexpressed in tetraconazole-treated CB45 matched to a strongly overexpressed gene (log 2 FC > 1.5) in tetraconazole-treated CT1 ( Figure 1B); however, seven of the overexpressed genes in Che did match Cps genes whose expressions were underexpressed (log 2 FC < −1.5) or unaffected (log 2 FC ≈ 0) following azole treatment ( Figure 1B). In only one of these matches, the Che gene, the putative FAD-dependent oxidoreductase (g8987), had a log 2 FC exceeding 1.5, demonstrating that expression of this gene was particularly impacted by tetraconazole treatment in Che compared to the closest matching gene in Cps.   Eight of the overexpressed CB45 genes had log 2 FC > 1.50. The most overexpressed gene with a predicted function encoded a polyamine transporter 4 protein (log 2 FC ≈ 1.97), followed by a cytochrome P450 monooxygenase (log 2 FC ≈ 1.94) and a cytochrome-b5 reductase (log 2 FC ≈ 1.80) (Table 3; Figure 2C). Several of these genes were thus annotated to have functions generally associated with xenobiotic detoxification activity [20]. Among the underexpressed genes in Che, four had log 2 FC values below −1.5, and only two (a predicted tyrosinase, g7560, according to our conserved domain search and a proteinase) had predicted function (Table 3; Figure 2D). Though possessing a log 2 FC value of roughly −1.4 and failing to reach our thresholds, there was a gene predicted to be associated with aminotriazole resistance that was underexpressed in Che.
When focusing on the GO annotations for these differentially over-and underexpressed genes, we noted that overexpressed genes in tetraconazole-treated CB45 tended to be associated with oxidoreductase activity (GO label occurred seven times among 17 genes) and redox processes (occurred in 10 of the 17 genes) ( Figure 3C). A relatively high prevalence of overexpressed genes also tended to be linked to the cell membrane (occurred in five of 17 genes). For the underexpressed genes, the only GO names that prevalently re-occurred were those associated with the membrane or the integral component of the membrane (both being observed in seven of the 17 genes). Hence, tetraconazole exposure appeared to largely reduce expression of select genes associated with the cellular membrane.

CYP51A and CYP51B Expression Levels
In Cps isolate CT1, reads were mapped to the CYP51A pseudogene (g5325), though expression was low and statistically insignificant between the control and tetraconazoletreated cultures (CPM < 2; p-value > 0.05). The CYP51B gene (g5734) was considered significantly overexpressed in the treated culture based on p-value (<0.05), but the log 2 FC value failed to reach our threshold of 1.5 (~0.61).

Prediction of Ergosterol Biosynthesis Genes
There were nine genes in CB45 predicted to play a role in ergosterol biosynthesis using our GO name search (Table 4). Two pairs of these genes were annotated to fulfill the same function (cytochrome P450 61; NADPH-cytochrome P450 reductase) and matched to the same genes in Cps according to our multi-blast analyses (g3336 and g9924, respectively). One of these genes, the predicted cytochrome P450 61 g133, had 100% identity to Che g1530 and Cps g3336 but was halted by a putative premature stop codon, upstream of predicted protein domains in the amino acid sequence of g1530. Although reads mapped to this gene, it is possible that Che g133 encodes a pseudogene unable to contribute to ergosterol biosynthesis, which may be complemented by g1530 in Che. Using data from [17] and NCBI blastp searches against Fusarium species, we identified an additional 17 predicted ergosterol biosynthesis genes in CB45 (Table 5). Notably, only one gene in the sensitive species, Cps g12006 (Table 4), showed significant overexpression in tetraconazole-treated CT1 (Table 1; Figure 1A, Figure 2A). No predicted ergosterol biosynthesis gene in CB45 show significant differences in gene expression.

Discussion
Boxwood blight is an economically important disease caused by two sister fungal species, Cps and Che, in landscape and plant nurseries across the globe [1,2]. As documented in a previous study [3] and reaffirmed in the analysis presented herein (Table 1), these two species differentially respond to select azole fungicides, the latter demonstrating a five-fold greater EC 50 tolerance to tetraconazole, and significantly reduced sensitivity to propiconazole, both members of the azole class (Table 1). In this study, we applied transcriptomic and bioinformatic analyses to gain insight into possible genetic determinants for the observed differences in azole fungicide sensitivity.
Results from our previous study showed that a CYP51A pseudogene in Cps may explain differential sensitivity to azoles [16]. As expected, few reads and no statistically significant difference in gene expression were detected for the CYP51A pseudogene in Cps. Given the stringent mapping parameters applied here, and the observation that the few reads that were mapped to the CYP51A pseudogene mapped to the region upstream of the premature stop codon, we can assume that a lack of functional CYP51A protein expression could result in less compensatory ergosterol synthesis. By contrast, in Che, the expression of both CYP51A and CYP51B gene paralogs could result in greater ergosterol synthesis, as both CYP51A and CYP51B proteins are known to catalyze the same step in the ergosterol pathway. In Candida albicans, chromosomal duplications of CYP51 are linked to fluconazole-resistance [29]. Numerous studies demonstrate overexpression of CYP51 as a mechanism for resistance to azoles (reviewed in [30]). Future complementation studies with a functional CYP51A in Cps are warranted to carefully assess the impact of the CYP51A pseudogene.
Although CYP51 and its paralogs receive the majority of focus when it comes to azole fungicide resistance, the synthesis of ergosterol is a complex metabolic process driven by over 20 genes from as many as three separate biochemical pathways [17,18]. Here we identified potential ergosterol biosynthesis genes in Cps and Che, though, intriguingly, only one biosynthesis gene, observed only in our sensitive Cps isolate (g12006, Tables 2 and 4), showed significant overexpression following one hour after treatment. This finding sug-gests that Che is not reliant upon upregulation of ergosterol biosynthesis following exposure to offset the inhibitory effects of tetraconazole. Interestingly, our observations differ from those observed by [31] in Cochliobolus sativus. This pathogen, which exhibits emerging fungicide resistance, displayed gene expression patterns whereby select ergosterol biosynthesis pathway genes were upregulated within hours of exposure. Our results in Che may instead support the hypothesis that the presence of an additional, functional CYP51 paralog in Che could confer the increased tolerance towards azoles relative to Cps [16]. In Schizosaccharomyces pombe, it was demonstrated that the addition of a gene encoding a C-5 sterol desaturase participating in ergosterol biosynthesis increased synthesis of ergosterol and was correlated with improved tolerance towards exposure to extreme temperatures and the presence of ethanol [32]. In light of the differential thermotolerance previously reported in Cps and Che, comparative analyses of lipid profiles and ergosterol levels are warranted in these two species [3].
Analysis of RNA-Seq expression profiles of tetraconazole-treated Cps and Che revealed overexpressed genes predominantly associated with oxidoreductase activity, oxidationreduction processes, and monooxygenase activity, somewhat parallel to findings from previous transcriptomic studies of azole sensitivity in fungi (Tables 3 and 4) [17,21,31]. The overexpression of cytochrome P450 monooxygenases is implicated in phase one of xenobiotic detoxification, two of which have been functionally shown to play a role in multidrug resistance in the turfgrass fungal pathogen, Clarireedia jacksonii (synonym = Sclerotinia homoeocarpa) [20,21,33]. Furthermore, our candidate detoxification genes and their differential expression patterns in Che and Cps were consistent with findings from previous studies in both Cercospora beticola and Clarireedia jacksonii [17,21]. The putative detoxification genes identified in the current study, including subsets of predicted transferases (Tables S5 and S6), could provide a more comprehensive list of putative targets for functional characterization and/or serve as genetic markers for monitoring the evolution of resistance and potentially guiding future fungicide development. It should be noted that gene deletion experiments in Cercospora beticola to verify the role of candidate detoxification genes were inconclusive, while in Clarireedia jacksonii experiments validated the function of two novel cytochrome P450 monooxygenases [17,21]. The amenability to genetic transformation of Cps and Che remains unknown, but findings presented here serve to build a knowledge base to enable such studies.
As reported by [21], non-synonymous polymorphism in a zinc transcription factor led to gain-of-function resulting in overexpression of detoxification genes from the aforementioned gene families of the fungal xenome. The direct comparisons of gene expression patterns relied on the apparent clonal nature of the sensitive and resistant fungal isolates of Clarireedia jacksonii, which also facilitated screening for transcription factor polymorphism. However, in the current study, the genetic divergence and the higher levels of genetic polymorphism between Cps and Che (warranting their placement in two separate species [3]) hampered efforts to readily identify polymorphism that could be linked to differential expression of an assortment of genes representing multiple gene families. It should be noted that a homolog of a negative regulator of nitrate metabolism genes, nmrA (g4897), was underexpressed in tetraconazole-treated Cps (Table 2; Figure 1A, Figure 2B). Though the linkage of this gene to azole sensitivity has not been previously reported, nmrA deletion mutants have been shown to display reduced sensitivity to rapamycin and methylmethanesulfonate, presumably through rapamycin signaling and DNA-damage [34,35]. Another candidate gene from Cps that could serve to influence expression levels of downregulated genes shows similarity to RNQ1 (Table 2; Figure 1A, Figure 2B), a prion encoding gene that is thought to have a role in epigenetics and phenotypic changes in Saccharomyces cerevisiae, and which also mediates polyQ aggregation in concert with heat shock proteins in a Huntington's Disease yeast model [36,37]. The significance of these genes in the response to an inhibitory concentration of tetraconazole by Cps remains to be determined.
Numerous genes were downregulated in Cps in response to tetraconazole treatment (Table 2; Figure 1A, Figure 2B), which may reflect the shunting of energy from housekeeping and maintenance to xenobiotic detoxification, analogous to findings from [21], which found a several-fold greater abundance of genes downregulated in the sensitive isolates following exposure to azole exposure, including genes encoding proteins with predicted association with cell membranes. In both [21] and the current study, the sensitive isolate/species were exposed to an azole concentration of five times the EC 50 or above, presumably a lethal dose, which could explain the much greater number of downregulated genes and the fact that equivalent genes in Che tended to be unaffected/upregulated ( Figure 1A). Further, several of these genes were predicted to encode efflux transporter proteins (Table 2; Figure 1A, Figure 2B) [13,20,[25][26][27][28]. Given the sensitivity of Cps CT1 to tetraconazole and its aforementioned exposure to an apparent lethal dosage of the chemical, fungal cells may have in response triggered reduction in expression of suites of genes unrelated to azole detoxification.
In Che, several overexpressed genes are associated with the first and third phases of detoxification, with GO terms suggesting association with the cell membrane and oxidoreductase activity (Table 3; Figure 1B, Figure 2C, Figure 3C). A similar pattern is reported in Cochliobolus sativus, Cercospora beticola and Clarireedia jacksonii which, following azole treatment, display significant increases in gene expression for phase I cytochrome P450 proteins and phase III ABC-G transporters [17,21,31]. Of particular interest is a cytochrome P450 gene described as a benzoate-4-monooxygenase, which showed significant overexpression (g8537). A previous study on antifungal response in Fusarium oxysporum and Neocosmospora solani found a benzoate-4-monooxygenase to be overexpressed following treatment with amphotericin b and posoconazole [38]. This gene, denoted CYP53, is a member of the family of CYP450 enzymes exclusive to fungi that are involved in phenolic detoxification. A study looking into natural compounds to target Cochliobolus lunatus found that natural phenolic compounds inhibiting CYP53A15 activity showed antifungal activity [39]. Perhaps this gene target may serve to inform novel fungicide treatments in boxwood blight.
Few genes were underexpressed in Che following treatment, and, of note, included several genes with predicted peptidase activity (Table 3; Figure 1B, Figure 2D). Aspartyl peptidase activity has been shown to vary in response to azoles in a Candida albicans, though future studies are necessary to begin to elucidate how peptidase activity can influence azole sensitivity in fungi [40]. According to our conserved domain search, a predicted tyrosinase gene (g7560) showed the highest degree of underexpression following tetraconazole exposure, and though azole chemicals are known inhibitors of tyrosinase activity, their impact on gene expression remains unclear [41]. A third gene of interest, a predicted short-chain dehydrogenase, showed a minimal degree of underexpression following tetraconazole exposure. In Aspergillus fumigatus, an important human pathogen, deletion mutants for HorA (a mitochondrial short-chain dehydrogenase) showed increased tolerance to azoles, suggesting to be indicative of a metabolic state linked to mitochondrial impairment [42]. Taken together, diverse metabolic pathways appear to be downregulated in Che to facilitate azole tolerance.

Estimating Effective Fungicide Concentrations for 50% and 85% Inhibition of Growth In Vitro of Cps and Che
To identify fungicides to which Cps and Che are differentially tolerant, in vitro radialgrowth assays were used with a range of concentrations of fungicides, as listed in Table S1. All assays were performed using Cps specimen CpsCT1s (hereafter: CT1) and Che type specimen CB045 (hereafter: CB45); for tetraconazole and propiconazole, an additional isolate of each species was used, as indicated in Table S1. All assays used commercial formulations of fungicides (Table S1), added to autoclaved half-strength potato-dextrose agar (half-PDA) in amounts producing the concentration ranges specified in Table S1, based on labeled percentages of active ingredient. After thorough mixing, exactly 20 mL of the amended media (or unamended control media) were then poured into 100 × 15 mm plastic Petri dishes, which were then used the following day for inoculation with test isolates. Inoculations of amended media used 10-day-old cultures growing on Spezieller Nährstoffarmer agar (SNA) [43], using 5.0-cm diameter plugs taken with a cork borer from culture margins and placed in the center of the dish. Two replicates were used for each isolate x treatment/control combination. After sealing plates with Parafilm M ® (Amcor, Zurich, Switzerland), cultures were incubated at room temperature in the dark until control cultures (i.e., those not amended with fungicide) reached no less than 1.0 cm from the dish perimeter, for the number of days as indicated in Table S1. Probit analysis was used to determine the effective concentrations of fungicide necessary for 50% and 85% inhibition of growth, as follows. Diameters were measured using digital calipers (Mahr GmbH, Esslingen am Neckar, Germany), with measurements automatically uploaded into a spreadsheet. After subtracting the inoculation-plug diameter, growth rates were calculated by dividing adjusted diameters by the number of growth days, then normalized by dividing by the growth rate of the unamended control. These normalized growth rates, expressed as percentages, were then transformed to "probits" (probability units) [44]. The linearized regression of probits against log-transformed fungicide concentrations was then used to calculate EC 50 and EC 85 values, as listed in Table S1.

Genome Sequencing and Assembly
The genome sequences and predicted genes of Cps isolate CT1 (also known as isolate CBS 139,707 and cpsCT1) and Che isolate CB45 (also known as CBS 138,102 and CB045) reported in [16] were used to guide our differential expression analysis.

RNA Sequencing of Control versus Tetraconazole-Treated CT1 and CB45
Using 10 mL of half-strength PDB, conidia were washed from ten-day-old cultures of CT1 and CB45, growing on half-strength PDA that had been scored four days previously with a sterile scalpel to promote sporulation. After estimating conidia concentrations with a hemocytometer, 1000 conidia were plated onto half-strength PDA and measured the following day for viability. Based on those estimates, 10,000 putatively viable conidia were used to inoculate each of the sixteen 50-mL tubes containing 20-mL of half-strength PDB, placed at a 45 • angle to maximize surface area, and allowed to grow at room temperature without shaking for six days. The resulting mycelial mats were then transferred with sterile spatulas to sixteen 50-mL tubes containing 20 mL of half-strength PDB, half of which had been amended with 1.0 ppm tetraconazole (Minerva). After incubation on a rocking platform for one hour, the mycelial mats were blotted on Kimwipes over ice, then transferred to 2-mL tubes, plunged in liquid N 2 , and lyophilized. The resulting lyophilized mycelium was stored at −80 • C for one day before RNA extraction using the Ambion Trizol Plus kit (ThermoFisher Scientific, Waltham, MA, USA). Resulting RNAs were then processed at the Yale Center for Genome Analysis (West Haven, CT, USA), where they were first analyzed for quality and quantity on a Bioanalyzer (Agilent Technologies, Inc., Santa Clara, CA, USA), confirmed to have RIN scores of 9.0 or greater, and then processed for 2 × 75-bp next-generation sequencing on an Illumina HiSeq-2500 (Illumina, Inc., San Diego, CA, USA).

Genome Functional Annotation and Multi-Blast
Unless specified otherwise, the following procedures were performed using OmicsBox software [45]. For both Calonectria genomes, the protein sequences of each predicted gene were extracted using our gene annotation files. A blastp-fast search was performed in September 2019 for each gene using the OmicsBox CloudBlast feature and the nr_v5 nonredundant protein sequence database using the following parameters: searches limited to fungi (taxa: 4751); e-value ≤ 10 −10 ; word size = 6; retaining only the top 10 blast hits for each amino acid sequence. A blast description annotator and a filter removing matches with regions of low complexity were applied. After performing the CloudBlast search, an InterProScan analysis (CloudIPS) acquired protein domain information for each protein in our files [46]. All available databases other than PatternScan, Coils, and MobiDBLite were selected for this analysis. Finally, Gene Ontology (GO) mapping was performed on each protein sequence in the ensuing file using the most recently updated GO database available [46].
The amino acid sequences for predicted genes in CB45 and CT1 were extracted using OmicsBox [45]. Two multi-blasts (blastp) were performed using the CLC Genomics Workbench 11.0 (QIAGEN, 2017), the first querying CB45 genes against CT1 and the second querying CT1 genes against those of CB45. The same parameters were used for both BLAST searches (match = 2, mismatch = -3; gap existence cost = 5, extension cost = 2; e-value ≤ 10 −10 ; word size = 11; number of threads = 10; filter low complexity matches). The two multi-blasts allowed for verification of whether we could find reliable matches of equivalent CB45 genes in the CT1 genome.

Identification of Candidate Detoxification Genes
Genes with protein domains associated with the three phases of detoxification (oxidation, conjugation, and efflux) were mined from the CT1 and CB45 annotation files [20]. Amino acid sequences possessing any of the following four InterPro IDs were extracted from the annotation files for both species: Cytochrome P450 (IPR001128; oxidation/ breakdown phase); Major Facilitator Superfamily (IPR020846; MFS, efflux phase); PDR-like Subfamily G, Domain 1 (IPR034001; ABC-transporter, efflux phase); and PDR-like Subfamily G, Domain 2 (IPR034003; ABC-transporter, efflux phase). Transferases were extracted by filtering the GO-annotated CT1 and CB45 genomes for the GO function name transferase activity. The full lists of these genes are stored as Supplementary Tables S3-S6.

RNA-Seq Analysis of Control versus Tetraconazole-Treated Isolates
Three biological replicates were sequenced for control and tetraconazole-treated CT1 and CB45. Each biological replicate had three technical replicates. In total, there were 36 experimental units, each with a minimum of 4,000,000 total RNA-Seq reads.
All of the following procedures in this section were performed using the OmicsBox software [45]. RNA-Seq reads were converted into FASTQ file format and read file quality was assessed and poor reads were filtered out using the default parameters of FastQC version 0.11.8 [47]. RNA-Seq reads meeting the quality threshold for CT1 and CB45 isolates were, respectively, aligned to the CT1 and CB45 reference genomes using the STAR aligner software in the OmicsBox environment [48], with read overhang set to one base pair lower than the maximum read length (75 base pairs) following the manufacturer's protocol. The annotation files of the two genomes were used to help guide the alignment. The assigned minimum and maximum intron lengths were 20 base pairs and 2000 base pairs, respectively [49]. Only one mismatch between the RNA-Seq read and the reference genome was permitted, and reads were allowed to match with only one region in the genome. These stringent mapping parameters were used to diminish the likelihood of false matching between RNA reads and predicted genes.
Using the HTSeq 0.9.0 package, two read count tables were generated, one comparing control and treated CT1 RNA-Seq reads and the other comparing control and treated CB45 reads [50]. The CT1 and CB45 structural annotation files were respectively used to acquire the RNA-Seq read counts per gene. The counts were made according to read overlap with the exon region of each gene (overlap mode: union; strand specificity: non strand specific). Since Illumina RNA-Seq procedures generate high-quality reads with low error rates (less than one base in 500), a high minimum mapping quality value of 25 (or an error rate of one in 500 bases) was specified for this step.
Using the edgeR software package within the OmicsBox environment and our count tables, separate differential expression analyses were run for CT1 and CB45 [51]. In both instances, a count per million (CPM) filter was set to one and at least one sample needed to reach this criterion. These parameters indicated that, for every million reads in a technical replicate, at least one of the reads must map to a particular gene for that gene to be included in the analysis. On average, there were between three and four million reads which mapped to a reference genome per replicate. Thus, for a given gene to be included in this analysis, at least three reads were needed to map to that gene for a minimum of one replicate. Since there were many genes whose samples possessed two or fewer matches, this criterion was considered suitable to acquire reliable calculations of gene differential expression. Additionally, the strict RNA-Seq alignment conditions decreased the likelihood that any of the 75-bp reads were mapped to similar regions of incorrect genes, thereby increasing confidence that only genes with reliably mapped reads were included in this analysis. However, we had also run each differentially over-and underexpressed gene against the NCBI Conserved Domain Search tool to further check the predicted functions of these genes. Non-matching genes were described as encoding a hypothetical protein. GraphPad Prism version 9.1.2 for Windows (GraphPad Software; San Diego, CA, USA, www.graphpad.com) was used to generate plots displaying normalized gene expression values and GO annotation counts for the differentially expressed genes in Cps CT1 and Che CB45 treated with tetraconazole.

Identification of Ergosterol Biosynthesis Genes
A list of putative ergosterol biosynthesis genes were extracted from our annotated Che CB45 genome using the GO name ergosterol biosynthetic process. A protein blast was also performed against the CB45 genome (gap existence cost = 11; gap extension cost = 1; e-value = 10 −10 , word size = 3; number of threads = 4; filter low complexity matches; BLO-SUM62) using various NCBI reference ergosterol biosynthesis genes from Fusarium species (Table S6). This step was performed to try and capture putative ergosterol biosynthesis genes missing from our GO name search. We used the pathway discussed by [17] to help us guide our NCBI gene search. Using the two aforementioned multi-blasts (see Genome functional annotation and multi-blast), we documented the closest matching genes in Cps CT1.

Conclusions
The findings presented here contribute to the body of knowledge on azole resistance determinants in filamentous fungal plant pathogens, and though previous transcriptomics studies of this nature have undertaken similar methodologies, the novelty of this study lies in part in the phylogenetic relationship of Cps and Che, two distinct sister species, and the presence of a CYP51A pseudogene in Cps. Follow up studies may have relevance to novel azole resistance mechanisms in human fungal pathogens [6]. Many medically important fungal pathogens, including Aspergillus fumigatus and Candida albicans, have shown reduced sensitivity and/or resistance towards azole antifungals. Candida auris, a species formally described in 2009 [52,53] is comprised of distinct clades, and treatment for disease from this pathogen is complicated by azole resistance and resistance to other major groups of medical antifungals through target gene mutations, gene duplications, and increased expression of efflux proteins. As new reports of antifungal resistance continue to emerge, an increasing need for novel drug targets arises. Future transcriptomic studies of this nature will provide a roadmap to facilitate development of improved fungicides and antifungal drugs for medicine, while serving to uncover cross-species patterns and novel mechanisms underpinning multidrug resistance in fungi.

Conflicts of Interest:
The authors declare that the research was conducted in the absence of any commercial/financial relationships that could be construed as a conflict of interest.