Dual RNA-Seq Analysis of the Pine-Fusarium circinatum Interaction in Resistant (Pinus tecunumanii) and Susceptible (Pinus patula) Hosts

Fusarium circinatum poses a serious threat to many pine species in both commercial and natural pine forests. Knowledge regarding the molecular basis of pine-F. circinatum host-pathogen interactions could assist efforts to produce more resistant planting stock. This study aimed to identify molecular responses underlying resistance against F. circinatum. A dual RNA-seq approach was used to investigate host and pathogen expression in F. circinatum challenged Pinus tecunumanii (resistant) and Pinus patula (susceptible), at three- and seven-days post inoculation. RNA-seq reads were mapped to combined host-pathogen references for both pine species to identify differentially expressed genes (DEGs). F. circinatum genes expressed during infection showed decreased ergosterol biosynthesis in P. tecunumanii relative to P. patula. For P. tecunumanii, enriched gene ontologies and DEGs indicated roles for auxin-, ethylene-, jasmonate- and salicylate-mediated phytohormone signalling. Correspondingly, key phytohormone signaling components were down-regulated in P. patula. Key F. circinatum ergosterol biosynthesis genes were expressed at lower levels during infection of the resistant relative to the susceptible host. This study further suggests that coordination of phytohormone signaling is required for F. circinatum resistance in P. tecunumanii, while a comparatively delayed response and impaired phytohormone signaling contributes to susceptibility in P. patula.


Introduction
One of the most important pathogens to natural and industrial pine forests is the pitch canker fungus, Fusarium circinatum Nirenberg and O'Donnell [1][2][3]. Since the first report of the pathogen, and the disease it causes, in the southeastern US, F. circinatum has been identified in more than ten countries world-wide, resulting in significant losses in both nurseries and plantations [2][3][4]. The pathogen was originally identified due to the visible symptom development on pines and classified as a necrotrophic pathogen of pine, capable of infecting douglas fir (Pseudostuga menziesii) and up to 60 different species of pine, including economically important species such as Pinus pinaster, Pinus taeda and Pinus radiata [3]. Alarmingly, the last natural stands of P. radiata are under threat of extirpation due to pitch canker [5]. Multiple studies, however, have since shown a broader range of potential ecological activities [6]. The pathogen has been shown to endophytically infect pine [7], as well as certain grass species [8], maize [9] and some dicots [10]. Thus, there is a high likelihood for F. circinatum to spread into unaffected regions through naturally occurring inoculum reservoirs on unmonitored species [2,[7][8][9]. from inoculated samples. Comparably for P. tecunumanii, an average of 0.02 ± 0.01% and 0.30 ± 0.17% reads mapped to the F. circinatum transcriptome from mock-inoculated samples, with 0.16 ± 0.04% and 1.62 ± 0.81% reads mapped from inoculated samples.
Filtering of the expression data identified 25,000, 20,614 and 5,003 expressed genes for P. patula, P. tecunumanii and F. circinatum, respectively, across all samples. Subsequent filtering of F. circinatum expression resulted in 4,354 high-confidence expressed genes ( Table 1, Table S5). Differential expression (DE) analysis of inoculated samples (P. tecunumanii versus P. patula) identified 132 and 470 significant F. circinatum differentially expressed genes (DEGs) for 3-and 7-dpi, respectively (Table 1, Table S5). P. patula DE analysis identified 323 and 7453 significant DEGs (inoculated versus mock-inoculated) at 3-and 7-dpi, while 735 and 2499 significant DEGs were identified for P. tecunumanii (Table 1, Tables S6  and S7).  b Host genes differentially expressed in inoculated relative to mock-inoculated host expression data. c F. circinatum genes differentially expressed in inoculated relative to mock-inoculated samples in the full expression data set (including both host and pathogen mapped reads) for each host. Up-regulated genes represent high confidence F. circinatum expressed genes. Down-regulated genes were excluded from downstream analysis. d F. circinatum genes differentially expressed in P. tecunumanii relative to P. patula inoculated samples from pathogen expression data.

