Differential gene expression in banana roots in response to Fusarium wilt

Abstract

Banana wilt disease, caused by Fusarium oxysporum f. sp. cubense tropical race 4 (Foc TR4), is a devastating disease in the banana industry. The molecular mechanism underlying the wilt resistance of AAA-type banana cultivars has not been fully characterized. Our objective was to analyse the genes differentially expressed in roots during the time of flower bud differentiation in plants grown in fields infested with Foc TR4. We compared the Foc TR4-tolerant ‘Formosana’ banana and the susceptible ‘Brazil’ banana. cDNA libraries were constructed from RNA isolated from the roots of both cultivars, and generated approximately 40 million and 35 million high-quality reads, respectively. We mapped 28 353 810 ‘Formosana’ reads and 26 917 421 ‘Brazil’ reads to the banana genome. A search of the NR (non-redundant) database resulted in the annotation of 34 408 genes. Additionally, 107 genes were differentially expressed between the two cultivars. Among the 48 differentially expressed genes with known functions, 41 were expressed more highly in ‘Formosana’ than in ‘Brazil’. These genes were divided into several categories, including specific resistance mechanism-related enzyme genes and phytoalexin synthesis-related genes. The expression levels of eight genes were validated by a quantitative real-time polymerase chain reaction. Our results indicate that ‘Formosana’ may enhance its tolerance to Foc TR4 by increasing the expression of defence-related genes (possibly constitutively). This cultivar also exhibits up-regulated expression of non-specific stress-related genes, which may enhance its overall disease tolerance. This study provides important transcript-level details potentially useful for clarifying the molecular mechanism underlying the disease resistance of Cavendish banana cultivars.

Résumé

Le flétrissement de la banane, causé par la race tropicale 4 de Fusarium oxysporum f. sp. cubense (Foc TR4), est une maladie dévastatrice qui frappe l’industrie de la banane. Le mécanisme moléculaire qui est à la base de la résistance au flétrissement des cultivars de banane de type AAA n’a pas été entièrement défini. Notre objectif était d’analyser les gènes différentiellement exprimés dans les racines durant la période de différenciation des boutons floraux chez les plants cultivés dans des champs infestés par Foc TR4. Nous avons comparé le cultivar ‘Formosana’, tolérant à Foc TR4, et le cultivar ‘Brazil’ réceptif à son égard. Des bibliothèques d’ADNc ont été construites à partir d’ARN obtenu des racines des deux cultivars, et elles ont engendré environ 40 millions et 35 millions de lectures, réciproquement. En ce qui a trait au génome de la banane, nous avons cartographié 28 353 810 lectures de ‘Formosana’ et 26 917 421 lectures de ‘Brazil’. Une recherche effectuée dans la base de données non redondante a entraîné l’annotation de 34 408 gènes. En outre, 107 gènes ont été différentiellement exprimés dans les deux cultivars. Parmi les 48 gènes différentiellement exprimés possédant des fonctions connues, 41 étaient plus fortement exprimés dans ‘Formosana’ que dans ‘Brazil’. Ces gènes ont été divisés en plusieurs catégories, y compris en gènes d’enzymes associés aux mécanismes de résistance spécifique et en gènes associés à la synthèse de phytoalexines. Les niveaux d’expression de huit gènes ont été validés par un test quantitatif de PCR en temps réel. Nos résultats indiquent que ‘Formosana’ peut accroître sa tolérance à Foc TR4 en augmentant l’expression des gènes associés aux mécanismes de défense (probablement de manière constitutive). Ce cultivar affiche aussi une expression régulée positivement de gènes non spécifiques associés au stress, susceptible d’améliorer sa tolérance globale à la maladie. Cette étude fournit d’importants détails quant au degré de transcription, détails possiblement utiles pour préciser le mécanisme moléculaire qui est à la base de la résistance aux maladies chez les cultivars de banane Cavendish.

Keywords:

• banana
• differentially expressed genes
• flower bud differentiation period
• Fusarium oxysporum f. sp. cubense
• transcriptome sequencing

Mots clés:

• Banane
• gènes différentiellement exprimés
• Fusarium oxysporum f. sp. cubense
• période de différenciation des boutons floraux
• séquençage du transcriptome

Introduction

Fusarium wilt of banana (Musa spp.), caused by Fusarium oxysporum Schlechtend.:Fr. f. sp. cubense (E.F. Smith) Snyder & Hansen (Foc), was first reported from Panama as early as 1890 and is commonly known as Panama disease. The pathogenic fungus invades plants mainly through the radicle and damaged corms, after which it migrates through the vascular bundles to the pseudostem and leaves, ultimately leading to plant death (Guo & Yang 2009). This disease occurs worldwide in the major banana-producing areas, including Panama, the Philippines, Taiwan, and in southern China (Lin 2010). Thus, it is a serious threat to a sustainable global banana industry. Some non-AAA varieties, such as plantains (AAB), are highly resistant to Fusarium wilt disease, but because of their relatively small-scale production, they are not as economically important as the AAA varieties. The tribute banana (AA) is also somewhat resistant to wilt disease (Smith et al. 1993). In the current large-scale commercially cultivated AAA type bananas, Taiwan reported a few resistant somaclonal-derived varieties such as GCTCV-53, GCTCV-119, GCTCV-44–1, GCTCV-215 and GCTCV-218 (Hwang & Ko 2004). Furthermore, three cultivars from Guangdong, China are considered tolerant to wilt disease, namely ‘Yueyoukang 1ʹ, ‘Nongke No. 1ʹ and ‘Kangku No. 5ʹ. Moreover, there are ongoing efforts to improve the resistance of cultivated banana species to Foc through the introduction of transgenic resistance genes (Pei et al. 2005; Mahdavi et al. 2012).

