Analysis of microsatellites from the transcriptome of downy mildew pathogens and their application for characterization of Pseudoperonospora populations

In silico identification and analyses of SSRs in predicted downy mildew transcriptomes

Transcriptomes predicted from genome assemblies for two downy mildew pathogens, Pseudoperonospora cubensis and Hyaloperonospora arabidopsidis that are publicly available were used for SSR identification. The FASTA file for the P. cubensis genome sequences and the .gff3 file with predicted gene coordinates were downloaded from Savory et al. (2012a) and Bioperl was used to generate a predicted transcriptome (Fig. S1). The genome assembly of H. arabidopsidis is described in Baxter et al. (2010), and is located in EMBL/Genbank/DDBJ databases under accession GCA_000173235.2. The transcripts predicted from the genome assembly were downloaded from EnsmblProtists database under HyaAraEmoy2_2.0.

The Microsatellite Identification Tool (MISA) (Thiel et al., 2003) was used to search the P. cubensis and H. arabidopsidis transcriptomes for the presence of microsatellites. MISA reported perfect and compound microsatellites ranging from one to six base pair units with a specified minimum number of repeats for each motif, specifically, 1/20, 2/5, 3/4, 4/3, 5/3, 6/3 (unit size/minimum number of repeats). Comparisons between microsatellite abundance, frequency, and motif types in the two downy mildew transcriptomes were calculated with a proportion test as conducted in Garnica et al. (2006).

Primer design

Bioperl programing was used to parse the MISA output file so the program Primer3 (Rozen & Skaletsky, 2000) could use sequence coordinates reported by MISA to design primers flanking the identified microsatellites. The program designed primers that would amplify products between 100 and 300 bp. The primers were to be between 18 and 27 bp with optimum length of 20 bp, with GC content between 20% and 80% with optimum GC content of 50%, and with a melting temperature between 57 and 63 °C with the optimum melting temperature of 60 °C. Primers representative of the different motif groups (100 total) were ordered from and manufactured by Integrated DNA Technologies, Inc. (Coralville, IA, USA). Forward primers were designed to include an M13 tail for fluorescent labeling of products and later fragment analyses (Schuelke, 2000).

Tissue collection and DNA extraction

A total of 11 Pseudoperonospora cubensis isolates and two P. humuli isolates were used to evaluate 100 Primer3-generated microsatellite primers for amplification, transferability, and polymorphism via gel electrophoresis. For analysis with fragment analysis, 38 P. cubensis isolates and 22 P. humuli isolates were screened (Table 1). Cucurbit and hop leaves infected with downy mildew were collected from several locations throughout North Carolina (NC) in 2013 and 2014 (Lenoir, Haywood, Johnston, and Rowan counties). The presence of the pathogen on leaves was confirmed by the observation of sporulation using a dissecting microscope. Leaf lesions were excised by sterile scalpel, placed in a microfuge tube, and stored at −80 °C until the time of DNA extraction. The isolate used in the sequencing of the P. cubensis genome (MSU-1) was included as a positive control (Savory et al., 2012a; Savory et al., 2012b). Collaborators provided isolates from other geographic regions and DNA was extracted from pelleted sporangia (Table 1). Tissue was disrupted and DNA was extracted and purified via phenol-chloroform extractions adapted from previous work (Ahmed et al., 2009). DNA was purified with ethanol washes then suspended in 1× TE buffer and quantified using a NanoDrop ND 1000 spectrophotometer and NanoDrop 2.4.7c software (NanoDrop Technologies Inc., Wilmington, DE, USA). Integrity of the DNA was confirmed by gel electrophoresis with the presence of a >12,000 bp band. Each isolate was amplified in the nrITS and mitochondrial Nad1 and Nad5 (Quesada-Ocampo, Granke & Olsen, 2012) where presence of a band of appropriate size (700 bp, 500 bp, 300 bp, respectively) confirmed the presence of P. cubensis or P. humuli DNA in lesion tissue and a lack of a band confirmed negative controls (uninfected cucumber leaf tissue and water control).

Primer evaluation with gel electrophoresis