Over-Represented GO Terms within Pathogen Datasets
GO enrichment analysis in the high confidence expressed F. circinatum genes showed shared biological process (BP), cellular compartment (CC) and molecular function (MF) terms between all data sets related to ribosomes, translation and lipid metabolism, indicative of growth, as well as BP terms related to responses to farnesol ( Figure S3, Table S8). BP terms related to pectin hydrolysis were only enriched in the 3-dpi P. patula data set. Most enriched terms in the P. tecunumanii 3-dpi data set were also enriched in both the P. tecunumanii and P. patula 7-dpi data sets. This included CC terms related to membranes and mitochondria, MF terms related to hydrolysis, lipid binding and oxidoreductases, and BP terms related to responses to oxidative stress and respiration ( Figure S3, Table S8). For both hosts at 7-dpi, there were also many enriched BP terms related to localization, cell wall organization, alcohol biosynthesis and sterol biosynthesis.
Further GO enrichment analysis of F. circinatum DEGs between P. tecunumanii and P. patula at both timepoints was performed to identify differences in pathogen responses between hosts (Table S9). Few terms were enriched for F. circinatum genes up-regulated during infection of P. tecunumanii relative to P. patula. No terms were enriched in the 3-dpi data set, while there was enrichment for the CC term, extracellular region, and the BP terms, carbohydrate metabolic process, and, oxidationreduction process, in the 7-dpi data set. In contrast, many terms were enriched for F. circinatum genes down-regulated in inoculated P. tecunumanii relative to P. patula samples. At both time points, cytoplasmic translation terms were enriched in all three GO categories (Table S9). At 3-dpi, BP terms related to glycolysis and energy production were enriched, while at 7-dpi MF and BP terms related to sterol and alcohol biosynthesis were enriched. Due to the enrichment of terms related to ergosterol biosynthesis in the high confidence 7-dpi data set from both hosts, as well as in the down-regulated 7-dpi data set, this pathway was investigated further.

Transcriptional Responses Related to Ergosterol Biosynthesis in the Pathogen
Candidate genes for all ergosterol biosynthesis steps could be identified in the F. circinatum transcriptome except HMG-CoA synthase (Table S10). Due to the physiological importance of HMG-CoA synthase, this likely indicates incompleteness of the genome rather than absence of the gene. ERG10 and ERG20 showed significantly lower expression during infection of P. patula than P. tecunumanii at 7-dpi ( Figure 1). While this indicates that F. circinatum produces more farnesyl diphosphate (FDP) during infection of P. tecunumanii, FDP is the precursor for a wide array of metabolites. Conversely, five genes involved in late ergosterol biosynthesis, CYP51/sterol-14α-demethylase (ERG11), C-4 methylsterol oxidase (ERG25), sterol 24-C-methyltransferase (ERG6), C-8 sterol isomerase (ERG2) and δ7-sterol 5-desaturase (ERG3), were expressed at lower levels during infection of P. tecunumanii than P. patula at 7-dpi ( Figure 1).

Over-Represented GO Terms within Host Datasets
Few GO terms were enriched for P. patula DEGs at 3-dpi (Table S11). Analysis of up-regulated DEGs identified the CC terms, cell wall and external encapsulating structure, as well as the BP terms, syncytium formation, cytoplasmic translation, response to oxygen-containing compound and response to chitin. All DEGs underlying the enriched syncytium formation term were predicted expansins, proteins involved in cell wall relaxation. The only enriched GO term in the down-regulated DEGs at 3-dpi was the MF term xyloglucan:xyloglucosyl transferase activity, an enzymatic activity involved in cell wall reinforcement.
For up-regulated DEGs at 7-dpi, most enriched CC terms were related to the nucleus, mitochondria, cytoplasm and ribosomes (Table S11). Enriched MF terms included many terms related to ribosomal activity and control of reactive oxygen species (ROS) during a defence response. Many BP terms were enriched, including a wide array of known defence related terms. Among these were terms related to: apoptosis, ROS production, response to oxidative stress, increased cytokine production, terpenoid biosynthesis, the lipoxygenase pathway, response to chitin (similar to the 3-dpi data set), MAPK signaling, responses to salicylic acid (SA), jasmonic acid (JA) and ethylene (ET), as well as JA-mediated induced systemic resistance (ISR) signaling. All BP terms enriched at 3-dpi were also enriched at 7-dpi except syncytium formation, however, a larger complement of expansins were up-regulated at 7than 3-dpi. For down-regulated genes, the enriched CC terms were all plastid related and the only enriched MF term was beta-amylase activity. The enriched BP terms included plastid organization and plastid processes, such as starch metabolism, apoptosis and hypersensitive response (HR), as well as SA responses and SAR.
A large number of BP terms were enriched in the P. tecunumanii up-regulated DEGs at both time points (Table S12). Potential defence related BP terms enriched at both timepoints were related to response to chitin, ubiquitin-mediated proteolysis (UMP) and vesicle mediated transport. Enriched terms unique to the 3-dpi dataset were related to UMP, ET signaling, JA mediated ISR and terpenoid biosynthesis. The enriched terms unique to 7-dpi were related to the SA response, SA-mediated SAR, JA responses, ET biosynthesis and the biosynthesis of various phytoalexins, including camalexin.
Enriched GO terms for down-regulated P. tecunumanii DEGs at both timepoints were mainly related to cellular growth and replication (Table S12). Most enriched BP terms for down-regulated DEGs at 3-dpi were also enriched for down-regulated DEGs at 7-dpi.