To develop an effective screening procedure for the breeding of wilt-resistant banana varieties, the mechanism underlying wilt disease resistance should be characterized. The organ structure, physiology, and biochemistry of infected plants undergo a series of changes in response to the invading pathogen. Additionally, several studies have been conducted regarding the mechanism regulating wilt disease resistance in banana. When studying the invasion of roots by Foc, Beckman et al. (1961) observed that in comparisons between susceptible and moderately or highly wilt-resistant banana varieties, the space between vessel perforation plates was smaller in the disease-resistant varieties, as was the vessel lumen diameter. Thus, in the wilt-resistant banana varieties, the migration of fungal spores through the vascular bundles is somewhat inhibited, which minimizes the pathogen infection rate. Zuo et al. (2010) reported that the organic components of the root exudates from highly disease-resistant banana varieties were lethal to fungal spores. The organic components in the root exudates of moderately disease-resistant varieties significantly inhibited spore germination and mycelial growth, while the root exudates of susceptible varieties induced spore germination and mycelial growth. De Ascensao & Dubery (2000) confirmed that peroxidases, polyphenoloxidases, phenylalanine ammonia lyases and other related enzymes were involved in the lignification of the host cell wall in banana roots treated with an elicitor from the pathogen cell wall. Thangavelu et al. (2003) infected disease-resistant Cavendish banana plants with Foc and Pseudomonas fluorescens (pf10 strain). They observed that the abundance of leaf phenolic compounds, peroxidases, phenylalanine ammonia lyases, chitinase and β-1,3-glucanase significantly increased within 3–6 days after inoculations. Similar changes were not detected in susceptible banana plants. Li et al. (2010) reported that phenylalanine ammonia lyase, β-1,3-glucanase and exo-chitinase activities are important biochemical indicators for identifying wilt-resistant banana varieties. Moreover, using suppression subtractive hybridization, Van Den Berg et al. (2007) compared the expression of defence-related genes between GCTCV-218 and Cavendish ‘Williams’ (i.e. susceptible) plants infected by Foc TR4. Several defence-associated genes were identified, including those encoding PR-1, pectin acetylesterase, xylanase inhibitor, peroxidase, catalase 2 and metallothionein. They also revealed that phenolics were more abundant in GCTCV-218 plants than in ‘Williams’ plants. This suggested that constitutively produced defence compounds are present in higher concentrations in GCTCV-218 plants, even prior to pathogen invasion.

High-throughput sequencing techniques represent new and effective ways to study the wilt disease resistance mechanism in bananas (Schmutz et al. 2010; Wang X et al. 2010; Wang Z et al. 2010; Liu et al. 2011; Tirumalaraju et al. 2011). These techniques have been applied to study the transcriptome of banana plants infected by Foc TR4 (Li et al. 2012, 2013; Wang et al. 2012; Bai et al. 2013). Li et al. (2012) and Bai et al. (2013) used tissue-cultured plantlets to analyse the transcriptomes of wilt-resistant Cavendish banana cultivars (‘Nongke No. 1ʹ and ‘Yueyoukang No. 1ʹ) and a wilt-susceptible wild-type ‘Brazil’ cultivar 1–10 days after the roots were infected with Foc TR4. However, characterizing the transcript-level differences between banana cultivars infested with Foc TR4 under field conditions may be more useful for elucidating the disease resistance mechanism. Additionally, the leaves of infested banana plants start to turn yellow and wilt during the flower bud differentiation stage (Li et al. 2007). In Hainan, China, the predominant ‘Brazil’ banana cultivar is susceptible to Foc TR4, whereas ‘Formosana’ banana, which was derived from the Giant Cavendish somaclonal variant GCTCV-218, is tolerant to Foc TR4. Thus, ‘Formosana’ banana plants are increasingly being grown in Hainan. We compared the root transcriptomes of ‘Formosana’ and ‘Brazil’ plantlets infested with Foc TR4 under field conditions. We analysed the flower bud differentiation period, and identified the differentially expressed genes. The gene expression level differences between ‘Formosana’ and ‘Brazil’ banana plants may be due to multiple factors related to variety characteristics.

Materials and methods

Plant materials, treatments and RNA extraction

Experiments were conducted using ‘Formosana’ (Musa acuminate L. AAA Cavendish cv. ‘Formosana’) and ‘Brazil’ (Musa acuminate L. AAA Cavendish cv. ‘Brazil’). Cavendish cultivars (AAA genome) are highly susceptible to Foc TR4. ‘Formosana’ banana was derived from the Fusarium wilt-tolerant Cavendish selection GCTCV-218 screened at the Taiwan Banana Research Institute. In 2004, it was introduced into Hainan, China by the transfer centre of the International Network for the Improvement of Banana and Plantain in Belgium. In 2012, the ‘Formosana’ banana was approved by the Crop Variety Certification Committee of Hainan province (variety certification number: Qiongrenxiangjiao2012003). During field trials in Hainan, ‘Formosana’ banana plants were observed to be moderately resistant to Fusarium wilt, which is consistent with similar field tests in Guangdong (Huang et al. 2005). ‘Formosana’ bananas are superior to ‘Brazil’ bananas regarding horticultural characteristics (e.g. a more robust pseudostem, stronger petioles, thicker leaves, better hand formation, more uniform hand size and higher yield). However, ‘Formosana’ plants grow more slowly than the ‘Brazil’ plants, requiring an extra month to complete the growth period (Hwang & Ko 2004). The two banana cultivars used in this study were obtained from the Environment and Plant Protection Institute, Chinese Academy of Tropical Agricultural Sciences (Hainan, China), and planted at the Hongguang Farm in Fushan (Hainan, China).

At the experimental field site, ‘Formosana’ and ‘Brazil’ plants were planted next to each other in soil (sticky red loam) infested with Foc TR4. Within each row, individual plants were separated by 2 m. Adjacent rows were also separated by 2 m. The Foc TR4 strain was cultured on Wheat Bran agar medium at 28°C for 10 days. The spores were scraped into sterile water and the spore suspension was adjusted with sterile water to a concentration of about 1.5 × 10⁶ condia mL⁻¹. To confirm the viable spore concentration, the suspension was diluted with distilled water to concentrations of 10^–1, 10^–2, 10^–3, 10^–4, 10^–5 and plated on potato dextrose agar for 30 h at 28°C, then the colonies were counted. Before planting the banana, Foc TR4 spore suspension was added to the soil (5 L m⁻²). To ensure Foc TR4 was maintained at a suitable density in the soil, spore suspension was added three times, once every 20 days. Before banana planting, the soil was sampled and then PCR used to identify the Foc TR4. Root samples were collected when the plants reached the flower bud differentiation stage (i.e. after approximately 4.5 months), which coincided with when the ‘Brazil’ plants began to exhibit widespread wilt symptoms. The roots (150 g; including the main roots, lateral roots and root hairs with consistent morphology) collected from five plants for each banana variety were combined, immediately frozen in liquid nitrogen, and stored at −80°C until analysed.

