Introduction
The dependence of parasites on their hosts and the high virulence of many parasites make the host–parasite relationship one of the strongest biotic interactions. Due to the diverse physiologies and molecular phenotypes of different hosts, these intimate relationships pose a challenge for generalist parasite species that attack multiple host species. Parasites with broad host ranges often display high levels of intraspecific variation (Archie and Ezenwa Reference Archie and Ezenwa2011; Kaci-Chaouch et al. Reference Kaci-Chaouch, Verneau and Desdevises2008). This may be partly maintained by balancing selection within populations or local adaptation across populations on different hosts (Frank Reference Frank1993). Within-host competitive interaction can also act as a maintaining force of intraspecific genetic diversity (Bashey Reference Bashey2015). Alternatively, in generalist parasites, plastic responses to hosts can be key to successful host generalism (Lajeunesse and Forbes Reference Lajeunesse and Forbes2002; Leggett et al. Reference Leggett, Benmayor, Hodgson and Buckling2013). While gene-for-gene and simple genetic architectures of host–parasite interactions have been frequently studied (Guo et al. Reference Guo, Tian, Wamboldt and Alfano2009; Watson Reference Watson1970), we know less about host–parasite interactions with a complex genetic basis, such as those between parasitic plants and their hosts. As many parasitic plants are major agricultural pests, an important goal is to understand the genetics of parasitism across diverse host crops (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013; Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015).
Purple witchweed [Striga hermonthica (Delile) Bent.] (Orobanchaceae) is an obligate root hemiparasite of grasses in sub-Saharan Africa causing billions of dollars in crop loss annually (Spallek et al. Reference Spallek, Mutuku and Shirasu2013). Striga hermonthica commonly attacks grain sorghum [Sorghum bicolor (L.) Moench], maize (Zea mays L.), and pearl millet [Pennisetum glaucum (L.) R. Br.], with increasing incidence of parasitism now being observed in sugarcane (Saccharum officinarum L.) and upland rice (Oryza sativa L.) fields (Spallek et al. Reference Spallek, Mutuku and Shirasu2013). Crop infestations of parasitic weeds have devastating effects on crop yield around the world, especially in Africa (Rodenburg et al. Reference Rodenburg, Demont, Zwart and Bastiaans2016). Farmers have reported losses between 20% and 80%, and high-density infestations might lead to the total loss of the crop and abandonment of a field (Atera et al. Reference Atera, Itoh, Azuma and Ishii2012). Particularly for Striga spp., weed control is challenging, because 1) the severe reaction of hosts to parasitism initiates almost immediately after haustorial attachment, negating the value of weeding of aboveground parasites as a crop-loss mitigation strategy; 2) mature plants produce tens to hundreds of thousands of tiny seeds that are readily dispersed by wind, water, contaminated grain, livestock, and people; and 3) Striga seeds are long-lived, capable of establishing a persistent seedbank. Striga hermonthica exhibits genetic variation in interactions both among and within crops with variability in host preference and virulence, hindering the development of resistant cultivars (Bozkurt et al. Reference Bozkurt, Muth, Parzies and Haussmann2015; Koyama Reference Koyama, Rattunde, Geiger, Haussman, Hess, Koyama and Grivet2000; Olivier et al. Reference Olivier, Glaszmann, Lanaud and Leroux1998; Unachukwu et al. Reference Unachukwu, Menkir, Rabbi, Oluoch, Muranaka, Elzein, Odhiambo, Farombi and Gedil2017). Reports of tolerant and/or resistant genotypes vary among crop species. For example, several grain sorghum genotypes have been shown to be tolerant or resistant to Striga (Haussmann et al. Reference Haussmann, Hess, Omanya, Folkertsma, Reddy, Kayentao, Welz and Geiger2004; Hess et al. Reference Hess, Ejeta and Butler1992; Maiti et al. Reference Maiti, Ramaiah, Bisen and Chidley1984; Vogler et al. Reference Vogler, Ejeta and Butler1996). Meanwhile, tolerance or resistance in maize has rarely been reported (Amusan et al. Reference Amusan, Patrick, Menkir, Housley and Ejeta2008; Lane et al. Reference Lane, Child, Moore, Arnold and Bailey1997; Mutinda et al. Reference Mutinda, Masanga, Mutuku, Runo and Alakonya2018).
Genetic and phenotypic responses to host variation in parasitic plants of the Orobanchaceae family are likely complex, given the multiple developmental stages to interactions and mechanisms of host defenses that also vary among host species. In parasitic plants like S. hermonthica, researchers are only beginning to uncover the genetic and phenotypic basis of host-specific interactions (Ejeta Reference Ejeta2007; Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013; Swarbrick et al. Reference Swarbrick, Huang, Liu, Slate, Press and Scholes2008). Expression profiling of genes expressed at the host–parasite interface suggest that one possible strategy for parasitic plants is to express a set of shared parasitic genes irrespective of the host and a different set of genes that are expressed on specific hosts (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013). Yet, other studies have found limited evidence of genetic differentiation between S. hermonthica from different hosts of origin (Unachukwu et al. Reference Unachukwu, Menkir, Rabbi, Oluoch, Muranaka, Elzein, Odhiambo, Farombi and Gedil2017; Welsh and Mohamed Reference Welsh and Mohamed2011). The dynamic and complex nature of variation in host-specific interactions in natural S. hermonthica populations challenges our understanding of its molecular basis.
Transcriptome profiling offers the opportunity to rapidly study functional variation in the genome and to connect this variation to specific genetic loci, evolutionary processes, and environmental conditions. By studying transcriptomic variation in individuals from natural populations, we can build a more mechanistic understanding of host–parasite interactions and their impacts on populations and communities (McGill et al. Reference McGill, Enquist, Weiher and Westoby2006). Recent examples demonstrate the utility of transcriptomic studies of wild organisms in natural systems (Swenson et al. Reference Swenson, Iida, Howe, Wolf, Umaña, Petprakob, Turner and Ma2017). Furthermore, gene expression analyses can help reveal complex mechanisms of host–parasite interactions (Ichihashi et al. Reference Ichihashi, Mutuku, Yoshida and Shirasu2015). Gene expression analysis in the parasite yellowbeak owl’s-clover (Triphysaria versicolor Fisch. & C.A. Mey.) showed how this parasite relies on a generalist response in which functionally overlapping but distinct gene sets are expressed based on the parasitized host (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013). In a comprehensive comparative transcriptomic study of individuals in the Orobanchaceae family, Yang et al. (Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015) identified a set of core parasitism genes used by members of three parasite species with different extents of nutritional dependence. The core parasitism genes belonged to specific members of gene families encoding cell wall–modifying enzymes and peroxidases, proteins known to be part of the parasitic invasion process (Antonova and Terborg Reference Antonova and Terborg1996; Losner-Goshen et al. Reference Losner-Goshen, Portnoy, Mayer and Joel1998; Pérez-de-Luque Reference Pérez-de-Luque, Joel, Gressel and Musselman2013; Singh and Singh Reference Singh and Singh1993).
Despite prior efforts aimed at characterizing host-specific responses of S. hermonthica, genomic studies have been limited by the lack of functional genomic data (i.e., beyond reduced-representation sequencing). Transcriptomes of S. hermonthica have been studied across developmental stages and tissues on a single host under laboratory conditions (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015; Parasitic Plant Genome Project, http://ppgp.huck.psu.edu). However, parasite transcriptomic profiling across hosts may be used to build a functional understanding of host specificity. Furthermore, growth under controlled laboratory conditions may not fully represent functional diversity of parasite response in natural environments. Here, we leverage transcriptomic data obtained from natural populations of S. hermonthica parasitizing different host species to investigate gene expression variation across hosts. Specifically, we analyze transcriptomes from S. hermonthica growing in the field on maize and grain sorghum to address two central questions about host–parasite interactions: 1) to what extent are S. hermonthica populations on different hosts divergent in genotype and expression phenotypes, and 2) which, if any, parasite genes exhibit host-specific expression?
Materials and methods
Study site and sampling
Field populations of S. hermonthica were sampled from farms located around the town of Mokwa in Niger State (Supplementary Figure 1; Supplementary Table 1). Mokwa is located at an elevation of 160 m above sea level and has a tropical savanna climate (Aw on the Koppen climate classification) (Geiger Reference Geiger, Landolt and Börnstein1954). The average temperature and precipitation are 27.4 C and 1,056 mm yr−1, respectively (Karger et al. Reference Karger, Conrad, Böhner, Kawohl, Kreft, Soria-Auza, Zimmermann, Linder and Kessler2017). More than 95% of annual rain tends to fall during a rainy season from April to October, which is followed by a dry period from November to March.
Heat map showing the expression pattern and clustering of the 38 differentially expressed transcripts (P < 0.01 and log fold change > 2) with at least 60% of their variation explained by host in Striga hermonthica. The dendogram shows that the two main branches cluster samples by host (sorghum vs. maize). Colors correspond to log2CPM (counts per million).
Top 1% of Striga hermonthica single-nucleotide polymorphisms (SNPs) with highest FST between hosts
a Transcripts with nonconservative amino acid changes are indicated in bold.
b Position of the SNPs within the transcript. Positions with a nonsynonymous amino acid change are indicated with an asterisk (*).
c RefAvFS, averaged frequency of the reference allele in a population with sorghum host.
d RefAgFM, averaged frequency of the reference allele in population with maize host.
e FST is a measure of structure in natural populations based on expected heterozygosity for each host-specific population relative to the total population.
In November 2009, we sampled 14 populations of S. hermonthica from 14 fields: 3 growing on maize and 11 growing on grain sorghum (Supplementary Figure 1; Supplementary Table 1). From each population, we pooled young leaves from 10 to 20 individuals and, separately, we pooled flower buds from the same individuals. The tissues were immediately flash frozen in liquid nitrogen in the field.
RNA extraction, library preparation, and sequencing
For each pooled collection (leaf or flower bud tissue), the entire collection was emptied into a prechilled mortar and thoroughly mixed in a pregrinding step before sampling (and further grinding) the amount needed for single RNA isolations. RNA isolations were carried out using Ambion’s RNAqueous™ Midi Kit with Plant Isolation Aid (Thermo Fisher Scientific, Waltham, MA) according to the manufacturer’s recommendations. The quality of the resulting RNA was assessed using an Agilent Bioanalyzer (RNA integrity number [RIN] values ≥ 7.4, based on a scale of 1 to 10, with 1 being most degraded and 10 most intact). The RNA was precipitated (with 0.1 volumes sodium acetate and 3 volumes of 100% ethanol) to concentrate the RNA. Samples were DNase treated, column purified, and evaluated on the Agilent Bioanalyzer (Agilent Technologies, Santa Clara, CA) for quality and quantity. Finally, leaf and flower bud RNA samples were pooled in equimolar amounts for each of the 14 field populations. To make a cDNA library, 10 μg of each RNA pool was used.
Indexed cDNA libraries for paired-end sequencing were made using Illumina’s mRNA Seq Sample Prep Kit and Multiplexing Sample Prep Oligonucleotide Kit (Illumina, San Diego, CA), with the size selection targeting an approximately 185-bp insert size. Because library concentrations were low after the final gel purification step, 3 to 12 additional cycles of amplification (with effective increases in library concentrations ranging from 10X to 100X) were performed to reach the target concentration of 10 nM before sequencing. Finally, libraries were validated by Agilent Bioanalyzer (Agilent Technologies, Santa Clara, CA) assay and qPCR. The qPCR validation was performed using KAPA SYBR® FAST qPCR Master Mix (F. Hoffmann-La Roche AG, Basel, Switzerland) and Illumina PhiX control DNA (PhiX control kit v. 2.0, CT-901-2001, Illumina, San Diego, CA) with a size-correction calculation according to the Kapa Biosystem manual. Libraries were sequenced as 75-bp PE (paired-end sequencing lanes) in two lanes (eight libraries per lane) on the Illumina GAIIx platform (Illumina, San Diego, CA). Each library was sequenced at the same sequencing depth in two lanes, and the resulting files were combined.
Data preprocessing and de novo assembly
Raw paired-end reads were trimmed and filtered using Trimmomatic v. 0.36 to remove adapter sequences, short reads, and low-quality reads from raw sequence data (Bolger et al. Reference Bolger, Lohse and Usadel2014). To remove adapter sequences, we looked for seed matches with a maximum of three mismatches and an initial length of 10. These seeds were extended and clipped when paired-end reads had a score of 30 or lower. We removed leading and trailing bases of low quality (Phred quality score <20) or N’s. Likewise, we cut and removed the 3′ end of reads if the quality of a 4-base-wide sliding window dropped below 15. Finally, we discarded any reads remaining that were shorter than 50 bases long. Filtered reads were inspected using FastQC v. 0.11.5 (http://www.bioinformatics.babraham.ac.uk/projects) to confirm that they met our quality standards.
We used the S. hermonthica de novo transcriptome assembly (StHeBC4) as a reference for both expression and single-nucleotide polymorphism (SNP) discovery analyses from the Parasitic Plant Genome Project (PPGP) website (http://ppgp.huck.psu.edu/download.php; http://bigdata.bx.psu.edu/PPGP_II_data). In brief, the reference de novo assembly was performed with Trinity (Grabherr et al. Reference Grabherr, Brian, Yassour, Levin, Thompson, Amit, Adiconis, Fan, Raychowdhury, Zeng, Chen, Mauceli, Hacohen, Gnirke and Rhind2011) using RNA-Seq data from the following parasite developmental stages: roots of germinated seedlings after exposure to GR24, roots of germinated seedlings after exposure to host roots, expression in haustoria attached to host roots (prevascular connection), haustoria attached to host roots after vascular connection, late postattachment from belowground plants, stem and leaf tissue, and flowers (Westwood et al. Reference Westwood, dePamphilis, Das, Fernández-Aparicio, Honaas, Timko, Wafula, Wickett and Yoder2012). The assembled contigs were postprocessed into nonredundant sets including predicted coding sequences and their corresponding translations using the PlantTribes pipeline (https://github.com/dePamphilis/PlantTribes). They were also screened against the grain sorghum genome and the National Center for Biotechnology Information (NCBI) nonredundant proteins database (nr) to remove host and other contaminant sequences. Detailed information about the PPGP parasitic plants reference transcriptome assemblies is provided on the PPGP website and described in Yang et al. (Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015).
SNP discovery and allelic frequency analyses
We further sought to understand host-specific differences in S. hermonthica parasitizing different hosts by characterizing SNPs in aboveground expressed transcripts. Trimmed reads were mapped to the de novo assembled transcriptome using aln and sampe commands in BWA, recommended for 75-bp paired-end reads (Li and Durbin Reference Li and Durbin2009). We converted the SAM files containing the mapped reads to BAM with SAMtools v. 1.5 (Li et al. Reference Li, Handsaker, Wysoker, Fennell, Ruan, Homer, Marth, Abecasis and Durbin2009) and used Picard Tools v. 2.8.2 (https://github.com/broadinstitute/picard) to sort and remove duplicate mappings. We filtered the sorted and deduplicated BAM files with SAMtools v. 1.5 and kept only proper pairs with map quality greater than 30. We used the BAM files containing high-quality pair-ended mappings to produce an mpileup file in SAMtools with minimum base quality of 20.
We next identified host-associated allelic variation. Allelic variation in transcriptomes can be due to host-associated allelic frequency or allele-specific expression. Both processes are of interest for identifying loci associated with host specificity. It has been shown that SNPs can be accurately called from pooled transcriptome data using conservative parameters to select sites for allelic frequency estimations (Konczal et al. Reference Konczal, Koteja, Stuglik, Radwan and Babik2013). Thus, we called genotypes in BCFtools for biallelic SNPs with a minimum coverage of 5X and removed SNPs with coverage less than 5X in more than one sample. To investigate how the individual samples clustered based on their genotype, we performed a principal components analysis (PCA) with the R package SNPRelate using the default options (Zheng et al. Reference Zheng, Levine, Shen, Gogarten, Laurie and Weir2012a). For allelic frequency estimation, we only used SNPs that were present in all samples. Allele frequencies at the remaining polymorphic sites were calculated for each library with VCFtools based on counts of reference and alternate alleles. We grouped libraries into two populations: a population for host maize (n = 3) and a population for host grain sorghum (n = 11). For each SNP, we computed fixation index (F ST) based on expected heterozygosity for each host-specific population relative to the total population using a custom perl script (Nei Reference Nei1977). We selected SNPs in the upper 1% tail (outliers) of F ST between host populations. To gauge whether differentiation of these “outlier” SNPs exceeded that expected by chance, we calculated F ST values for three grain sorghum populations located nearest the maize populations (Sorghum_1, Sorghum_5, and Sorghum_11) versus the eight remaining sorghum populations. The sorghum versus sorghum comparison average F ST was 0.017, and the maize versus sorghum comparison average F ST was 0.016. The distribution of F ST values for this sorghum versus sorghum analysis was compared with the distribution of F ST values obtained for the maize populations versus the eight grain sorghum populations (excluding Sorghum_1, Sorghum_5, and Sorghum_11).
We annotated candidate outlier SNPs following the same steps used for the differentially expressed transcripts and further characterized them with respect to their predicted impact on coding sequence. We aligned the transcripts containing the outlier SNPs to their predicted coding sequence (CDS) using Geneious 9.1.8 to determine whether the SNP was located within the coding region and whether it coded for a synonymous or nonsynonymous mutation. Nonsynonymous changes were categorized as conservative or nonconservative (see “Results and Discussion”), and we used the score from the substitution matrix BLOSUM-62 to indicate the structural similarity between the amino acids.
Differential expression analysis
We used the alignment-free package Salmon to perform differential expression (DE) analysis (Patro et al. Reference Patro, Duggal, Love, Irizarry and Kingsford2017). Salmon is a tool for expression quantification that couples the concept of quasi-mapping with a two-phase inference procedure. To characterize host-related expression patterns, the output matrix containing count estimates from Salmon were analyzed with edgeR (Robinson et al. Reference Robinson, McCarthy and Smyth2010). edgeR identifies differentially expressed transcripts based on the assumption that the number of reads produced by each transcript is proportional to its abundance. Because Salmon produces estimates instead of total counts, we followed the guidelines and rounded the estimates for downstream analysis with edgeR. Populations of S. hermonthica parasitizing the same host crop (maize or grain sorghum) were considered biological replicates, with 3 replicates for maize and 11 for grain sorghum. Before performing the DE analysis, we excluded transcripts with very low expression values using a cut-off of 1 CPM (counts per million) in at least three samples (number of biological replicates for libraries with maize as host crop), and we used the trimmed mean of M-values method for normalizing the counts. We used an exact test implemented in edgeR, using host as group, in the transcripts remaining after filtering. We focused on those transcripts that were differentially expressed with a log fold change (logFC) significantly greater than 2 and a P-value < 0.01.
To estimate the importance of host-associated expression genome-wide, we implemented variation partitioning in the R package Vegan 2.3-5 (Oksanen et al. Reference Oksanen, Blanchet, Kindt, Legendre, Minchin, O’Hara, Simpson, Solymos, Stevens and Wagner2013). Variation partitioning can be used to estimate the proportion of variance explained by covariates (host) with multivariate responses (transcriptomes). We also tested whether the effect of host (proportion of transcriptomic variation explained) was significant genome-wide using permutation tests that randomized host identity across libraries. Finally, we used the R package VariancePartition (Hoffman and Schadt Reference Hoffman and Schadt2016) to quantify the contribution of host to variance explained at each level of individual transcript.
To focus on transcripts with a significant and large effect size of host on expression, we combined the results from the DE and VariancePartition analysis. We selected significantly differentially expressed transcripts (P < 0.01 and logFC > 2) with at least 60% of the transcript expression variation explained by host for annotation with Blast2GO (Conesa and Götz Reference Conesa and Götz2008; Conesa et al. Reference Conesa, Götz, García-Gómez, Terol, Talón and Robles2005). We annotated these focal transcripts using the Basic Local Alignment Search Tool (blastx) (e-value cutoff 1.0 × 10−5) against the nr database from the NCBI (Pruitt et al. Reference Pruitt, Tatusova and Maglott2007) and InterPro from the European Bioinformatics Institute (Finn et al. Reference Finn, Attwood, Babbitt, Bateman, Bork, Bridge, Chang, Dosztányi, El-Gebali, Fraser, Gough, Haft, Holliday, Huang and Huang2017). We ran an enzyme code search and a Kyoto Encyclopedia of Genes and Genomes pathway analysis (https://www.genome.jp/kegg; Kanehisa and Goto Reference Kanehisa and Goto2000). Then, we retrieved gene ontology (GO) terms (http://www.geneontology.org; Ashburner et al. Reference Ashburner, Ball, Blake, Botstein, Butler, Cherry, Davis, Dolinski, Dwight, Eppig, Harris, Hill, Issel-Tarver, Kasarskis and Lewis2000; Gene Ontology Consortium 2017) from the blastx, InterPro, and enzyme hits. We extracted functional gene categories for each transcript, comparing its translated sequence with the Embryophyta orthologue database in OrthoDB v. 9.1 (Zdobnov et al. Reference Zdobnov, Tegenfeldt, Kuznetsov, Waterhouse, Simão, Ioannidis, Seppey, Loetscher and Kriventseva2017). Finally, we searched for noncoding RNA families and other structured RNA elements using the Rfam database (Griffiths-Jones et al. Reference Griffiths-Jones, Bateman, Marshall, Khanna and Eddy2003).
Our libraries captured differentially expressed transcripts in aboveground tissue of parasites but did not contain information on other tissues relevant for host–parasite interactions. To investigate tissue specificity of the focal DE transcripts, we used S. hermonthica RNA-Seq expression data across parasite developmental stages obtained from the PPGP v. 2 (PPGP II), a large database of replicated expression data for species of Orobanchaceae, including S. hermonthica collected in Borno State (Nigeria) that was growing on grain sorghum (http://bigdata.bx.psu.edu/PPGP_II_data). The developmental stages represented in this data are: imbibed seeds (stage 0), roots of germinated seedlings after exposure to GR24 (stage 1), roots of germinated seedlings after exposure to host roots (stage 2), expression in haustoria attached to host roots (prevascular connection) (stage 3), haustoria attached to host roots after vascular connection (stage 4), late postattachment from belowground plants (stage 5), stem and leaf tissue (stage 6.1), and flower tissue (stage 6.2) (Westwood et al. Reference Westwood, dePamphilis, Das, Fernández-Aparicio, Honaas, Timko, Wafula, Wickett and Yoder2012).
Results and discussion
Data preprocessing and de novo assembly
The total number of raw reads across the 14 libraries of S. hermonthica ranged from 2,837,888 in sample St_He_NG_12 to 7,749,848 in St_He_NG_2 (Supplementary Table 1) (raw sequencing data available at http://bigdata.bx.psu.edu/PPGP_II_data/StHeBC4 and in the Dryad repository https://doi.org/10.5061/dryad.s98q453). Likewise, the minimum and maximum number of high-quality clean reads were found in St_He_NG_12 and St_He_NG_2 (2,694,360 and 7,381,916, respectively) (Supplementary Table 1). Quality filtering of individual reads retained an average of 93.11% of read pairs.
After assembly of the quality trimmed reads with Trinity (Grabherr et al. Reference Grabherr, Brian, Yassour, Levin, Thompson, Amit, Adiconis, Fan, Raychowdhury, Zeng, Chen, Mauceli, Hacohen, Gnirke and Rhind2011), a series of postprocessing steps were used to reduce the assembly to a single transcript per putative locus (see Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015; Supplementary Table 2). Out of 94,969 postprocessed S. hermonthica de novo assembly transcripts, we selected 59,243 putative locus representatives to use as the reference transcriptome for subsequent SNP discovery and expression analyses.
Annotation of Striga hermonthica transcripts with highest FST between hosts containing nonsynonymous single-nucleotide polymorphisms
a Transcripts with nonconservative amino acid changes are indicated in bold, and transcripts with gene functions that match the core parasitism gene categories (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015) are noted with an asterisk (*).
b Gene description obtained from the best match (based on the e-value) in the Basic Local Alignment Search Tool (blastx).
c Functional gene categories obtained from the orthologue comparison in the OrthoDB.
Host-associated Allelic differentiation analysis: low genome-wide divergence in expressed allelic frequency between hosts
Before SNP calling, we mapped clean reads to the reference transcriptome with an average mapping rate of 83.67% (ranging from 79.66% to 90.19%). The initial SNP calling contained 401,544 SNPs, but after filtering for biallelic SNPs and minimum coverage of 5X, we retained 4,806 SNPs distributed across 899 transcripts. A low SNP retention rate primarily resulted from coverage of transcripts expressed at low levels. The set of filtered SNPs was used to investigate how the individual samples clustered based on their genotype using a PCA. The first six eigenvectors captured 51.4% of the data variation, with each eigenvector explaining a similar amount of variation (9.7%, 8.8%, 8.4%, 8.3%, 8.1%, and 7.9%). The PCA clustering did not mirror geographic distribution of S. hermonthica samples, and clustering with respect to host was not observed (Supplementary Figure 2). After removing non-polymorphic SNPs and SNPs with missing values, we calculated pairwise FST values for each SNP. We focused on the top 1% of sites with highest differentiation between groups from different host species (Supplementary Figure 3), which yielded 56 SNPs distributed in 50 transcripts (Table 1). For the 56 SNPs, difference in reference allelic frequency between populations from maize versus grain sorghum ranged from 12.5% to 45.4%.
Heat map showing expression in transcripts per million (TPM) from tissue-specific data in different developmental stages of Striga hermonthica. Data were extracted from the Parasitic Plant Genome Project II (http://ppgp.huck.psu.edu/). Stage 0, expression in imbibed seeds; Stage 1, expression in roots of germinated seedlings after exposure to GR24; Stage 2, expression in roots of germinated seedlings after exposure to host roots; Stage 3, expression in haustoria attached to host roots (prevascular connection); Stage 4, expression in haustoria attached to host roots after vascular connection; Stage 5, late postattachment from belowground plants; Stage 6.1, expression in stem and leaf tissue; Stage 6.2, expression in flower tissue. Values correspond to the averaged relative expression from three replicates in each stage.
The low FST between hosts averaged across SNPs found in expressed transcripts (0.012) indicated minimal genome-wide divergence across the study region (~180 km2). Insect-pollinated species like S. hermonthica (Safa et al. Reference Safa, Jones and Musselman1984) can display high gene flow among populations (Ahmed et al. Reference Ahmed, Stephen, Roger and Philip2009). This and the high spatial heterogeneity in available hosts (adjacent farm fields with different crops) are likely to promote extensive gene flow between S. hermonthica parasites on different hosts. Seed distribution from neighboring fields via wind and water and distribution from near and more distant locations via contaminated crop seed, animals, equipment, and farmers would also contribute to lower population differentiation. Moreover, in the study region, farmers frequently rotate crop species between growing seasons, facilitating gene flow of S. hermonthica across hosts. Eight of 14 populations in this study were growing with another crop species (intercrop) when sampled (Supplementary Table 1). Of these eight, two populations were from fields where maize and grain sorghum were intercropped, and only one population was from a field planted to the same crop (grain sorghum) for the 3 yr before sampling. Altogether, we anticipated low genome-wide divergence between S. hermonthica populations from different hosts. Although FST values calculated for individual SNPs between host groups were generally low, we also found SNPs with moderate to high differentiation (i.e., exceeding ~0.2). We cannot rule out the role of neutral processes underlying high levels of differentiation among these transcripts (Supplementary Figure 4). However, high host-associated allelic differentiation could indicate host-specific selection over generations, causing strong differentiation at loci important for host specificity. Because parasites initiate the infection process belowground, it could be expected that selection on host-specific responses would be strongest during initial stages of parasite germination, development, and attachment (Conn et al. Reference Conn, Bythell-Douglas, Neumann, Yoshida, Whittington, Westwood, Shirasu, Bond, Dyer and Nelson2015). Additionally, allelic differentiation could result if seeds of S. hermonthica genotypes germinate in response to certain host genotypes. Following attachment, different physiological/metabolic environments presented by different hosts may also impose selection on loci important for physiological integration at the vegetative stage.
To understand the role that the transcripts containing FST outlier loci had in S. hermonthica, we predicted the impact that the change of allele had on the coding sequence and annotated the coding sequences. Six out of the 56 outlier SNPs were located outside a predicted CDS (Table 1). For the 50 SNPs located within a CDS, we identified 32 synonymous and 18 nonsynonymous substitutions (Table 1). Among the nonsynonymous substitutions, 9 were conservative and 9 were nonconservative changes. While all nonsynonymous mutations lead to an amino acid change in the resulting protein, the properties of the amino acid remain unchanged in conservative substitutions, while the properties of the resulting amino acid differ in nonconservative substitutions, potentially resulting in protein functional changes. We retrieved annotation for 45 out of the 50 transcripts containing highly differentiated SNPs (Table 2; Supplementary Table 3). To identify relevant genes for the host–parasite interaction, we focused on the functional gene categories and the GO terms extracted from the annotation. The three most common GO terms (level 3) identified during the annotation process for the category biological processes were cellular metabolic process, organic substance metabolic process, and nitrogen compound metabolic process. In the category molecular function, we found the terms protein binding, ion binding, and organic cyclic compound binding. In the cellular component category, the three most common GO terms were intracellular part, intracellular organelle, and membrane-bounded organelle. The functional annotations showed two core parasitic gene categories previously reported as enriched among core parasitism genes and formation of haustoria (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015): transporters and cell wall biogenesis (Supplementary Table 3). The most abundant category was transporters, with five transcripts containing nonconservative substitutions (Table 2). Finally, the enzyme code search indicated that the most common categories were oxidoreductases and lyases (Supplementary Table 3). When we checked pathway annotations for transcripts from the FST analyses, we found 34 pathways that included 22 enzymes, the most common ones were purine and thiamine metabolism (Supplementary Table 4).
Transcripts of Striga hermonthica detected as differentially expressed with P < 0.01 and log fold change (logFC) > 2 and variation explained by host > 60%
a Transcripts not found in the tissue-specific data are highlighted in bold, and those belonging to core parasitism gene sets (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015) are denoted with an asterisk (*).
b Host in which the transcript was identified as upregulated.
c Gene description obtained from the best match (based on the e-value) in the Basic Local Alignment Search Tool (blastx).
d Functional gene categories obtained from the orthologue comparison in the OrthoDB.
e Developmental stages where this transcript was found to have the highest expression in the data from the Parasitic Plant Genome Project (PPGP II).
Genes related to nutrient transport, metabolism, and defense may be among those most important for regulating physiological integration at the vegetative stage. Among the top 1% of sites most genetically differentiated between groups were genes that reside in the functional category nutrient transport (Table 2). Two candidate transcripts with nonconservative changes encoded a predicted serine acetyltransferase and a dnaJ homologue subfamily B protein. Besides their role in nutrient transport, these transcripts are also associated with energy production and conversion. An ABC transporter gene associated with ABA-activated signaling and active transport was detected by the expression analysis. Parasitic Orobanchaceae have been shown to recruit plant regulatory machinery associated with multiple phytohormones to perform parasite-specific functions (Tomilov et al. Reference Tomilov, Tomilova, Abdallah and Yoder2005). Moreover, in a comparative transcriptomic study of three Orobanchaceae parasites, an ABC transporter was upregulated (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). Finally, ABA has been shown to inhibit the germination of the obligate parasite branched broomrape [Phelipanche ramosa (L.) Pomel] (Lechat et al. Reference Lechat, Pouvreau, Péron, Gauthier, Montiel, Véronési, Todoroki, Le Bizec, Monteau, Macherel, Simier, Thoiron and Delavault2012; Pouvreau et al. Reference Pouvreau, Gaudin, Auger, Lechat, Gauthier, Delavault and Simier2013). A second ABA-associated transcript containing a conservative substitution is a homologue of a primary amine oxidase with quinone-binding activity. In other root parasitic plants, a quinone oxidoreductase was found to be necessary for haustorial initiation (Bandaranayake et al. Reference Bandaranayake, Filappova, Tomilov, Tomilova, Jamison-McClung, Ngo, Inoue and Yoder2010; Ngo et al. Reference Ngo, Albrecht, Tsuchimatsu and Grossniklaus2013). Another transcript containing a nonconservative substitution was a WRKY transcription factor associated with defense responses. This transcript encodes a promising candidate for host-specific interactions, because a member of this family, WRKY45, has been found to be associated with resistance to S. hermonthica infection in rice (Mutuku et al. Reference Mutuku, Yoshida, Shimizu, Ichihashi, Wakatake, Takahashi, Seo and Shirasu2015).
We also identified a set of candidate transcripts with synonymous mutations between plants parasitizing maize and grain sorghum (Table 1). Traditionally, synonymous mutations were thought to have no consequence for phenotype or fitness. However, there is growing evidence that synonymous mutations can affect the phenotype and fitness by altering mRNA stability and transcript splicing (Chamary and Hurst Reference Chamary and Hurst2005; Kudla et al. Reference Kudla, Murray, Tollervey and Plotkin2009) and therefore can be subject to selection (Bailey et al. Reference Bailey, Hinz and Kassen2014; Hunt et al. Reference Hunt, Simhadri, Iandoli, Sauna and Kimchi-Sarfaty2014; Ingvarsson Reference Ingvarsson2010). One transcript of potential interest was an ABA-associated, synaptosomal SNAP25-like homologous to SNAP33. In Arabidopsis thaliana, loss of function in SNAP25-type proteins leads to large necrotic lesions on leaves, resulting in a lethal dwarf phenotype (Heese et al. Reference Heese, Gansel, Sticher, Wick, Grebe, Granier and Jürgens2001). Moreover, SNAP25-like proteins are associated with defense responses, and the expression of SNAP33 increases response to pathogenic infection and mechanical stimulation (Wick et al. Reference Wick, Gansel, Oulevey, Page, Studer, Dürst and Sticher2003). We also found two transcripts with synonymous mutations that were associated with response to salicylic acid (SA): salicylate 3-hydroxylase DMR6-like and a cell wall kinase receptor. When attacked by a pathogen, plants require SA for the initiation of a hypersensitive response and to trigger systemic acquired resistance (Durrant and Dong Reference Durrant and Dong2004; Fu and Dong Reference Fu and Dong2013). Moreover, nonparasitic plants pathogens have been show to modify SA concentrations in the host for their own benefit (Cui et al. Reference Cui, Bahrami, Pringle, Hernandez-Guzman, Bender, Pierce and Ausubel2005; Zheng et al. Reference Zheng, Spivey, Zeng, Liu, Fu, Klessig, He and Dong2012b). The SA-related cell wall kinase may also be an interesting candidate, because cell wall–associated receptor kinases are transmembrane proteins with an extracellular domain that binds cell wall–associated pectins. Members of this family are known to be involved in disease resistance, hormone signaling, legume–rhizobium symbiosis, senescence, abiotic stress, and wounding stress responses (Bethke et al. Reference Bethke, Thao, Xiong, Li, Soltis, Hatsugai, Hillmer, Katagiri, Kliebenstein, Pauly and Glazebrook2016; Nishiguchi et al. Reference Nishiguchi, Yoshida, Sumizono and Tazaki2002).
Host-specific, differentially expressed transcripts in aboveground tissue—transporters, cell wall modification, and defense pathogenesis–related and oxidation–reduction processes
After we removed weakly expressed transcripts and normalized read counts, the count matrix retained 45,733 of the 59,243 transcripts in the de novo reference assembly. Across all 45,733 transcripts, host species explained 8.43% of the total variation in expression, although permutations indicated this was nonsignificant (P = 0.379, 999 permutations). This absence of significance for the whole set of transcripts is not surprising, because most genes were not host-specific in expression. However, the percentage of variation explained by host for each individual transcript ranged from 0% to 97.09%, indicating that there are individual transcripts with strong host-specific expression.
To reduce the influence of outlier samples within groups of small size (e.g., n = 3 for maize), we selected a subset of differentially expressed transcripts that were highly associated with host species (P-value < 0.01, log fold change > 2, variation explained by host of at least 60%). Using these criteria, we identified 38 transcripts, of which 33 were downregulated in grain sorghum and 5 upregulated in grain sorghum compared with maize (Figure 1). The P-values obtained in the DE analysis ranged from 5.68 × 10−9 to 0.0083, and the percentage of the variation explained by host ranged from 97.09 to 60.09 (Supplementary Table 5). Samples clustered based on host, indicating that host is the main factor associated with expression of these transcripts (Figure 1).
To further investigate the levels of these differentially expressed transcripts in other parts of the plant, we obtained expression data (transcripts per million [TPM]) for the following developmental stages in S. hermonthica from the PPGP II: imbibed seeds (stage 0), roots of germinated seedlings after exposure to GR24 (stage 1), roots of germinated seedlings after exposure to host roots (stage 2), expression in haustoria attached to host roots (prevascular connection) (stage 3), haustoria attached to host roots after vascular connection (stage 4), late postattachment from belowground plants (stage 5), stem and leaf tissue (stage 6.1), and flower tissue (stage 6.2). The transcripts from all of these stages were isolated from parasites growing on grain sorghum as the host. Six of the differentially expressed transcripts were absent in the libraries from the tissue-specific experiments (Figure 2). Five of those six transcripts were found to be upregulated in the maize samples but downregulated in the grain sorghum samples (Figure 1). Moreover, the tissue-specific expression data (from grain sorghum hosts) contained two transcripts highly expressed in all stages except in the aboveground tissues, and those were also found to be upregulated in the maize samples and downregulated in grain sorghum.
We annotated the 38 focal transcripts, finding 34 of them had a significant hit in the blastx search (Table 3; Supplementary Table 5). Once again, we focused on the functional gene categories and the GO terms extracted from the annotation to find important transcripts for the host–parasite interaction. One-third of the focal transcripts matched gene categories reported as enriched among core parasitism genes and formation of haustoria (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). The most common group in the data was transporters, followed by cell wall modification (Table 3). Most of the transporters belonged to carbohydrate and amino acid transporters and metabolic functional gene categories. The three most common GO terms (level 3) identified during the annotation process for the category biological processes were all related to metabolism: cellular metabolic process, organic substance metabolic process, and primary metabolic process (Supplementary Table 5). In the category molecular function, we found protein binding, organic cyclic compound binding, and heterocyclic compound binding (Supplementary Table 5). Finally, for the category cellular component, the three most common GO terms were intracellular part, intracellular organelle, and membrane-bounded organelle (Supplementary Table 5). The most common enzyme categories were hydrolases, followed by oxidoreductases, transferases, and lyases (Supplementary Table 5). Because host-specific responses may be linked to important physiological mechanisms, we identified pathways associated with the focal transcripts. We found 26 pathways that included nine enzymes (Supplementary Table 6). The most common pathways were thiamine metabolism and purine metabolism, consistent with their annotation in the amino acid and metabolic functional gene category. For those transcripts with no annotation, we searched the Rfam database but found no match. Finally, although minimum coverage criteria for the FST analysis precluded detection of differentially expressed transcripts among the set of high FST transcripts, the most common pathways detected for the transcripts containing the top 1% of highly differentiated SNPs were consistent with the ones found in the DE transcripts.
If the high gene flow in S. hermonthica detected by our allelic frequency analyses limits adaptation to specific hosts, plastic host-specific gene expression may be key to successful parasitism. So far, expression studies looking for host-specific interactions with root parasites have focused on the haustorium and host–root interface, overlooking parasite aboveground responses to host (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013; Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). Here we focused on expression of aboveground tissue in an attempt to identify a group of focal transcripts associated with host-specific, aboveground processes in this obligate hemiparasite. Of the total variation in expression, 9% was attributable to host. High genetic diversity among populations, even at limited spatial scales (Bozkurt et al. Reference Bozkurt, Muth, Parzies and Haussmann2015; Unachukwu et al. Reference Unachukwu, Menkir, Rabbi, Oluoch, Muranaka, Elzein, Odhiambo, Farombi and Gedil2017; Welsh and Mohamed Reference Welsh and Mohamed2011), together with temporal and environmental variation impossible to control in a natural setting, are both likely contribute to this variation. Still, we identified a set of transcripts exhibiting strong associations with host, supporting differential expression of a subset of genes involved in host-specific interactions and suggesting that aboveground organs also play a role in the parasitic process (Table 3).
Among differentially expressed focal transcripts, we detected several previously reported and relevant gene categories for host–parasite interactions (core parasitic genes) (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). Common functional categories detected in the analysis were genes encoding transporters, cell wall modification, oxidation–reduction process, and defense mechanisms (Table 3; Supplementary Table 5). The most abundant functional categories were amino acid and carbohydrate transport and metabolism (Table 3). Nutrient transport was also the main predicted functional gene category detected in transcripts containing FST outliers, suggesting that differences in nutrient transport in aboveground tissues may be a key characteristic of parasitism on distinct hosts. In the stem parasite field dodder (Cuscuta pentagona Engelm.), transcriptomic analyses revealed that after the establishment of vascular connections, amino acid and sugar transporters are upregulated (Ranjan et al. Reference Ranjan, Ichihashi, Farhi, Zumstein, Townsley, David-Schwartz and Sinha2014). Differences in expression in these transcripts between hosts suggests that as the parasitic process progresses, the host–parasite interaction relies increasingly on the transfer of nutrients and solutes from host to parasite. In agreement, Yang et al. (Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015) found transcripts for several amino acid and sugar transporters consistently upregulated across several Orobanchaceae parasites. Moreover, S. hermonthica infection is associated with increased amino acid levels in the host xylem sap (Pageau et al. Reference Pageau, Simier, Le Bizec, Robins and Fer2003). One strongly differentially expressed transcript between hosts (upregulated in sorghum) encoded an organic cation carnitine transmembrane transporter. Carnitine is associated with fatty acid metabolism in plants (Bourdin et al. Reference Bourdin, Adenier and Perrin2007), and it is present in many seeds, where it acts as a store of energy and has a critical role in postgermination growth (Lawand et al. Reference Lawand, Dorne, Long, Coupland, Mache and Carol2002; Panter and Mudd Reference Panter and Mudd1969). Another pivotal gene category for parasitic functions is cell wall modification (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013; Ranjan et al. Reference Ranjan, Ichihashi, Farhi, Zumstein, Townsley, David-Schwartz and Sinha2014; Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). To establish a successful host–parasite connection, parasites have to loosen or degrade the host cell wall and establish connections with the host vascular tissue (Lawand et al. Reference Lawand, Dorne, Long, Coupland, Mache and Carol2002; Saucet and Shirasu Reference Saucet and Shirasu2016). A cell wall modification gene of interest as a candidate for host specificity was a pectinesterase inhibitor. This transcript was strongly differentially expressed between hosts in the field data (upregulated in maize) and was absent in the tissue-specific expression data from laboratory-grown plants. Pectinases mediate the deposition of Casparian strips in the endodermis by recruiting the lignin polymerization machinery (Roppolo et al. Reference Roppolo, Boeckmann, Pfister, Boutet, Rubio, Dénervaud-Tendon, Vermeer, Gheyselinck, Xenarios and Geldner2014). Pectinesterase inhibitors have been linked to enhanced resistance to fungal infection in plants (Liu et al. Reference Liu, Sun, Pei, Zhang, Wang, Li, Li and Hou2018) and to symbiotic nodulation (Young et al. Reference Young, Debellé, Oldroyd, Geurts, Cannon, Udvardi, Benedito, Mayer, Gouzy, Schoof, Van de Peer, Proost, Cook, Meyers and Spannagl2011; Zouari et al. Reference Zouari, Salvioli, Chialva, Novero, Miozzi, Tenore, Bagnaresi and Bonfante2014).
Another functional gene category shown to be upregulated in parasitic plants is defense response and pathogenesis (Honaas et al. Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013; Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). Striga manipulates host immunity for its own benefit (Saucet and Shirasu Reference Saucet and Shirasu2016) and can alter ABA, gibberellic acid, and cytokinin concentrations during host infection (Frost et al. Reference Frost, Gurney, Press and Scholes1997). In addition, it has been suggested that Striga can produce a set of effectors that subdue host defense (Runo and Kuria Reference Runo and Kuria2018). In this study, we found several strongly differentially expressed transcripts associated with defense and pathogenesis. One of them was a gene encoding an auxin-signaling F-box that was upregulated in the maize samples but absent from the grain sorghum host tissue-specific expression data. Finally, another relevant category is oxidation–reduction process. The latter was found to be enriched in orthogroups, with greater dN/dS in parasitic lineages compared with nonparasitic lineages (Yang et al. Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015). The most strongly differentially expressed transcript between hosts belongs to this category, a phenylcoumaran benzylic ether reductase-like protein. This protein is associated with oxidoreductase activity and flavonoids and has been previously associated with witchweed [Striga asiatica (L.) Kuntze] haustorial formation (Liu Reference Liu2012). Certain flavonoids act as haustorium-inducing factors (Albrecht et al. Reference Albrecht, Yoder and Phillips1999; Riopel and Timko Reference Riopel, Timko and Verma1992), while oxidoreductase activity triggers haustorial formation (Bandaranayake et al. Reference Bandaranayake, Filappova, Tomilov, Tomilova, Jamison-McClung, Ngo, Inoue and Yoder2010). Finally, three transcripts had no hit in the blastx search, but we found a functional gene category in the orthologue search. This was also reported in Yang et al. (Reference Yang, Wafula, Honaas, Zhang, Das, Fernandez-Aparicio, Huang, Bandaranayake, Wu, Der, Clarke, Ralph, Landherr, Altman and Timko2015) and Honaas et al. (Reference Honaas, Wafula, Yang, Der, Wickett, Altman, Taylor, Yoder, Timko, Westwood and dePamphilis2013): several parasitic-specific transcripts had nonsignificant blast alignments or had a hit to genes encoding proteins of unknown function, suggesting that many parasite-specific transcripts have a yet to be discovered function, even when evidence is based on homology comparisons.
Most of the host-associated, differentially expressed transcripts were upregulated in maize compared with grain sorghum (33 vs. 5, respectively) (Table 3). The majority of the upregulated functional gene categories and pathways were associated with nutrient metabolism (i.e., sugars, amino acids) (Table 3; Supplementary Table 6). Because nutrient uptake from the host is essential for the parasite’s growth and reproduction, differences in physiological integration at the vegetative stage can be inferred by studying expression patterns of nutrient transport and metabolic transcripts. These differences in expression are likely to be caused by selection pressures imposed by a different growing environment (i.e., host crop). We hypothesize that our results might be associated with host tolerance level to the parasite infection. Tolerance to S. hermonthica varies greatly within and among crops (Amusan et al. Reference Amusan, Patrick, Menkir, Housley and Ejeta2008; Haussmann et al. Reference Haussmann, Hess, Omanya, Folkertsma, Reddy, Kayentao, Welz and Geiger2004; Hess et al. Reference Hess, Ejeta and Butler1992; Lane et al. Reference Lane, Child, Moore, Arnold and Bailey1997; Maiti et al. Reference Maiti, Ramaiah, Bisen and Chidley1984; Vogler et al. Reference Vogler, Ejeta and Butler1996). In this regard, tolerance in grain sorghum to infection with S. hermonthica has been reported in many cases, but tolerance in maize is rare (Amusan et al. Reference Amusan, Patrick, Menkir, Housley and Ejeta2008; Lane et al. Reference Lane, Child, Moore, Arnold and Bailey1997; Mutinda et al. Reference Mutinda, Masanga, Mutuku, Runo and Alakonya2018). This variation could be attributed to the evolutionary history of these crops. Grain sorghum is a native African crop that has presumably coevolved with S. hermonthica and other Striga species, while maize evolution took place in the absence of the parasite. The almost absent resistance in maize might benefit the host–parasite physiological integration, facilitating nutrient uptake from the host at the root level, with the nutrients then metabolized in the parasite’s aboveground tissue. Further work is required to test this hypothesis. Finally, because our data come from aboveground tissue, finding haustorium activity–related transcripts could indicate a global pattern of host-specificity expression (i.e., across all tissues). Alternatively, aboveground organs could be a source of gene expression after the parasite’s germination and establishment of the haustorial connection, contributing to the maintenance of the connection with the host. Overall, the DE analysis revealed a set of predicted functional gene categories with candidate transcripts relevant to host-specific interactions (Table 3). The importance of several of these candidate transcripts in host-specific interactions was further supported by tissue-specific expression data from the PPGP II. The tissue-specific data came from S. hermonthica plants parasitizing grain sorghum, and interestingly, almost all of these transcripts (five out of six) were upregulated in maize and downregulated in grain sorghum (Table 3). Nevertheless, we have to consider the possibility that the absence of these transcripts in the tissue-specific expression data might be explained by methodological causes (not captured in libraries or expressed at very low levels).
In summary, the agriculturally important parasitic weed S. hermonthica infects multiple host species, making it a valuable system to study the complex genetic mechanisms of parasite responses to diverse hosts. Parasite interactions with different hosts likely require and result in physiological changes throughout the parasitic plant, especially in obligate parasitic weeds such as S. hermonthica. In addition, as host recognition acts at several levels (i.e., germination, haustorial initiation, attachment and penetration, and physiological integration), the parasite response to the host is likely to vary throughout development. After attachment, obligate parasites become a metabolic sink for carbohydrates and nitrogen from the host (Musselman Reference Musselman1980). It therefore remains critical to understand variation in dynamics of the host–parasite interaction beyond the more well-studied initial stages of host detection, development, and attachment of the haustorium, and establishment of vascular connections. Host transcriptomic studies have previously provided insights into Striga genotype–host genotype/species interactions, but so far, no study has addressed how genetic polymorphism may contribute to plasticity and/or host-specific genotypes. Here, we leveraged transcriptomic data from aboveground tissue to identify a set of host-associated transcripts in natural populations of S. hermonthica and evaluate levels of genetic polymorphism in transcripts expressed among populations parasitizing different hosts, illuminating the role that expression variation and genetic polymorphism may play in mediating host–parasite dynamics during later stages of the interaction. Despite high gene flow among populations, we identified a set of host-associated transcripts using differential expression and population genomic analysis of expressed genes. We found several functional gene categories relevant to host–parasite interactions supported by previous studies in parasitic plants. More specifically, our results highlight the role of nutrient transporters and the probable importance of gene regulation in aboveground organs in the specificity of the parasitic process. Those parasite genotypes that are not removed by selection during germination and haustorial connection may be selected later by pressures that the host, as an environment, imposes on the parasite. In addition, we found several defense- and pathogenesis- related genes together with plant hormone-response genes, supporting the idea that pathways involved in plant defense/pathogenesis and hormone response can be used by the parasite during the parasitic process. Overall, we provide a set of candidate transcripts that show host-specific interaction in aboveground tissue in the parasitic plant S. hermonthica. This suggests that signals of host-specific processes can be detected aboveground, expanding the focus of host–parasite interactions beyond the haustorial connection.
Supplementary material
To view supplementary material for this article, please visit https://doi.org/10.1017/wsc.2019.20.
Author ORCID
Lua Lopez https://orcid.org/0000-0002-3619-763X.
Acknowledgments
No conflicts of interest have been declared. This work was supported by research funds from the International Institute of Tropical Agriculture (IITA) Nigeria to CWD and the National Science Foundation (NSF IOS-1238057) to CWD and MPT, and by start-up funds from Penn State’s Eberly College of Science to JRL. Emily S. Bellis was supported by an NSF Postdoctoral Research Fellowship in Biology (grant number 1711950). The S. hermonthica field sampling and sequencing was funded under a USAID linkage grant from IITA. We thank Joshua P. Der for helpful ideas and comments on the project design. We thank Rafael Guerrero for his advice with the allelic frequency analysis. Likewise, we thank two anonymous reviewers for their comments on the article.