Transcriptional Responses Related to Host Phytohormone Signalling
Phytohormones play crucial roles during growth and development. Interactions between these hormones also serve to regulate gene expression during stress responses. Hormone biosynthesis and signaling related DEGs were investigated to identify putative pathways involved in the pine-Fusarium interaction due to the enrichment of different phytohormone related GO terms in the up-regulated data at 3-and 7-dpi in P. tecunumanii, as well as in both the up-and down-regulated data at 7-dpi in P. patula ( Figure 2, Tables S13 and S14).

Cytokinin
At both time points in P. tecunumanii, there was up-regulation of UDP-glycosyl transferase (UGT) and cytokinin oxidase/dehydrogenase (CKX) genes, enzymes related to cytokinin (CK) degradation, with more UGTs up-regulated at 7-dpi ( Figure 2, Table S14). There was also down-regulation of cytochrome P450s (CYP) related to CK biosynthesis and histidine kinase (HK) receptor genes at 7-dpi. Increased degradation and down-regulation of biosynthesis indicate suppression of CK at both time points. In P. patula, the only CK related DEGs at 3-dpi were up-regulated UGT genes, while at 7-dpi, similar to P. tecunumanii, there was up-regulation of UGT and CKX genes as well as down-regulation of HK genes ( Figure 2, Table S13). Furthermore, two A-ARR genes, negative regulators of CK signaling, showed up-regulation and two B-ARR genes, positive regulators of CK signaling and A-ARRs, showed down-regulation in P. patula at 7-dpi. This, together with the up-regulation of superoxide dismutase (SOD) genes, indicate suppression of CK signaling.