Total RNA was extracted from each sample using the TRIzol^® reagent (Invitrogen, CA, USA) following the manufacturer’s recommended procedure. Additionally, total DNA was extracted from 5 g root tissue for each sampled banana plant using a published CTAB method (Murray & Thompson 1980). The purified DNA was used to assess the infection of banana plants by Foc TR4 according to a polymerase chain reaction (PCR)-based method previously described by Lin et al. (2009).

To compare the effects of Foc TR4 on the roots of ‘Formosana’ and ‘Brazil’, a few hydroponically grown plantlets of each banana variety were treated with green fluorescent protein (GFP)-labelled Foc TR4 (10⁶ spore mL⁻¹) for 2 h. Uninoculated ‘Brazil’ plantlets were used as control samples. The plantlets were then returned to the hydroponic culture for 48 h. The banana roots were collected, and the fluorescence of freehand sections was analysed using a FV10-ASW laser scanning confocal microscope (Olympus, Tokyo, Japan) after the sections were exposed to a 488-nm laser line to excite the GFP. The recombinant GFP-labelled Foc TR4 strain used in this study was generated according to a polyethylene glycol-mediated transformation method (Zhang et al. 2007). The pCT74 plasmid used to transform Foc TR4 cells included the strong constitutive ToxA promoter from Pyrenophora tritici-repentis to regulate GFP expression. The recombinant fungal strain was stored at the Environment and Plant Protection Institute, Chinese Academy of Tropical Agricultural Sciences in Haikou, China.

cDNA library construction and sequencing

High-throughput Illumina sequencing using the HiSeq2000 platform was completed by Shanghai Majorbio Bio-Pharm Technology Co., Ltd according to the manufacturer’s recommended procedure. The TruSeq™ RNA Sample Prep Kit (Illumina) was used to construct libraries. The whole process involved the preparation of the library and flow cell as well as the high-throughput sequencing. Details regarding the procedures are as follows: (1) 5 μg total RNA was used as the template for constructing libraries; (2) mRNA was isolated using magnetic beads (TruSeq™ RNA Sample Prep Kit); (3) double-stranded cDNA was synthesized, after which a single A base was added to the 3′ end and the index joints were linked (TruSeq™ RNA Sample Prep Kit); (4) the library was enriched and a PCR was completed using 15 cycles; (5) the amplification products were analysed by 2% agarose gel electrophoresis (Certified Low Range Ultra Agarose) and the target product was recovered; (6) the purified amplification product was quantified using TBS380 (Pico green), mixed, and loaded onto the machine; (7) clusters were generated with a bridge PCR amplification on cBot; (8) Paired-end 2 × 100-bp sequencing was completed using the HiSeq2000 sequencing platform.

Data filtering, mapping and gene annotation