The microsatellite primers report from the Primer3 output was divided by microsatellite motif types (tri-, tetra-, penta-, hexa-nucleotide repeats) and arranged from highest number of repeating units in descending order, as microsatellites with more repeats tend to be more polymorphic (Guichoux et al., 2011; Ellegren, 2000). The top fifteen primers of each motif type were selected to be validated for amplification and tested for consistency and transferability to P. humuli in a screen against the eleven Pseudoperonospora cubensis isolates from different cucurbit hosts (Cucumis sativus, Cucumis melo, Cucurbita pepo, Cucurbita maxima, Cucurbita moschata, Citrullus lanatus, Momordica balsamina, and Momordica charantia) and two P. humuli isolates (Table 1). Forward primers were labeled with a partial M13 tail (GACGGCCAGT) on the 5′ end so they could be used in fragment analysis downstream. Preliminary results suggested microsatellites with certain motifs were more likely to be polymorphic, so more markers with tri- and hexa-nucleotide repeats were evaluated. PCR reactions were performed with 10 µL of 2xGoTaq^® Hot Start Green Master Mix (Promega, Madison, WI, USA), 1 µL of 10 µM forward primer, 1 µL of 10 µM reverse primer, 10 ng of DNA, and sterile water. The thermal cycler (Bio-Rad, Hercules, CA) program, CDMSSR1, was set to have an initial denature of 94 °C for 3 m, and 35 repeating cycles consisting of denaturing at 94 °C, annealing at 53 °C, and extension at 72 °C, each step having a duration of 30 s. This program concluded with a 5 m final extension at 72 °C.

P. cubensis isolates from watermelon, bitter melon, and balsam apple, and the positive and negative controls were amplified by touchdown PCR (TDSSR1) because infected leaves with low levels of sporulation did not amplify reliably with standard PCR settings. These isolates underwent PCR with thermal cycler settings with an initial denaturing step of 94 °C for 5 m, then 20 cycles of a 30 s denaturing step at 94 °C, a 45 s at an annealing temperature starting at 62 °C, and an extension step at 72 °C for 2 m. At each cycle, the annealing temperature would decrease by 0.5 °C. Then the reaction continued with another twenty cycles of a 30 s denature step at 94 °C, a 45 s annealing step at 55 °C, and a 2 m extension step at 72 °C. The reaction ended with a 5 m final extension at 72 °C (Korbie & Mattick, 2008). These 13 Pseudoperonospora isolates were screened against a total of 100 SSR markers and evaluated via gel electrophoresis. PCR products were run on 4% ultrapure agarose gels at 40 volts for approximately five hours to evaluate amplification and determine product size.

Fragment analysis of polymorphic SSR primers

From the initial screening of 100 microsatellite markers with agarose gel electrophoresis, a subset of 17 primers with consistent amplification across samples and appearance of polymorphism across P. cubensis isolates were further analyzed. These 17 primers were applied to 38 P. cubensis isolates and 22 P. humuli isolates and evaluated for polymorphism and consistency via fragment analysis (Table 1). PCR products of downy mildew isolates amplified with the polymorphic SSR primers were subjected to a second PCR, CDMSSR2, in order to attach fluorescent dyes to the amplified products. Reactions were carried out in 10 µL volume consisting of the same reagents and concentrations as above with the exception of an M13 primer (TGTAAAACGACGGCCAGT) tagged with a fluorescent dye in place of the site-specific forward primer (Schuelke, 2000). Thermal cycler settings for the program CDMSSR2 were the same as the CDMSSR1 program, but with 15 repeating cycles of denaturing, annealing, and extension steps. Products from the CDMSSR2 were diluted fifty to twenty-five-fold and pool-plexed, combining multiple PCR products of different microsatellite primers labeled with different fluorescent dyes (VIC and 6FAM) (Applied Biosystems, Foster City, CA, USA). A genotyping reaction was performed where HiDi Formamide, LIZ600 size standard, and the diluted, pool-plexed sample were combined then submitted to the NCSU Genomic Science Laboratory (GSL, Raleigh, NC, USA) for genotyping with a 3730xl DNA Analyzer (Applied Biosystems, Foster City, CA, USA).

Results were individually analyzed and binned with the microsatellite plugin for the program Geneious version 8.1 (Kearse et al., 2012). On a given isolate fragment analysis run, peaks occurring one length of the repeat motif away from the peak with the highest signal and were less than 15% of the height of the larger peak were removed to decrease the risk of genotyping stutter peaks. Two alleles were assumed to be present at each locus because Pseudoperonospora spp., belonging to Oomycota, are diploid organisms. If one peak was observed at any given locus, homozygosity was assumed. Six of the 17 makers evaluated with fragment analysis proved to be monomorphic across the evaluated P. cubensis isolates and were removed from further analysis.

Data analysis

Basic summary statistics were calculated for 11 reliable polymorphic primers across 38 P. cubensis isolates from diverse hosts, years, and locations, and 22 P. humuli isolates from diverse locations (Table 1). The R package Poppr version 2.2.0 was used to calculate descriptive population statistics such as the genotype accumulation curve, heterozygosity, and evenness, through the “gac” and “poppr” function (Nei, 1978; Pielou, 1975; Grünwald et al., 2003; Kamvar & Grunwald, 2014). To assess genotypic richness while accounting for sample size, a rarefaction curve was generated with the “vegan” package and “rarecurve” function. The dataset was clone-corrected with the “clonecorrect” function, and then, to determine if P. cubensis and P. humuli were sexual populations, the index of association was calculated with the “ia” function (Agapow & Burt, 2001). To evaluate the genetic distance between individuals, an UPGMA dendrogram based on Bruvo’s distance was created with 1,000 bootstrap replicates (R Core Team, 2016; Kamvar & Grunwald, 2014; Kamvar, Brooks & Grünwald, 2015).