Gibberellic Acid
At 3-dpi in P. patula, the only gibberellic acid (GA) related DEG was an up-regulated GA methyl transferase (GAMT, Figure 2, Table S13). Unlike other phytohormones, methylation results in GA degradation (Eckardt, 2007). GA 3-oxidase (GA3ox) and putative ent-copalyl diphosphate synthase (CPS) genes were up-regulated at both time-points in P. tecunumanii, as well as at 7-dpi in P. patula indicating possible GA biosynthesis at 7-dpi in both hosts. However, there were also down-regulated CPS genes, with more CPS genes down-regulated in P. patula than P. tecunumanii. GID1, GA2ox, LBD and AMY genes were up-regulated in both hosts at 7-dpi. A LFY gene was up-regulated at 7-dpi in P. patula and phytochrome-interacting factor 3 (PIF3) genes were down. While up-regulation of LFY and AMY genes indicate the presence of GA signaling between 3 and 7-dpi, up-regulation of DELLA responsive GA2ox, GID1 and LBD indicate suppression of GA signaling at 7-dpi in both hosts.
Microorganisms 2019, 7, x FOR PEER REVIEW 8 of 20 Figure 2. Summary of phytohormone related host DEGs during F. circinatum challenge. Up-and down-regulation (inoculated relative to mock-inoculated) of genes related to cytokinin, abscisic acid, gibberellic acid, brassinosteroid, ethylene, jasmonic acid, salicylic acid and auxin signalling, in P. patula (blue) and P. tecunumanii (green) at three (dark colours) and seven (light colours) days post inoculation, are indicated by arrows (↑ and ↓ respectively; Additional file 5:

Abscisic Acid
In P. tecunumanii a 9-cis-epoxycarotenoid dioxygenase (NCED) was down-regulated at both time points, however, at 7-dpi another NCED, as well as a putative xanthoxin dehydrogenase (SDR) was up-regulated and an abscisic acid (ABA) hydroxylase was down-regulated ( Figure 2, Table S14). This could indicate suppression of ABA biosynthesis at 3-dpi and induction at 7-dpi in the resistant host. In P. patula, an NCED and an ABA hydroxylase were up-regulated at 3-dpi. At 7-dpi, although an SDR1 was up-regulated, other ABA biosynthesis genes, zeaxanthin epoxidase (ZEP) and abscisic-aldehyde oxidase (AAO3), were down-regulated while an ABA hydroxylase and a carotenoid cleavage dioxygenase 8 (CCD8) were up-regulated. CCD8 diverts carotenoid metabolism away from Zeaxanthin. This indicates suppression of ABA biosynthesis and increased degradation at 7-dpi in the susceptible host.
ABA receptors were up-regulated at both time points in P. tecunumanii and at 7-dpi in P. patula. A type 2C protein phosphatase (PP2C) gene was down-regulated at 7-dpi in P. patula. The up-regulation of receptors and a CIPK20 could indicate ABA signaling at both time points in P. tecunumanii. At 7-dpi in P. patula, down-regulation of ABA biosynthesis indicates suppression of ABA levels, while up-regulation of receptor and CIPK20 genes with down-regulation of a PP2C suggest ABA signaling.
Cross talk between CK and ABA signaling is mediated by ABA-insensitive 4 (ABI4) and ABI5, positive regulators of ABA signaling, and A-ARRs. ABI5 activity is attenuated by interaction with A-ARRs, allowing A-ARRs to negatively regulate both CK and ABA signaling [32]. ABI4 positively regulates A-ARR5, resulting in suppression of CK responses by ABA signaling [33]. An ABI4 gene was down-regulated at 7-dpi in P. tecunumanii, which could indicate a lack of ABA signaling despite the up-regulation of receptors. As ABA has been implicated in suppression of GA responses by stabilizing DELLA proteins, preventing their degradation [34], the up-regulation of ABA biosynthesis in P. tecunumanii could be related to GA suppression. 2.6.5. Ethylene 1-aminocyclopropane-1-carboxylic acid (ACC) synthase (ACS), ACC oxidase (ACO) and a large amount of ET response factor (ERF) genes were up-regulated at both time points in P. tecunumanii and at 7-dpi in P. patula (Figure 2, Table S13). More ERFs were up-regulated at 3-dpi in P. tecunumanii than at 7-dpi. There was also a S-adenosyl-methionine (SAM) synthetase up-and an ET insensitive 2 (EIN2) gene down-regulated at 7-dpi in P. patula. The up-regulation of ET biosynthesis genes and ERFs indicates active ET signaling at both timepoints in P. tecunumanii and at 7-dpi in P. patula. However, down-regulation of EIN2 and up-regulation of a predicted reversion-to-ethylene sensitivity 1 (RTE1), a known negative regulator of ET signaling [35], in P. patula could interfere with ET signaling. ET signaling has been implicated in the suppression of ABA biosynthesis as well as negative regulation of ABA signaling [36]. Thus, despite up-regulation of ABA receptors at both time points in P. tecunumanii, up-regulation of ET biosynthesis could indicate suppression of ABA signaling by ET.

Salicylic Acid
There are two main routes for SA production in plants, the isochorismate (IC) and phenylalanine ammonia-lyase (PAL) pathways, both originating from the shikimate pathway product chorismate. In P. tecunumanii, shikimate pathway and PAL genes were up-regulated at both time points, with more genes up-regulated at 7-dpi ( Figure 2, Table S14). Phytoalexin deficient 4 (PAD4) was also up-regulated at both time points, however, enhanced disease susceptibility 1 (EDS1) was only up-regulated at 7-dpi. Conversely, in P. patula, while shikimate pathway and PAL genes were also up-regulated at 7-dpi, EDS1 and the majority of PAD4 DEGs were down-regulated. Isochorismate synthase (ICS), the first step in the IC pathway [37], was also down-regulated at 7-dpi. Furthermore, isochorismatase, an enzyme that diverts isochorismate away from SA synthesis, was up-regulated. PR-1 genes, classic SA response genes [38], were up-regulated at both time points in P. tecunumanii, while in P. patula, they were up-regulated at 7-dpi and down-regulated at 3-dpi (Tables S13 and S14). Thus, while SA signaling appeared to be active at 7-dpi in P. tecunumanii, at 7-dpi in P. patula there were indications of SA biosynthesis and signaling suppression.

Auxin
In P. patula, the only auxin related DEGs at 3-dpi were an up-regulated auxin response factor (ARF) 19 and a down-regulated small auxin up-regulated (SAUR)-like gene ( Figure 2, Table S13). At 7dpi, an auxin biosynthesis gene (YUCCA) was down-regulated in P. patula. As YUCCA is only involved in one auxin biosynthesis pathway, auxin could still be produced via other routes [39]. ARF4, ARF6, ILL5 (IAA-amino acid hydrolase), CAND1 (Cullin-associated and neddylation dissociated 1) and three Aux/IAA auxin repressor genes were also down-regulated at 7-dpi in P. patula, indicating suppression of auxin signaling, despite the up-regulation of putative SAUR genes as well as GH3 (Gretchen Hagen 3 acyl acid amido synthetase family proteins), IAMT1 (IAA carboxyl methyltransferase 1) and ILL6 genes. In P. tecunumanii, no auxin biosynthesis genes were differentially expressed, however, auxin influx carrier genes were up-regulated at both time points and an auxin efflux carrier gene was up-regulated at 3-dpi. Furthermore, Aux/IAA, GH3 and SAUR genes were up-regulated at both time points and there was up-regulation of IAMT, ILL6 and CAND1 genes as well as down-regulation of some Aux/IAA and SAUR genes at 7-dpi. The up-regulation of ILL6 and influx proteins indicate an increase of auxin levels in P. tecunumanii at both time points, however, the efflux carrier indicates lower auxin levels at 3-than 7-dpi, which is reflected in the amount of auxin response DEGs.
A MES1 gene, a methyl esterase capable of hydrolyzing MeSA, MeJA and MeIAA, was up-regulated in both hosts at 7-dpi. There was also up-regulation of MES17, a MeIAA specific methyl esterase, at 7-dpi in P. tecunumanii. Systemic signaling through hormone methyl esters requires demethylation to activate the hormone [40][41][42]. Thus, up-regulation of IAMT1 and MES genes could indicate systemic signaling at 7-dpi in both hosts.
At 7-dpi, cullin genes were up-regulated in P. tecunumanii but down-regulated in P. patula. Cullins are critical structural proteins of Skp-Cullin-F-box (SCF) complexes. This could indicate a decreased ability to activate defence signaling in P. patula as SCF complex catalyzed ubiquitination is an important component of GA, JA and auxin signaling [43][44][45].

Discussion
The disparate F. circinatum resistance phenotypes of P. patula and P. tecunumanii [15] provided a pathosystem to study resistant and susceptible host-pathogen interactions between Pinus spp. and F. circinatum. The susceptibility of P. patula and resistance of P. tecunumanii to F. circinatum challenge was previously confirmed by the significant difference in lesion development between species, as well as mortality of P. patula seedlings, while P. tecunumanii seedlings showed signs of recovery [28]. Significantly higher read mapping to the F. circinatum transcriptome from inoculated relative to mock-inoculated samples, and up-regulation of most genes in the F. circinatum high confidence expressed genes, for each sample set supported the presence of pathogen sequence reads, as expected from the reported infection [28]. This was corroborated by higher mapping from 7-relative to 3-dpi inoculated samples, which is indicative of fungal growth.
Expression of F. circinatum genes indicated compromised ergosterol biosynthesis during infection of the resistant host at 7-dpi, which could increase pathogen susceptibility to PR-1 proteins [14]. All detected F. circinatum orthologs of five ergosterol biosynthesis genes, ERG11, ERG25, ERG6, ERG2 and ERG3, were expressed at lower levels during infection of P. tecunumanii relative to P. patula. Azole group fungicides inhibit fungal growth by inhibiting ERG11, blocking ergosterol biosynthesis [46]. One known mechanism of azole resistance in Candida albicans results from a loss of function mutation of ERG3, without affecting virulence [47]. Conversely, in F. graminearum ERG3 mutation has been associated with decreased virulence [48]. Although many fungal pathogens exhibit resistance to azole fungicides, plant coumarins have also been associated with ERG11 inhibition in C. albicans [49], and the transformation of Arabidopsis and barley with an ERG11 targeting double-stranded RNA, to elicit host-induced gene silencing, resulted in complete immunity to F. graminearum [50]. Additionally, a tomato glycoalkaloid has been associated with suppression of ergosterol biosynthesis in Saccharomyces cerevisiae by inhibiting ERG6 [51] and treatment of C. albicans with the terpenoid farnesol resulted in down-regulation of ERG11, ERG25, ERG6 and ERG3 [52]. Thus, the lower relative expression of these genes during infection of the resistant host could aid in host resistance.
Enriched GO terms indicated a delayed and imprecise response to F. circinatum challenge by P. patula. The overrepresented GO terms for P. tecunumanii DEGs suggested active host defence responses at both time points, while defence related GO terms were only enriched at 7-dpi for P. patula. Additionally, in P. tecunumanii there was enrichment of ubiquitin-mediated proteolysis and cell cycle regulation terms in the up-regulated DEGs and DNA and cellular replication terms in the down-regulated data which were absent for P. patula. Ubiquitin-mediated proteolysis is a critical process for the activation and regulation of GA, JA, SA and auxin signaling pathways [43][44][45]53]. Several lines of evidence suggest the existence of a growth-defence trade-off in plants [54], thus decreasing replication could assist the host in mounting a successful defence. Furthermore, although ET, JA and SA related terms were enriched for the up-regulated 7-dpi P. patula DEGs, similar to P. tecunumanii, there was also enrichment for terms related to oxidative stress, apoptosis and ROS production. In plants, these responses are associated with HR [55], which has been linked to increased susceptibility against necrotrophic pathogens and has been shown to be promoted by necrotrophs to facilitate infection [56].
Phytohormone related DEGs at 3-dpi indicated roles for auxin and ET in P. tecunumanii defence responses. Plant defence against necrotrophic pathogens is usually associated with active signaling by both the ET and JA pathways [57][58][59][60]. However, while no JA biosynthesis genes were up-regulated, putative JAZ genes (repressors of JA responses) and a JA hydroxylase (involved in JA degradation) were up-regulated, indicating suppression of JA signaling. The jasmonate-insensitive 1 (JIN1, a.k.a. AtMYC2) MYC protein has been shown to negatively regulate EIN3 expression, inhibiting the expression of ERF1 [61]. Thus, JA suppression could allow for the activation of a larger repertoire of ET responses. Auxin has also been shown to have an antagonistic effect on JA signaling by stimulating the expression of JAZ proteins [62].
At 7-dpi, P. tecunumanii DEGs suggested the inclusion of JA and SA in host resistance. JA biosynthesis genes were up-regulated, indicative of JA signaling and, while more ET biosynthesis genes were up-regulated, markedly fewer ERFs were up-regulated, reflecting the expected JA/ET antagonism [61]. Additionally, a larger array of auxin response genes were differentially regulated relative to 3-dpi. A synergistic interaction between auxin and JA signaling has been associated with enhanced host resistance to necrotrophic pathogens [45,63]. There was also up-regulation of SA biosynthesis genes and the SA response marker genes PAD4 and EDS1 [64,65], as well as more PR-1 genes at 7-relative to 3-dpi. Therefore, while at 3-dpi P. tecunumanii seemed to induce ET and auxin while suppressing JA signaling, at 7-dpi, transcriptomic responses indicated a complex regulation of host defence using auxin, ET, JA and SA. Although SA signaling is usually classified as antagonistic to both auxin and JA signaling, synergistic interactions also exist [32,45,66]. SA signaling has been shown to induce JA biosynthesis and modulate JA defences in Arabidopsis effector-triggered immunity [67,68]. The early auxin responsive GH3 proteins play an important role in mediating crosstalk between auxin, JA and SA [45]. The core JA signaling component JAR1, which is required for production of the bioactive JA-Ile conjugate, is a GH3 [69], and increased GH3 expression has been shown to simultaneously induce the SA pathway and derepress the auxin pathway [70]. The SA conjugate salicyloyl-aspartate has also been implicated as a signaling molecule to induce systemic resistance [70,71]. Thus, auxin signaling could play an important role in coordinating and integrating phytohormone defence pathways, similar to the central role played in growth and development [45,[72][73][74].
Host responses at 3-dpi in P. patula indicated a lack of defence and increased membrane permeability. Less than 400 genes were differentially expressed at 3-dpi in the susceptible host. A previous study investigating P. patula responses to F. circinatum challenge at 1-dpi identified even fewer DEGs [75]. Although overrepresentation of the GO term, response to chitin, at 3-dpi indicated fungal perception [76,77], there was no enrichment for defence related terms and few phytohormone related or PR family DEGs. Furthermore, there was overrepresentation of the GO term, syncytium formation, in the up-regulated DEGs and the term, xyloglucan:xyloglucosyl transferase activity, in the down-regulated DEGs. Expansins are involved in loosening of the cell wall associated with growth as well as symbiotic interactions [78,79]. Suppression of these proteins has been associated with increased resistance to necrotrophic pathogens [79,80]. Xyloglucan:xyloglucosyl transferases are involved in covalent cross-linking of cell wall polymers, such as xylose, and have been associated with cell wall reinforcement through xyloglucan remodeling [81,82]. Increased levels of cell-wall bound xylose have been associated with resistance to necrotrophic pathogens, while decreased levels have been associated with susceptibility, in Arabidopsis [83]. Thus, the up-regulation of expansins and the down-regulation of xyloglucan:xyloglucosyl transferases could contribute to susceptibility in P. patula by increasing membrane permeability. Combined with the seeming lack of defence responses, this could indicate effector triggered susceptibility at 3-dpi.
Despite the enrichment of defence related GO terms, P. patula DEGs at 7-dpi suggested an impaired phytohormone defence response against F. circinatum. Regarding SA related DEGs, the biosynthesis gene ICS2, as well as the SA response marker genes EDS1 and PAD4 were down-regulated. Thus, unlike P. tecunumanii, the SA defence pathway appeared to be suppressed at 7-dpi. Biosynthesis genes for ET and JA, as well as many ERF and JAZ genes, were up-regulated. Auxin responsive SAUR and GH3 genes were up-regulated, while Aux/IAA and ARF genes were down-regulated. Although this could indicate signaling by these phytohormones, key JA and ET signaling genes, COI1, TPL and EIN2, were down-regulated and a negative regulator of ET signaling, RTE1, was up-regulated, indicating that these signaling pathways could be compromised. Auxin signaling has been shown to have an antagonistic effect on JA signaling by stimulating the expression of JAZ proteins [62,84,85] as well as inducing the expression of ERFs [86][87][88]. However, similar to JA and GA signaling, auxin signaling requires SCF-mediated repressor degradation [44,89]. Down-regulation of CAND1 and cullin genes suggests a decrease in available SCF complexes. Thus, despite the difference in the number of DEGs between time points, indicating a delayed response, immune signaling appeared to be compromised at 7-dpi. Nonetheless, the large number of up-regulated defence-related and PR-genes suggest some form of defence signaling. One of the down-regulated auxin response factors, ARF2, has been implicated in the negative regulation of COI1-independent defence responses in A. thaliana [90]. Similarly, the up-regulation of a rice GH3 has been associated with SA-and JA-independent immunity in rice [80]. Consequently, the differential regulation of defence related genes at 7-dpi could be the result of phytohormone-independent signaling.
In summary, this dual RNA-seq study suggested that ergosterol could be required for F. circinatum virulence in pine and identified phytohormone signaling pathways potentially involved in host resistance and susceptibility. This study purports that the observed reduction of F. circinatum ergosterol biosynthetic gene expression could compromise this pathway, combined with the up-regulation of host PR-1 genes, this could be a key factor contributing to host resistance. Future work to determine the effect of supressing F. circinatum ergosterol biosynthesis on host susceptibility and pathogen response to host metabolite treatment is required to provide further support for the observed responses. Furthermore, DEGs in P. tecunumanii indicated the integration and coordination of auxin, ET, JA and SA mediated defence responses could be required for resistance, while the absence of defence responses at 3-dpi and the down-regulation of phytohormone signaling components at 7-dpi in P. patula suggested pathogen inhibition of host responses. To the authors' knowledge, these results represent the first comprehensive investigation of F. circinatum gene expression during pathogenesis of pine, as well as the first comparison of host responses to F. circinatum between two different pine species. Although an ideal comparison would have been to compare responses between resistant and susceptible genotypes of each species, no susceptible P. tecunumanii LE genotypes are known and even the most tolerant P. patula genotypes are still susceptible. Despite this limitiation, the current approach improves knowledge on the pine-F. circinatum host pathogen interaction, as well as adding to the limited knowledge of defence responses in conifers. While this study investigated how these resistant and susceptible hosts respond to F. circinatum challenge on a transcriptomic level, under the assumption that changes in expression are implicated in defence, it is possible that differences in basal defences between the species would contribute to the resistant or susceptible outcome. Future work is required to investigate the role of physiological differences between these hosts.

Read Data From F. circinatum Inoculation Trial
RNA-sequencing libraries for P. patula and P. tecunumanii, generated from a F. circinatum inoculation trial, were obtained from the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA). In brief, six-month-old seedlings, from open-pollinated families, were inoculated with F. circinatum isolate FCC3579 and mock-inoculated using sterile glycerol [28]. The top 1 cm of shoot tissue was harvested at 3-and 7-days post inoculation (dpi) for three biological replicates per treatment group. A biological replicate consisted of tissue pooled from 16 seedlings. The complete disease progression was reported previously [28]. At 3-and 7-dpi, disease symptoms were not observed on either host species, however, at 14-dpi P. patula showed marked lesions from the point of inoculation. The expected difference in host resistance, represented by the difference in lesion development rates [15,28], was clearly visible by 21-dpi, with pronounced lesion development on inoculated P. patula and only mild discolouration on inoculated P. tecunumanii ( Figure S1). RNA was extracted using the Plant/Fungi RNA Purification Kit (Norgen Biotek Corp., Thorold, ON, Canada) and sent to Novogene (Novogene Corporation Inc, Chula Vista, CA, USA) for strand specific RNA-sequencing on an Illumina HiSeq2500 (Illumina, San Diego, CA, USA), PE125 for 3-dpi samples and PE150 for 7-dpi samples.