High-quality reads were obtained from the original sequencing data after the following quality control procedures using SeqPrep (https://github.com/jstjohn/SeqPrep) and Sickle (https://github.com/najoshi/sickle): (1) removal of low-quality sequencing reads (mass < 25); (2) removal of adaptor sequences; (3) removal of reads with several Ns; (4) removal of short sequences (< 25 nucleotides) after completing steps 1–4. The high-quality reads were mapped using TopHat (http://tophat.cbcb.umd.edu/) and annotated based on the banana genome sequence (http://www.ncbi.nlm.nih.gov/genome/10976). The mapped reads were functionally annotated with gene ontology (GO) terms using the Blast2GO program (http://www.blast2go.com/b2ghome). Phylogenetically diverse domain families were analysed based on Cluster of Orthologous Groups (COG) annotations completed using the blastx 2.2.24+ program and STRING 9.0 database (http://string-db.org/). The Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway annotations were completed using the blastx/blastp 2.2.24+ program and the KEGG genes database (http://www.genome.jp/kegg/genes.html). Additionally, a pathway enrichment analysis for all differentially expressed genes was performed using KOBAS (http://kobas.cbi.pku.edu.cn/home.do).

Analysis of differentially expressed genes

Differentially expressed genes were identified by comparing the abundance of the same transcript in the roots of ‘Formosana’ and ‘Brazil’. The banana genome sequence (D’Hont et al. 2012) was used as a reference. The number of clean reads mapped to genomic regions was used to calculate gene expression levels. The fragments per kilobase of exon per million mapped reads (FPKM) value of each transcript was used to represent the gene expression level. Furthermore, the significance of the expression level differences was assessed to identify the differentially expressed genes. A visual analysis was also completed. The expression patterns of differentially expressed genes (false discovery rate <0.05) were analysed using a clustering algorithm [hcluster (complete algorithm)]. The Spearman and Pearson correlations were based on the distance between samples and the distance between genes, respectively. We used the following programs: TopHat (http://tophat.cbcb.umd.edu/), Cuffdiff (http://cufflinks.cbcb.umd.edu/) and cummerbund (for visualization) (http://compbio.mit.edu/cummeRbund/).

qRT-PCR validation

To validate the Illumina sequencing results, eight genes related to a specific resistance mechanism were analysed by quantitative real-time (qRT)-PCR. Primers (Table 1) were designed using Primer 5.0, and the banana actin gene (ACT1; GenBank number AF285176) was used as the internal reference gene. The total RNA extracted for sequencing was reverse transcribed to cDNA using the Reverse Transcription System (Promega, USA). The cDNA samples were used as the templates for qRT-PCR experiments, which were conducted with a QuantStudio™ 6 Flex Real-Time PCR System (Applied Biosystems, USA) and the SYBR Green Master ROX (Takara, Japan) following the manufacturer’s recommended procedure. The qRT-PCR samples consisted of 0.8 µL (each) 10 µM forward and reverse primers, 10 µL 2× SYBR Premix Ex Taq II (Takara, Japan), 0.4 µL 50× ROX Reference Dye II, 2 µL cDNA products, and PCR-grade water for a final volume of 20 µL. The PCR conditions were as follows: 95°C for 30 s, 40 cycles of 95°C for 5 s, 56°C for 30 s and 72°C for 30 s. All assays were performed in triplicate. The relative gene expression levels were normalized against the ACT1 levels according to the 2^−ΔΔCt method. The relative expression levels were used to validate the RNA-sequencing data.

Table 1. Details regarding the primers used in this study.

Primer name	Primer sequence (5ʹ – 3ʹ)	Gene
Actin-F	CGAGGCTCAATCAAAGA	Actin (endogenous reference gene)
Actin-R	ACCAGCAAGGTCCAAAC
Cellulose-F	GCTCAACCAGGTGCTACG	Cellulose synthase catalytic subunit 12
Cellulose-R	ACTTGTCGGTGAGGAGGC
Glucosidase-F	AGGCTACGGCAAACTCAG	Periplasmic beta-glucosidase precursor, putative
Glucosidase-R	GACAGGGCCACAAAGAAA
peroxidase N-F	GCTGCTGCTGTTGCTGTAGA	PREDICTED: peroxidase N
peroxidase N-R	TCAAGAACCTGCTGAGCCA
laccase-4-F	GGTTATCCGCCCTGAATC	PREDICTED: laccase-4
laccase-4-R	AACGGCACGAGGCTCTAT
Chitinase-F	GGATGGGCTGAAGGTTGA	Putative chitinase
Chitinase-R	CAGGTCGAGGTAATACTGGTAA
T-phylloplanin-F	CATTTGCGATGCTACTGA	Similar to T-phylloplanin
T-phylloplanin-R	GCTACGGTGAACTTGGAC
TPA-F	ATCTACAACGAGACGAGCAT	TPA: class III peroxidase 136 precursor
TPA-R	CGAAGTCGCTGAAGAACC
Cla47-F	TTGGGCAGAAGAGCGAGAT	Glutathione-S-transferase Cla47
Cla47-R	AGGGAGGATGGACGGGTT

Results and discussion

Comparison of the wilt disease resistance traits between 'Formosana’ and ‘Brazil’

‘Formosana’ banana plants are moderately resistant to Fusarium wilt, but can be susceptible in certain situations. Thus, it is considered to be a Fusarium wilt-tolerant variety. After being grown for 4.5 months in soil infested with Foc TR4, ‘Brazil’ plants died because of the disease. In contrast, 90% of’Formosana’ plants continued to grow normally, with no obvious wilt disease symptoms (Fig. 1). The differences between the roots of 'Formosana' and 'Brazil' grown hydroponically and infected by Foc TR4 are presented in Fig. 2. The 'Brazil' roots had considerably more fungal hyphae attached (Fig. 2c) than the 'Formosana’ roots (Fig. 2b). The PCR results (Fig. 3) confirmed that the roots of 'Formosana’ and 'Brazil' were infected by Foc TR4.

Fig. 1 (Colour online) Comparison of disease response of ‘Brazil’ and ‘Formosana’ banana. (Left) ‘Formosana’ banana; (right): ‘Brazil’ banana.

Fig. 2 (Colour online) Comparison between the ‘Formosana’ and ‘Brazil’ roots infected by Foc TR4. A, Uninoculated ‘Brazil’ plantlet roots. B, ‘Formosana’ roots. C, ‘Brazil’ roots. Scale bar = 20 µm.

Fig. 3 Confirmation of the infection of banana plantlets by Foc TR4 based on a PCR assay. M: DL2000 marker; K+: Foc TR4 genomic DNA; K−: DNA extracted from the roots of sterilized ‘Brazil’ plantlets; Lanes 1–5: DNA extracted from the roots of five infected ‘Brazil’ plantlets; Lanes 6–10: DNA extracted from the roots of five infected ‘Formosana’ plantlets; K1: H₂O.

Transcriptome sequencing, mapping and annotation

We obtained 2 × 24 007 117 and 2 × 21 815 650 reads for 'Formosana’ and ‘Brazil’, respectively. A quality control step resulted in 2 × 20 083 619 and 2 × 17 773 064 high-quality clean data for 'Formosana’ and ‘Brazil’, respectively (Table 2). These included 28 353 810 ‘Formosana’ reads (70.6%) and 26 917 421 ‘Brazil’ reads (75.7%) that were mapped to the banana genome. A subsequent bioinformatics analysis identified a total of 37 639 transcripts, including 34 408 isoforms that were highly similar to homologous sequences in the NR database. The banana genome includes 36 542 protein-coding genes. The BLAST searches completed in this study identified 34 408 unique protein accessions, indicating that the Illumina paired-end sequencing results represented a substantial proportion of the known genes from banana. Our sequencing data were deposited in the NCBI Sequence Read Archive (accession number SRA108989; http://www.ncbi.nlm.nih.gov/Traces/sra).

Table 2. Details regarding the raw and high-quality reads in bananas ‘Formosana’ and ‘Brazil’.

Samples	The total number of sequences (strip)	The number of bases (bp)	The number of high-quality sequences (strip)	The number of high-quality bases (bp)
‘Formosana’	24 007 117*2	4,849,437,634	20 083 619 *2	3 775 720 372
‘Brazil’	21  815 650*2	4,406,761,300	17 773 064*2	3 341 336 032

Note: Each read consisted of 101 bp.

A total of 25 029 isoforms were assigned at least one GO term for the biological process, molecular function or cellular component categories (Fig. 4). Of these GO terms, ‘cell’ (17 235 sequences, 68.86%), ‘cell part’ (17 235 sequences, 68.86%), and ‘organelle’ (14 333 sequences, 57.27%) were the dominant cellular component terms, while ‘cellular process’ (13 934 sequences, 55.67%), ‘metabolic process’ (13 169 sequences, 52.61%), and ‘response to stimulus’ (5670 sequences, 22.65%) were the most common biological process terms. Additionally, ‘binding’ (12 355 sequences, 49.34%), ‘catalytic activity’ (11 275 sequences, 45.04%), and ‘transporter activity’ (1404 sequences, 5.61%) were the most abundant molecular function terms.

Fig. 4 (Colour online) Gene ontology classifications. The horizontal axis presents the GO term sub-categories, the left vertical axis corresponds to the proportion of genes included in the sub-categories, and the right vertical axis represents the number of genes in the sub-categories.

To further evaluate the completeness of the transcriptomes and predict possible functions, all transcripts were compared against the sequences in the COG database. The results revealed 9461 homologous sequences, which were assigned to appropriate COG clusters. These COG classifications included 23 functional categories (Fig. 5). The five largest categories were ‘general function prediction only’ (30.60%), ‘transcription’ (19.46%), ‘signal transduction mechanisms’ (16.39%), ‘replication, recombination and repair’ (15.60%), and ‘posttranslational modification, protein turnover and chaperones’ (10.33%).

Fig. 5 (Colour online) Cluster of orthologous group classifications. The horizontal axis presents the functional categories, while the vertical axis corresponds to the number of unigenes in the categories.

Genes differentially expressed between 'Formosana’ and ‘Brazil’

There were no significant differences in the distribution of varying FPKM values in the sequenced samples for ‘Formosana’ and ‘Brazil’. There were also no significant differences in the transcript levels of more than 99% of the genes in ‘Formosana’ and ‘Brazil’. Because the two banana cultivars belong to the same genotype (AAA), the two root transcriptomes are comparable. Although the root transcriptomes of ‘Formosana’ and ‘Brazil’ were highly similar, a few differentially expressed genes were identified. Our data revealed that 107 genes were differentially expressed in the roots of the two cultivars (Table S1), with the expression levels of 81 and 26 genes being higher and lower in ‘Formosana’ than in ‘Brazil’, respectively. These genes included 59 that encode unknown proteins. Of the remaining 48 differentially expressed genes, the expression levels of 41 and seven genes were up- and down-regulated in ‘Formosana’. These 48 genes were divided into the following five categories based on their functions: (1) 13 genes related to a specific resistance mechanism (12 were up-regulated in ‘Formosana’); (2) seven phytoalexin biosynthesis-related genes (all were up-regulated in ‘Formosana’); (3) five stress resistance-related genes (non-specific for biotic and abiotic stresses) (three and two were up- and down-regulated in ‘Formosana’, respectively); (4) 16 genes related to nutrition and growth (15 were up-regulated in ‘Formosana’); (5) seven genes related to other functions (four were up-regulated in ‘Formosana’).

Enriched KEGG pathways for the differentially expressed genes

The identified differentially expressed genes were mainly defence genes, with a few stress resistance genes. These genes were related to several pathways, including those related to flavonoid biosynthesis, flavone and flavonol biosynthesis, alpha-linolenic acid metabolism, starch and sucrose metabolism, and phenylpropanoid biosynthesis (Table 3). We detected two type-III polyketide synthase genes, two flavonoid 3′,5′-hydroxylase genes, and one dihydroflavonol 4-reductase gene associated with flavonoid biosynthesis pathways. There were also two flavonoid 3′,5′-hydroxylase genes and one caffeic acid O-methyltransferase gene related to flavone and flavonol biosynthesis pathways. We identified three genes involved in the alpha-linolenic acid metabolic pathway. Additionally, four genes were related to the starch and sucrose metabolic pathway, while another three genes, including a peroxidase N gene, were associated with the phenylpropanoid biosynthesis pathway.

Table 3. The most enriched pathways for the differentially expressed genes in banana.

Term	Id	Sample number	Background number	P-Value	Corrected P-Value	Genes
Flavonoid biosynthesis	ko00941	5	67	8.18E-07	3.19E-05	CUFF.21668.1|CUFF.29055.1|CUFF.26496.1|CUFF.4245.1|CUFF.16156.1
Flavone and flavonol biosynthesis	ko00944	3	29	5.76E-05	0.001123	CUFF.4241.1|CUFF.29055.1|CUFF.26496.1
alpha-Linolenic acid metabolism	ko00592	3	63	0.000588	0.007648	CUFF.26760.1|CUFF.26768.1|CUFF.26742.1
Starch and sucrose metabolism	ko00500	4	243	0.003623	0.035325	CUFF.24142.1|CUFF.4717.1|CUFF.5013.1|CUFF.2143.1
Phenylpropanoid biosynthesis	ko00940	3	144	0.006235	0.048633	CUFF.24142.1|CUFF.1865.1|CUFF.2143.1

Note: Data were analysed using a hypergeometric test and false discovery rate correction method (Benjamini & Hochberg 1995).

Previous studies have examined the transcriptome sequences of banana plants infected by Foc TR4 (Li et al. 2012, 2013; Wang et al. 2012; Bai et al. 2013). Similar to the findings of these earlier studies, we observed the up-regulated expression of defence genes, cell wall modification genes, and phytoalexin biosynthesis genes. However, there were some inconsistencies between our results and those reported in the earlier related studies (Table 4). For example, we did not detect any pathogen-associated molecular pattern (PAMP) recognition genes, signalling/regulatory genes, and transcription factors among the genes differentially expressed between 'Formosana’ and ‘Brazil’. These differences were probably because the previous studies focused on gene expression levels during the early stages of the Foc TR4 infection (i.e. within the first 10 days). In contrast, we examined the plants about 4.5 months after they were initially infected (i.e. adaptation stage). Thus, the expression levels of genes related to signal transductions and transcriptional regulations of early plant defence responses (e.g. synthesis of signalling molecules) were not up-regulated. Instead, many genes associated with terminal plant defence responses were being highly expressed. Our data suggest that the mechanism underlying the tolerance of ‘Formosana’ to Foc TR4 during the flower bud differentiation stage is a little different from the corresponding mechanism of banana seedlings during the early stages of infection. Our results also imply that the up-regulated expression of defence-related genes in ‘Formosana’ may be constitutive, and a distinguishing characteristic of ‘Formosana’ banana plants.

Table 4. Previous transcriptomics studies involving banana plants infected by Foc TR4.

	Banana variety	Plant growth stage	Pathogen	Sampling dates (after Foc infection)	Enriched pathways of differentially expressed genes
This study	‘Brazil’ banana and tolerant ‘Formosana’ banana	4.5 months	Foc TR4	4.5 months	Flavonoid biosynthesis, flavone and flavonol biosynthesis, alpha-linolenic acid metabolism, starch and sucrose metabolism, and phenylpropanoid biosynthesis
Wang et al. (2012)	‘Brazil’ banana	60 d	Foc TR4	0, 2, 4, and 6 d	Phenylalanine metabolism, Phenylpropanoid biosynthesis, Drug metabolism – cytochrome P450, alpha-linolenic acid metabolism, Amino sugar and nucleotide sugar metabolism
Li et al. (2013)	‘Brazil’ banana	50 d	Foc1 and Foc TR4	3, 27, and 51 h	PR proteins, Synthesis of phytoalexins and phenolpropanoids, Cell wall modification, Ethylene biosynthesis and signalling pathway
Li et al. (2012)	‘Brazil’ banana and resistant ‘Nongke No 1’	Plants with four or five leaves	Foc TR4	48 and 96 h	Perception of PAMPs by PRRs, Plant hormones biosynthesis and signalling, Transcription factors, Cell wall modification, Flavonoid biosynthesis, Programmed cell death, PR proteins
Bai et al. (2013)	‘Brazil’ banana and resistant ‘Yueyoukang 1'	8 weeks (plants with five leaves)	Foc TR4	0.5, 1, 3, 5, and 10 d	PR proteins, Transcription factors, Cell wall modification, Phenylpropanoid biosynthesis, Plant hormone signal transduction

In our study, the genes related to the flavonoid biosynthesis and flavone and flavonol biosynthesis pathways were more highly expressed in ‘Formosana’. These pathways influence the biosynthesis of phytoalexins, which are important for plant defences against pathogenic fungi and other microbes. This suggests that the up-regulated expression of these defence-related genes may contribute to the Foc TR4-tolerance of ‘Formosana’. However, the identified genes associated with the alpha-linolenic acid metabolic pathways were more highly expressed in ‘Brazil’. Alpha-linolenic acid affects other metabolic pathways related to the production of oxilipins and waxes associated with plant development and stress responses (Wang et al. 2004; Li et al. 2011). Additionally, the metabolism of alpha-linolenic acid leads to the production of jasmonic acid (JA), which contributes to pathogen resistance. A previous study concluded that biotic and abiotic stresses, such as pathogen infection, wounding and insect feeding, can induce JA biosynthesis by directly activating the expression of genes encoding specific biosynthetic enzymes (Creelman & Mullet 1995). Li et al. (2012) reported that the expression levels of genes related to JA biosynthesis and signalling were higher in resistant banana plants than in susceptible banana plants. However, our findings suggest that the expression levels of genes involved in alpha-linolenic acid metabolic pathways were lower in Fusarium wilt-tolerant banana plants than in susceptible banana plants, which is consistent with the results of a previous study by Bai (2013). Thus, the differences in the expression of JA biosynthesis genes may reflect the differences in the background expression levels between Fusarium wilt-tolerant and -susceptible banana plants.

Validation of gene expression patterns by qRT-PCR analysis

The qRT-PCR results for the eight selected candidate genes were generally consistent with the Illumina sequencing results (Table 5), suggesting the Illumina sequencing data were valid. The inconsistencies may be due to the differences in the algorithms and sensitivity of the two techniques.

Table 5. Relative expression levels of specific resistance mechanism-related genes according to qRT-PCR in banana.

		log2Fold change (A/B)
Gene ID	Description	qRT-PCR	RNA-Seq
CUFF.12751	Cellulose synthase catalytic subunit 12	4.62702	2.69907
CUFF.16474	Periplasmic beta-glucosidase precursor, putative	4.84333	3.56937
CUFF.1865	PREDICTED: peroxidase N	10.71611	3.63376
CUFF.29797	PREDICTED: laccase-4	5.76765	3.16864
CUFF.7118	Putative chitinase	7.60023	2.79053
CUFF.5976	Similar to T-phylloplanin	7.78922	3.87321
CUFF.8570	TPA: class III peroxidase 136 precursor	4.65535	3.31414
CUFF.13852	Glutathione-S-transferase Cla47	1.24489	3.20794

A: ‘Formosana’; B: ‘Brazil’.

Analysis of the molecular mechanism underlying the tolerance of ‘Formosana’ to Foc TR4

Because ‘Formosana’ and ‘Brazil’ were grown in the same field infested with Foc TR4, the differences in gene expression levels may be associated with phenotypic differences, including tolerance to disease. The genes involved in specific resistance mechanisms and phytoalexin biosynthesis probably contributed to the differences in disease tolerance between the ‘Formosana’ and ‘Brazil’ banana plants. The differentially expressed resistance mechanism-related genes are listed in Table 6. Except for the gene encoding snakin-1, all genes were more highly expressed in ‘Formosana’ than in ‘Brazil’, suggesting that these up-regulated genes are possibly involved in the tolerance of ‘Formosana’ to Foc TR4. The rationale for this is as follows. First, some of the genes are involved in the removal of excessive reactive oxygen species (ROS). Plants exposed to stresses, such as extreme temperatures and pathogen invasion, produce an excess of ROS, which is toxic to plants (Neill et al. 2002). Peroxidases help maintain ROS levels in plant tissues at acceptable levels. The glutathione-S-transferase Cla47 has antioxidant and detoxifying functions and can protect cells against oxidative damage (Reinemer et al. 1992). Second, the proteins encoded by some of the genes enhance the structural resistance of plant cells. For example, lignin is a plant cell wall polymer, and increasing the lignin content of a plant host is one of the most effective means of conferring resistance to pathogens. Peroxidases, which are abundant in plant cell walls, are widely believed to catalyse the last enzymatic step in lignin biosynthesis (Higuchi 1985). Studies have confirmed that laccase secreted in the secondary plant cell walls can catalyse lignin biosynthesis in the presence of oxygen (Sterjiades et al. 1992; Bao et al. 1993). Three laccase genes were observed to be more highly expressed in ‘Formosana’ than in ‘Brazil’. Caffeic acid O-methyltransferase is a key enzyme specific to the lignin pathway, and is involved in the synthesis of syringyl lignin (O’Malley et al. 1993; Ma and Xu 2008). Periplasmic beta-glucosidase is involved in the biosynthesis of phenylpropanoid compounds that contribute to lignin biosynthesis. Cellulose synthase and cellulose synthase catalytic subunit 12 are involved in the synthesis of cellulose, which is the main component of the cell walls of higher plants. Cellulose also helps regulate the structural resistance of plants. Our results revealed the importance of structural resistance in plant defences against Foc, and also verified the results of previous studies. Third, some of the identified genes encode pathogenesis-related (PR) proteins, which are produced when plants are infected by a pathogen or exposed to certain stressors. Chitinase, which can hydrolyse cell wall components of pathogenic fungi, is one of the major plant PR proteins. The observed higher expression levels of peroxidase and PR genes in ‘Formosana’ banana plants than in the ‘Brazil’ banana plants is consistent with the results of a previous study (Li et al. 2010).

Table 6. Specific resistance mechanism-related genes in banana.

Gene ID	A_FPKM	B_FPKM	log2fold_change	p_value	q_value	Gene description
CUFF.12751	48.4366	7.45885	−2.69907	2.07E-06	0.001876	cellulose synthase catalytic subunit 12 [Zea mays]
CUFF.16474	4.31267	0.363294	−3.56937	6.56E-05	0.024117	Periplasmic beta-glucosidase precursor, putative [Ricinus communis]
CUFF.1865	60.2792	4.8562	−3.63376	8.79E-08	0.000292	PREDICTED: peroxidase N [Vitis vinifera]
CUFF.20440	50.6218	6.64523	−2.92937	3.17E-07	0.000499	cellulose synthase [Betula luminifera]
CUFF.22126	13.7268	1.80391	−2.9278	3.70E-05	0.01672	laccase [Liriodendron tulipifera]
CUFF.25302	26.9838	4.76472	−2.50163	4.73E-05	0.019706	PREDICTED: laccase-4 [Vitis vinifera]
CUFF.29797	38.3655	4.26664	−3.16864	5.20E-07	0.000653	PREDICTED: laccase-4 [Vitis vinifera]
CUFF.4241	8.92429	0.425703	−4.38982	4.15E-05	0.018221	caffeic acid O-methyltransferase [Camellia sinensis]
CUFF.7118	69.1931	10.0007	−2.79053	5.08E-06	0.003759	putative chitinase [Epipremnum aureum]
CUFF.5976	46.4019	3.16654	−3.87321	2.30E-07	0.000429	similar to T-phylloplanin [Nicotiana tabacum]
CUFF.21873	171.312	1846.52	3.43011	5.02E-09	9.79E-05	RecName:Full = Snakin-1; Flags:Precursor snakin-1 [Solanum tuberosum]
CUFF.8570	30.5278	3.06931	−3.31414	3.42E-05	0.015689	TPA:class III peroxidase 136 precursor [Oryza sativa Japonica Group]
CUFF.13852	54.0644	5.85092	−3.20794	4.44E-06	0.003426	glutathione-S-transferase Cla47 [Elaeis guineensis]

A: ‘Formosana’; B: ‘Brazil’.

The expression levels of several phytoalexin synthesis-related genes were significantly higher in ‘Formosana’ than in ‘Brazil’ (Table 7). Phytoalexins, which are produced in response to diverse biotic and abiotic stress, are one of the most important physiological substances involved in plant defences. The most studied phytoalexins are the flavonoids and terpenoids. Anthranilate N-benzoyltransferase mediates the biosynthesis of phytoalexins (Reinhard & Matern 1989). Additionally, cytochrome P450 family proteins can catalyse the synthesis of defensive secondary metabolites and signalling molecules related to phytoalexin synthesis (Werck-Reichhart 1995; Morant et al. 2003). Flavonoid 3′,5′-hydroxylase and dihydroflavonol 4-reductase regulate the biosynthesis of flavonoids. Type-III polyketide synthase is a chalcone/stilbene synthase. Chalcone is a precursor of some of the flavonoid phytoalexins, and is important for the resistance of plants to diseases and parasites. Moreover, multiple flavonoid 3′,5′-hydroxylase and type-III polyketide synthase isoenzymes were higher in ‘Formosana’ than in ‘Brazil’. The phytoalexin synthesis-related genes whose expression levels were up-regulated in ‘Formosana’ may influence the observed tolerance to Fusarium wilt.

Table 7. Phytoalexin biosynthesis-related genes in banana.

Gene ID	A_FPKM	B_FPKM	log2fold_change	p_value	q_value	Gene description
CUFF.11375	16.2028	1.58274	−3.35575	6.78E-06	0.004599	Anthranilate N-benzoyltransferase protein, putative [Ricinus communis]
CUFF.24605	46.1247	5.77207	−2.99838	2.10E-05	0.011585	cytochrome P450 family protein [Musa acuminata]
CUFF.26496	31.277	5.34766	−2.54812	5.10E-05	0.020581	flavonoid 3ʹ,5ʹ-hydroxylase [Dendrobium hybrid cultivar]
CUFF.29055	125.182	22.5285	−2.4742	1.21E-05	0.007361	flavonoid 3ʹ,5ʹ-hydroxylase [Gossypium hirsutum]
CUFF.4245	21.3348	0.729475	−4.87021	1.64E-08	9.79E-05	dihydroflavonol 4-reductase [Anthurium andraeanum]
CUFF.16156	104.982	16.8295	−2.64109	3.66E-06	0.002952	type III polyketide synthase 4 [Musa acuminata AAA Group]
CUFF.21668	67.7461	10.6023	−2.67577	5.47E-06	0.003795	type III polyketide synthase 2 [Musa acuminata AAA Group]

A: ‘Formosana’; B: ‘Brazil’.

The expression levels of three stress resistance-related genes were significantly up-regulated in ‘Formosana’ (Table 8). Cysteine proteinase participates in responses to biotic and abiotic stresses, and is responsible for degrading stress-damaged proteins (Grudkowska and Zagdańska 2004). Nicotianamine synthase is critical for responses to iron deficiencies in monocots (Higuchi et al. 1996). The aluminium-activated malate transporter 10-like gene helps plants cope with aluminium stress (Hoekenga et al. 2006). The up-regulated expression of these genes in ‘Formosana’ may be related to the tolerance of these plants to Foc TR4.

Table 8. Stress resistance-related genes (non-specific responses to biotic and abiotic stresses) or genes related to environmental interactions in banana.

Gene ID	A_FPKM	B_FPKM	log2fold_change	P_value	q_value	Gene description
CUFF.20873	26.2729	2.6988	−3.28318	3.22E-06	0.002668	PREDICTED:nicotianamine synthase [Vitis vinifera]
CUFF.22946	9.15963	0.068897	−7.05471	1.95E-06	0.001823	cysteine proteinase precursor [Ipomoea batatas]
CUFF.25075	2.75965	0.044068	−5.96862	7.95E-05	0.027267	PREDICTED:aluminum-activated malate transporter 10-like [Brachypodium distachyon]
CUFF.29519	103.116	724.061	2.81184	2.63E-06	0.002313	arachidonic acid-induced DEA1 precursor [Solanum lycopersicum]
CUFF.97	10.4952	156.091	3.89459	1.25E-08	9.79E-05	PREDICTED:36.4 kDa proline-rich protein-like [Brachypodium distachyon]

A: ‘Formosana’; B: ‘Brazil’.

We determined that 16 of the genes differentially expressed between ‘Formosana’ and ‘Brazil’ affected plant growth and development (Table 9). These genes were functionally divided into three groups. The first group included some genes related to the assimilation, absorption, and transport of nitrogen (e.g. genes encoding a glutamine synthetase nodule isozyme and ammonium transporter 3 member 1-like). The second group included some genes involved in starch and sucrose metabolism (e.g. sugar transporter and cell wall invertase). Nitrogen and carbon metabolism form the basis of plant growth. The third group comprised genes encoding proteins mediating the growth and senescence of plants, including alpha-expansin 3, auxin-induced protein 5NG4 and a putative senescence-associated protein. These genes were expressed more highly in ‘Formosana’ than in ‘Brazil’, except for the gene encoding the senescence-associated protein. The differences in the gene expression levels between’Formosana’ and ‘Brazil’ were probably due to the differences in the development of Fusarium wilt symptoms and growth characteristics between the two varieties.

Table 9. Growth- and development-related genes in banana.

Gene ID	A_FPKM	B_FPKM	log2fold_change	P_value	q_value	Gene description
CUFF.10935	127.319	16.9563	−2.90856	4.62E-07	0.000627	putative thiamine biosynthesis protein [Oryza sativa Japonica Group]
CUFF.11605	66.3772	6.12059	−3.43894	1.99E-07	0.000423	Beta-expansin 3 precursor, putative [Ricinus communis]
CUFF.16771	108.372	21.125	−2.35897	2.94E-05	0.014471	Glutamine synthetase nodule isozyme [Vigna aconitifolia]
CUFF.18098	5.58978	0.499731	−3.48357	0.000104	0.03218	PREDICTED: ammonium transporter 3 member 1-like [Brachypodium distachyon]
CUFF.20530	3.21526	0.067103	−5.58242	9.99E-05	0.031406	PREDICTED: probable peptide/nitrate transporter At1g59740-like [Vitis vinifera]
CUFF.20569	14.2535	1.20206	−3.56774	6.57E-05	0.024117	PREDICTED: auxin-induced protein 5NG4 [Vitis vinifera]
CUFF.22363	5.4585	0.165955	−5.03964	9.36E-05	0.030066	PREDICTED: probable peptide/nitrate transporter At1g59740-like [Glycine max]
CUFF.2604	5.28788	0.31447	−4.07169	1.86E-05	0.010492	sugar transporter [Allium sativum]
CUFF.27654	23.7504	1.90324	−3.64142	6.01E-05	0.023312	germin-like protein [Musa acuminata AAA Group]
CUFF.415	14.6651	1.26488	−3.53531	5.25E-07	0.000653	inorganic pyrophosphatase [Nicotiana tabacum]
CUFF.597	14.3098	1.33667	−3.42028	2.31E-05	0.012322	alpha-expansin 3 [Populus tremula x Populus tremuloides]
CUFF.8264	19.4738	3.72506	−2.3862	0.000136	0.039947	Peptide transporter, putative [Ricinus communis]
CUFF.4717	11.5049	1.43439	−3.00375	6.62E-05	0.024117	cell wall invertase [Musa acuminata AAA Group]
CUFF.19956	8.65075	0.484782	−4.15742	5.04E-05	0.020581	PREDICTED: probable carboxylesterase 2-like [Vitis vinifera]
CUFF.29124	5537.25	39 630.5	2.83937	2.26E-05	0.012292	putative senescence-associated protein [Pisum sativum]
CUFF.5013	11.9841	1.12808	−3.40917	5.21E-06	0.003759	glycosyltransferase, CAZy family GT8 [Populus trichocarpa]

A: ‘Formosana’; B: ‘Brazil’.

According to our transcriptome analysis, there are several aspects to the mechanism underlying the resistance of 'Formosana’ to Foc TR4. First, 'Formosana’ exhibits up-regulated production of hydrolytic enzymes and PR proteins (e.g. chitinase). Consequently, there is increased hydrolysis of the cell wall of the invading Fusarium wilt pathogen and inhibition of fungal growth. Second, cell wall lignin production and the expression of cellulose synthesis-related genes are greater in ‘Formosana’ than in ‘Brazil’. These differences result in a stronger ‘Formosana’ cell wall that is better able to impede the invasion by Foc TR4. Third, the expression levels of multiple phytoalexin biosynthesis-related genes are up-regulated in 'Formosana’, implying that phytoalexins contribute to the tolerance of ‘Formosana’ to Fusarium wilt. Fourth, ‘Formosana’ exhibits up-regulated expression of genes encoding peroxidases and glutathione-S-transferase Cla47. The increased production of these enzymes enhances the removal of excessive ROS. Fifth, the expression levels of some stress-resistance genes are higher in ‘Formosana’ than in ‘Brazil’. Therefore, the tolerance of ‘Formosana’ to wilt disease may be the result of synergistic interactions among several stress-responsive and/or disease resistance genes. Our results may be useful for characterizing the molecular mechanism responsible for the disease resistance of AAA-type banana cultivars under natural conditions. Our data may also serve as a valuable transcriptomics resource for future functional genomics studies of economically important banana cultivars.

Acknowledgements

This work was supported by the National Natural Science Foundation of China [31471738]; the Earmarked Fund for China Agriculture Research System [CARS-32-04]; and the Chinese Agriculture Ministry Project for Tropical Crop Pest Monitoring and Control [15RZBC-16].

Supplemental material

Supplemental data for this article can be accessed online here: https://doi.org/10.1080/07060661.2017.1342693.