In silico identification and analysis of SSRs in predicted downy mildew transcriptomes

MISA analysis revealed that of the 23,522 P. cubensis sequences examined, 2,398 sequences contained microsatellites, with a total of 2,738 microsatellites identified. In H. arabidopsidis, 14,548 sequences were examined and 1,691 of the examined sequences contained microsatellites. A total of 2,119 microsatellites were identified in the H. arabidopsidis transcriptome. The number of SSR-containing sequences out of the total number of sequences examined were significantly different between species (p < 0.0002). A significant difference was also found between the percentage of total microsatellites identified out of the total number of sequences examined in P. cubensis (11.6%) and H. arabidopsidis (14.5%) (P < 0.0002). The total relative abundances of microsatellites in P. cubensis and H. arabidopsidis (101.79 and 152.01, respectively), show a greater difference between organisms compared to the percentage of SSR-containing sequences of the total number of sequences examined (10.2% and 11.6%, respectively). The total relative abundance (SSR/Mb) of microsatellites in P. cubensis and H. arabidopsidis revealed less of a difference compared to the total relative densities (bp/Mb) between the two species. H. arabidopsidis has a much higher density of microsatellites in the predicted transcriptome (2,322.89) compared to P. cubensis (1,421.66) (Table 2).

A majority of the identified microsatellites in both P. cubensis and H. arabidopsidis sequences were perfect (96% and 93%, respectively), meaning there are no interrupting sequences within the chain of repeating units. Also, 1.2% of the sequences examined from the P. cubensis transcriptome contained more than one microsatellite, whereas 2.2% of the H. arabidopsidis sequences examined contained more than one microsatellite.

MISA analysis indicated that tri-nucleotide repeats were the most frequently occurring motif found in both transcriptomes (Table 3). This motif-type made up 61% of the total microsatellites from the P. cubensis transcriptome and 71% of the total microsatellites from the H. arabidopsidis transcriptome (Table 3). In both species, di-nucleotide repeats were the second most frequently occurring motif, followed by tetra-, hexa-, and penta-nucleotide repeats. There were significant differences found between the percentage of di- and tri-nucleotide repeats between species (P < 0.0002), but there were no significant differences found between the percentage of tetra-, penta-, and hexa-nucleotide microsatellites (P > 0.1). In terms of relative abundance, the same trend held with tri-nucleotide repeats having the highest value in both species, followed by di-, tetra-, hexa-, and then penta-nucleotide repeats. Although the relative density values for each motif group in P. cubensis also kept this trend, tetra-nucleotide repeats had a higher density than di-, hexa-, and penta-nucleotide repeats in H. arabidopsidis. Also, in H. arabidopsidis higher relative abundance and relative density values were seen for each motif group except di-nucleotide repeats, in which P. cubensis held higher values (Table 3). It should also be noted that in the H. arabidopsidis transcriptome, there was only one monomer identified that met the specification stated in the MISA script of 20 repeating units. There were no monomers identified in the P. cubensis transcriptome.

Both transcriptomes had the same repeating sequences occur with high frequency. For example, the motif AGC/CTG was the most commonly occurring motif in both transcriptomes, being the repeating sequence in 451 of the microsatellites identified in the in the P. cubensis transcriptome and 495 of the microsatellites identified in the H. arabidopsidis transcriptome. Eleven out of the fifteen most common motif sequences were the same for both P. cubensis and H. arabidopsidis (Table 4). A majority of the microsatellites identified in each motif-type group fell toward the lower bound of the repeat range set for microsatellite identification using MISA (Fig. 1). MISA was to identify repeating sequences that exceeded five repeating motifs for di- nucleotide repeats, four repeating motifs for tri-nucleotide repeats, and three repeating units for tetra-, penta-, and hexa-nucleotide repeats. In both species, over 90% of tetra-, penta-, and hexa-nucleotide repeat microsatellites were made up of three repeating units. A higher percentage of longer chains of repeating units occurred in di- and tri- nucleotide repeats. For example, tri-nucleotide repeat microsatellites with more than four repeating units made up 19% of the tri-nucleotide repeat microsatellites in P. cubensis and approximately 23% of the tri-nucleotide repeat microsatellites in H. arabidopsidis (Fig. 1).

Primer design and lab validation