Mapping and Gene Expression Analysis
For both host species, the host and pathogen reference transcriptomes were combined ( Figure S2). Read mapping and expression quantification, against the combined references, was performed using Kallisto 0.42.4 [100], with sequence bias correction and 40 bootstrap samples. The count data was imported into R 3.4.2 [101], using tximport 1.4.0 [102], for differential expression (DE) analysis. Fungal genes expressed at equal levels in inoculated and mock-inoculated samples, when normalizing against the full read set, were likely due to endophyte expression. Therefore, a high-confidence pathogen expressed gene set was produced to exclude potential endophyte contamination by performing a DE analysis using the full count data set (including both host and pathogen mapped reads) for each host ( Figure S2). F. circinatum genes significantly up-regulated in inoculated relative to mock-inoculated samples were considered high confidence expressed genes for each time point. Significantly down-regulated F. circinatum genes, in inoculated relative to mock-inoculated samples, could represent potential endophyte genes and were discarded for downstream analysis. For pathogen DE analysis, F. circinatum count data from both host species was combined. Pathogen expression data was removed from the count data for host DE analysis. Transcripts with less than 20 reads from at least three read libraries were classified as low expression transcripts and filtered out. DE testing was performed with DESeq2 1.18.1 [103] using a Wald test with Benjamini and Hochberg (BH) false discovery rate (FDR) correction (p < 0.05, Abs|Log2(Fold-change)| ≥ 0.5). Biological pathways related to the differentially expressed genes (DEGs) were investigated using Mercator annotations with MapMan v3.5.1R2 [104], as well as KEGG orthology with KEGG mapper reconstruct pathway tool [105][106][107].

Gene Ontology Enrichment Analysis
DEGs for each comparison were divided into up-and down-regulated subsets (inoculated versus mock-inoculated for host transcripts; P. tecunumanii vs. P. patula for F. circinatum transcripts, Figure S2). Significant enrichment of GO terms (BH FDR, p < 0.10), relative to the transcriptome annotation for each species, was determined for each high confidence pathogen expressed gene set as well as the host and pathogen DEG subsets using GOSeq 1.28.0 [108].

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.