Primers were successfully designed by Primer3 to amplify 2,088 microsatellites out of the 2,738 microsatellites identified by MISA in the transcriptome of P. cubensis. A majority of these primers (97%) were for amplification of perfect microsatellites, meaning there are no sequences interrupting the repeating motif; however, 3% of the primers are predicted to amplify compound microsatellites, where non-SSR base pairs may be found within the repeating motif sequence. Of the 2,088 primers designed, 417 were di-nucleotide repeats, 1,125 were tri-nucleotide repeats, 248 were tetra-nucleotide repeats, 40 were penta-nucleotide repeats, and 91 were hexa-nucleotide repeats.

Experimental validation of 100 primer sets with gel electrophoresis (Table S1) revealed that 92% of the selected markers produced a product across the P. cubensis isolates used, and 90% of the 92 primer pairs that produced a product were the size predicted by Primer3 (83 total). The electrophoresis results also revealed a majority of the primers showed significant species transferability. Of the 85 primers that produced a product in more than just the positive control isolate (MSU-1), only one primer set, SSR94, did not produce a product in the P. humuli isolates screened. Of the 92 primers that produced a PCR product of the Primer3 prediction, seven primers only amplified the isolate that was used to sequence the P. cubensis genome (MSU-1) (Savory et al., 2012a; Savory et al., 2012b;). Of the initial 100 primer sets selected to be experimentally validated with gel electrophoresis, 17 primer sets had consistent amplification and appeared to be polymorphic across P. cubensis isolates. When the 17 primers were applied to a larger panel of P. cubensis isolates (n = 38), only 11 primer sets were identified to consistently amplify loci with polymorphic alleles within P. cubensis (Table 5). These 11 markers were determined to be polymorphic via the “informloci” function in poppr. When the loci were sampled 1,000 times without replacement, the genotype accumulation curve demonstrated that 7–8 markers were necessary to discriminate between 90% of the multilocus genotypes (Fig. 2). Five of these 11 polymorphic primers were tri-nucleotide repeats and four were hexa-nucleotide repeats, initially selected based on predicted number of repeats. The final two polymorphic primers were selected because they were predicted to be located in function-associated genes. These two markers, SSR97 and SSR92, were identified in putative Crinkler family proteins.

Marker characterization and summary statistics

When these 11 primers were applied to a diverse panel of P. cubensis isolates from all major commercial hosts and three non-commercial cucurbits spanning several years and geographic locations (n = 38) and P. humuli isolates from diverse geographic regions (n = 22), descriptive population statistics could be determined. Over eleven loci, 89 alleles were found in the isolates evaluated. Overall, more alleles were found at each locus across P. cubensis isolates, with the exception of SSR79, in which two alleles were found across the set of isolates for each species. The number of alleles per loci ranged from two to nine, with an average allele diversity of 5.18 alleles per locus in P. cubensis and 2.91 alleles per locus in P. humuli. Heterozygosity across all isolates had a mean value of 0.45, with the mean heterozygosity being 0.53 in P. cubensis and 0.38 in P. humuli. This higher genotypic richness seen in P. cubensis isolates was further confirmed by generating a rarefaction curve to account for differences in sample size (Fig. 3). The evenness of alleles at each locus ranged from 0.38 to 0.98, with the mean evenness of 0.69 and 0.70 in P. cubensis and P. humuli, respectively, showing a similar moderate distribution of multi locus genotypes (MLGs) in both species (Table 5).

The isolates were also tested to determine if either species were in linkage disequilibrium. In P. cubensis isolates, the ${\overline{r}}_{\mathrm{D}}$ value was 0.04, which lies within the distribution expected under linkage. The p value of 0.116 does not reject the null hypothesis that the alleles seen across the 11 loci are not linked, suggesting the P. cubensis isolates are sexually reproducing. On the other hand, the ${\overline{r}}_{\mathrm{D}}$ value for the P. humuli isolates was 0.28, which falls outside of the distribution expected under linkage. The p value of 0.001 supports the rejection of the null hypothesis that the alleles seen across the 11 loci are not linked, suggesting the P. humuli isolates are clonal (Fig. 4).

Distinct differences were seen between P. cubensis and P. humuli isolates when an UPGMA dendrogram was generated for Bruvo’s distance between individuals. Two clear groups could be seen based on comparison between individuals without taking population into consideration and is supported with 100% confidence, as the analysis had 1,000 bootstrap replications. Cluster 1 as seen in Fig. 5 includes only P. cubensis isolates. P. humuli isolates, on the other hand, group together in Cluster 2. There is some support for distinct clusters within the P. humuli group, the smaller of which contain three P. humuli isolates, Cas, Sant2-5, and Hdm501ba. These isolates are from North Carolina and Oregon, but so are several P. humuli isolates that fall into the larger P. humuli cluster. Cluster1 comprises of P. cubensis isolates, with no clear trend or support of sub-clustering by host or location. However, it is clear the diversity within P. cubensis isolates is much more varied and complex compared to that of P. humuli. Trends may be become apparent when a larger dataset with more isolates per host and location is evaluated.