Transcriptome-enabled marker discovery and mapping of plastochron-related genes in Petunia spp.
- Yufang Guo†1,
- Krystle E. Wiegert-Rininger†1,
- Veronica A. Vallejo1,
- Cornelius S. Barry1 and
- Ryan M. Warner1Email author
© Guo et al. 2015
Received: 5 May 2015
Accepted: 16 September 2015
Published: 24 September 2015
Petunia (Petunia × hybrida), derived from a hybrid between P. axillaris and P. integrifolia, is one of the most economically important bedding plant crops and Petunia spp. serve as model systems for investigating the mechanisms underlying diverse mating systems and pollination syndromes. In addition, we have previously described genetic variation and quantitative trait loci (QTL) related to petunia development rate and morphology, which represent important breeding targets for the floriculture industry to improve crop production and performance. Despite the importance of petunia as a crop, the floriculture industry has been slow to adopt marker assisted selection to facilitate breeding strategies and there remains a limited availability of sequences and molecular markers from the genus compared to other economically important members of the Solanaceae family such as tomato, potato and pepper.
Here we report the de novo assembly, annotation and characterization of transcriptomes from P. axillaris, P. exserta and P. integrifolia. Each transcriptome assembly was derived from five tissue libraries (callus, 3-week old seedlings, shoot apices, flowers of mixed developmental stages, and trichomes). A total of 74,573, 54,913, and 104,739 assembled transcripts were recovered from P. axillaris, P. exserta and P. integrifolia, respectively and following removal of multiple isoforms, 32,994 P. axillaris, 30,225 P. exserta, and 33,540 P. integrifolia high quality representative transcripts were extracted for annotation and expression analysis. The transcriptome data was mined for single nucleotide polymorphisms (SNP) and simple sequence repeat (SSR) markers, yielding 89,007 high quality SNPs and 2949 SSRs, respectively. 15,701 SNPs were computationally converted into user-friendly cleaved amplified polymorphic sequence (CAPS) markers and a subset of SNP and CAPS markers were experimentally verified. CAPS markers developed from plastochron-related homologous transcripts from P. axillaris were mapped in an interspecific Petunia population and evaluated for co-localization with QTL for development rate.
The high quality of the three Petunia spp. transcriptomes coupled with the utility of the SNP data will serve as a resource for further exploration of genetic diversity within the genus and will facilitate efforts to develop genetic and physical maps to aid the identification of QTL associated with traits of interest.
KeywordsPetunia Transcriptome sequencing Molecular markers Single nucleotide polymorphisms Development rate Floriculture
The genus Petunia resides within the Solanaceae family and contains 20 species and subspecies that are native to South America . Petunia × hybrida (petunia) is an important ornamental crop plant and represents a hybrid species derived in the nineteenth century from a cross between P. axillaris and P. integrifolia . Subsequent breeding has introgressed traits from additional Petunia spp. that, together with natural variation resulting from mutations in key genes, have contributed to the wide diversity in plant and floral morphology and flower color that exists within the pool of commercially available germplasm [2–7]. In cool climates in the Northern Hemisphere, petunia is often produced in greenhouses during the winter months for distribution to spring markets once it reaches an optimal size and begins to flower [8, 9]. Therefore, a high percentage of the cost of crop production is related to energy consumption and growers are often faced with the dilemma of either reducing greenhouse temperatures, thereby extending the growing time of the crop and incurring increased labor costs, or elevating the growing temperature and increasing energy costs but reducing the duration of crop growth [8, 9]. Thus, understanding the factors that impact crop timing traits may facilitate the selection of petunia varieties with an increased rate of vegetative node formation (development rate) at either optimal or suboptimal growing temperatures, or varieties that initiate flowering following emergence of fewer leaf nodes. We have previously documented that accessions of P. axillaris and P. integrifolia possess increased development rate when compared to a diverse pool of commercial petunia germplasm, suggesting genetic variation for this trait within the genus . This was confirmed in an interspecific F2 population of a cross between P. axillaris and P. integrifolia that identified three quantitative trait loci (QTL) on chromosomes 1, 2 and 5 that affected development rate and explained 34 % of the observed variation . The molecular basis underlying these QTL remains to be identified.
The genetic determinants of development rate, often referred to as plastochron, are multifaceted, complex and not fully understood but are, at least in part, linked to hormonal control of meristem size and activity. For example, transgenic tobacco (Nicotiana tabacum) plants with increased cytokinin oxidase activity and a concomitant reduction in cytokinin levels displayed reduced meristem size and delayed plastochron when compared to wild type . Similarly, characterization of Arabidopsis mutants with reduced auxin levels and disrupted auxin transport also influence plastochron [13–15]. In addition, mutations at the plastochron 1 and plastochron 2 loci of rice, which encode a cytochrome P450 of unknown function and a MEI2-like RNA binding protein homolog, respectively, also influence development rate but do so independently of each other [16, 17]. In Arabidopsis, the SQUAMOSA PROMOTER BINDING PROTEIN-LIKE (SPL) transcription factors SPL9 and SPL15 act redundantly to influence plastochron and over-expression of miR156, which targets multiple SPLs, shortens plastochron . Analogous to the relationship of the plastochron 1 and 2 loci of rice , the SPL/miR156 regulatory module acts independently of CYP78A5/KLUH, which is a putative ortholog of rice plastochron 1 . The involvement of the miRNA pathway in influencing plastochron is further supported by the characterization of the serrate and altered meristem program (amp1) mutants of Arabidopsis, which display reduced and increased rates of leaf initiation, respectively [19–22]. SERRATE encodes a zinc finger protein required for miRNA biogenesis and RNA splicing while AMP1 associates with ARGONAUTE1 at the endoplasmic reticulum and is required for translation inhibition through the exclusion of miRNA target mRNAs from polysomes [19, 20, 23]. Furthermore, mutations in AMP1 homologs in maize and rice confer similar pleiotropic phenotypes to those exhibited by Arabidopsis amp1 mutants, including altered plastochron [24, 25]. Together, these data suggest complex regulation of plastochron that involves different regulatory modules, including hormone and miRNA pathways.
The development of next generation sequencing technology has revolutionized biology and in particular, transcriptome sequencing provides a cost effective strategy for generating sequence and expression information from the gene space of non-model organisms or from species with large complex genomes [26, 27]. In plants, transcriptome sequencing has facilitated gene discovery, the development of molecular markers and large scale analyses of genetic variation [28–33]. Despite the economic and biological importance of petunia, genomic information and molecular marker resources for this genus are limited [11, 34–38], single nucleotide polymorphism (SNP) markers are currently unavailable and marker assisted selection is rarely utilized. In addition, although transcriptome resources are available for petunia, they are not extensive and most often are derived from the cultivated species P. hybrida or from specialized tissue types [39–41]. Herein, de novo assembly of reference transcriptomes of P. axillaris, P. integrifolia and P. exserta, derived from paired-end RNAseq analysis of five diverse tissues (callus, seedling, shoot apices, flowers and trichomes) is reported. Tissue types were selected to attempt to maximize the number of transcripts recovered while generating resources for studying traits of interest related to development and metabolism. These resources were utilized to develop a set of SNP, cleaved amplified polymorphic sequence (CAPS) and simple sequence repeat (SSR) markers that will facilitate QTL mapping and gene discovery for multiple traits within the genus, including those associated with development rate.
Results and discussion
Transcriptome assembly and annotation
Description of Petunia spp. tissues, libraries, and RNA-Seq data
Number of raw reads (millions)
Number of filtered reads (millions)a
54.6 (95.8 %)
52.3 (96.2 %)
67.9 (96.3 %)
62.8 (96.0 %)
44.6 (94.5 %)
34.8 (94.5 %)
63.0 (94.4 %)
53.6 (94.7 %)
41.9 (95.3 %)
41.6 (94.5 %)
54.7 (94.5 %)
61.2 (95.0 %)
58.2 (94.8 %)
47.5 (94.7 %)
45.7 (95.4 %)
Total number of de novo assembled transcripts in each of the Petunia axillaris, P. exserta, and P. integrifolia transcriptomes in relation to the number of representative transcripts and their total coverage of the transcriptome space
Representative transcripts (no.)
Assembled transcripts (no.)
Transcriptome coverage (bp)
Representative transcripts from each assembly were annotated using BLASTX searches against five publically available databases (Additional file 2: Dataset S1). Adopting quality thresholds of ≥30 % coverage and ≥70 % identity resulted in an annotation rate of approximately 70 % (Fig. 1c), with the largest number of annotations retrieved from the RefSeq and NCBI non-redundant databases. Similar rates of annotation were reported for transcriptome assemblies of chickpea and red clover [48, 49]. In addition, of the 1944 predicted Petunia proteins available in GenBank, ~89 % are present in the P. axillaris, P. exserta, and P. integrifolia transcriptome assemblies.
SNP detection and characterization
Summary of read depth coverage for Petunia axillaris, P. exserta, and P. integrifolia reads mapped to the P. axillaris transcriptome assembly after de-duplication
Granular Quartile Read Depths
% of Reference Transcriptome with Read Depth Coverage:
Total reads (millions)
Mean read depth
Third (25 %)
Median (50 %)
First (75 %)
Overall, there were 814,408 SNPs between the three species. Of these, 105,645 were located in 5ʹ untranslated regions (UTRs), 481,289 were located within the coding sequence (CDS), and 217,178 were located in 3ʹUTRs. SNP frequency was calculated by dividing the total length of the reference transcriptome by the total number of SNPs. When only considering the length of 5ʹUTRs, CDS, and 3ʹUTRs in transcripts (not the length of the genomic regions these transcripts spanned), the SNP frequencies were 1 SNP/69 bp, 1 SNP/61 bp and 1 SNP/62 bp in these regions, respectively. The overall SNP frequency was 1 SNP/63 bp.
KAAS (KEGG Automatic Annotation Server) analysis of super pathways involving SNP-containing transcripts
Number of genes
Number of SNPs
Amino acid metabolism
Biosynthesis of other secondary metabolites
Glycan biosynthesis and metabolism
Metabolism of cofactors and vitamins
Metabolism of other amino acids
Metabolism of terpenoids and polyketides
Xenobiotics biodegradation and metabolism
Genetic Information Processing
Folding, sorting and degradation
Replication and repair
Environmental Information Processing
Signaling molecules and interaction
Cell growth and death
Transport and catabolism
SNPs were further classified based on their location; CDS or UTR and their zygosity. Altogether, only 10 homozygous SNPs were polymorphic across all three species combinations (Fig. 4b). Between P. axillaris and P. exserta, 24,528 SNPs were homozygous, which comprised 94.9 % of all SNPs. There were 13,512 homozygous SNPs, 17 % of the total SNP positions, between P. exserta and P. integrifolia. Between P. integrifolia and P. axillaris, there were only 4539 homozygous SNPs, which constituted only 6.2 % of all polymorphisms, indicating that the self-incompatible P. integrifolia is highly heterozygous. Of all SNP loci, 60,701 were located within the CDS, of which 60,561 were biallelic; 9537 were located within 5ʹUTRs, of which 9512 were biallelic; 17,983 were located within 3ʹUTRs, of which 17,923 were biallelic. There is a small group of 786 SNPs whose location could not be determined.
SNP frequency, the number of unigenes which contain SNPs, and percentage of the total unigene length for all three species combinations
SNP-containing unigenes (no.)
Total unigene length containing SNPs (%)
Unigenes with ≥10 SNPs (no.)
SNP-containing unigenes with ≥10 SNPs (%)
Maximum SNP frequency per unigene
P. axillaris and P. exserta
P. axillaris and P. integrifolia
P. exserta and P. integrifolia
SNPs between P. axillaris and P. exserta were distributed across 12,060 unigenes (Table 5). There were only 110 unigenes with 10 or more SNPs. The maximum SNP frequency per unigene was 1 SNP/111 bp. Between P. axillaris and P. integrifolia, SNPs were identified in 17,949 unigenes. There were 1466 unigenes with 10 or more SNPs. The maximum SNP frequency per unigene was 1 SNP/96 bp. There were 18,032 unigenes with SNPs between P. integrifolia and P. exserta, and 1722 unigenes with more than 10 SNPs. The maximum SNP frequency per unigene was 1/89 bp.
Due to the stringent nature of the SNP discovery parameters employed (i.e., filtering out three SNPs occurring within 100 bp of each other and filtering out SNPs within the first and last of 30 bp of a transcript, etc.), the above SNP frequencies calculated after the filtering steps are likely underestimated. Additionally, for organisms without a fully sequenced genome, using a transcriptome de novo assembly followed by variant detection can result in underestimation of expressed variants . Even with the collection of different tissue types in our study, there are still likely to be undetected variants.
SNP distribution and corresponding transition to transversion ratio (Ts/Tv) in coding regions (CDS), 5′UTRs and 3′UTRs
All SNPs (biallelic)
SNPs in CDs (biallelic)
SNPs in 5’UTRs (biallelic)
SNPs in 3′UTRs (biallelic)
P. axillaris and P. exserta (biallelic)
SNPs in CDs (biallelic)
SNPs in 5′UTRs (biallelic)
SNPs in 3′UTRs (biallelic)
P. integrifolia and P. axillaris (biallelic)
SNPs in CDs (biallelic)
SNPs in 5′UTRs (biallelic)
SNPs in 3′UTRs (biallelic)
P. exserta and P. integrifolia (biallelic)
SNPs in CDs (biallelic)
SNPs in 5′UTRs (biallelic)
SNPs in 3′UTRs (biallelic)
A total of 55 primer pairs were selected for PCR amplification. Of these, 50 amplified single DNA fragments, containing a total of 51 SNP loci, matching the expected sizes in P. axillaris and were selected for further sequencing (Additional file 1: Table S1). The failed/unclear amplification might be from primers located across intron-exon boundaries, or from the amplification being across a large intron or in paralogous genes. Three of the loci were monomorphic, five yielded poor sequence quality, which was associated with either failure of the sequence reaction or ambiguous base calls that were likely due to heterozygosity of the parents that resulted in amplification of distinct alleles of a given locus. The remaining 43 loci were polymorphic between at least two species resulting in an overall SNP validation rate of 84.3 % (43/51). For the three monomorphic SNP loci, the minor allele sequence read count frequencies were 0.07, 0.26 and 0.33, indicating that low minor allele sequence read count frequency might be one factor that can introduce false positives in SNP identification. The validation rate obtained in this study is comparable with previously reported results for transcriptome derived SNPs. For example, the validation rate for pea (Pisum sativum) transcriptome derived SNPs was 84.5 % by selecting a 1920 SNPs set in a GoldenGate VeraCode assay on a large number of pea accessions . The validation rate for the SNP detected from Douglas-fir (Pseudotsuga menziesii) was 72.5 %  while, in potato, the validation rate was 85 % for set of 96 transcriptome derived SNPs .
CAPS marker design and validation
Twenty restriction enzymes were used to computationally predict CAPS markers from 15,701 SNP loci (Additional file 5: Dataset S4). Some SNP loci are predicted to be digested by two different restriction enzymes, with one recognizing the reference allele and the other recognizing the alternative allele, or by different restriction enzymes with overlapping recognition sequences. For example, for the restriction enzymes AluI and SacI, SacI is a 6 bp cutter recognizing/cutting sequences of GAGCTˇC and AluI is a 4 bp cutter that recognizes/cuts AGˇCT, thus all AluI recognition cites are also recognized by SacI. These analyses led to the identification of a total of 17,635 predicted CAPS markers (Supplemental dataset 4). Fourteen putative CAPS markers were randomly selected for validation. All primer pairs amplified products of the expected size, and 11 were polymorphic. All 14 amplicons for P. axillaris and P. exserta had the expected genotype, while four amplicons from P. integrifolia were homozygous while the genotype prediction was heterozygous. However, the CAPS markers generated from this study have some limitations. First, the CAPS marker design pipeline only considered possible restriction enzyme recognition differences between species, the number/size of the fragments generated from digestion was not taken into consideration. Thus, there might be multiple small fragments yielded from the digestion, or the digested fragments between genotypes may only vary by a few base pairs. Either situation could not be scored successfully on agarose gel genotyping systems. Second, there could be sequence differences between the P. axillaris reference genome and the actual parental genotype used for this study (a different P. axillaris accession), and some restriction enzyme recognition sites predicted in the pipeline might not exist when screening our samples.
Simple sequence repeat (SSR) marker development
Summary of SSR markers in the P. axillaris transcriptome
Total number of sequences examined
Total size of examined sequences (bp)
Total number of identified SSRs
Number of SSR containing sequences
Number of sequences containing more than 1 SSR
Number of SSRs present in compound formation
Summary of SSR repeat motif types and their corresponding repeat unit numbers for di- and tri- nucleotide repeats
A total of 2949 primer pairs were developed from 3027 unigenes (Additional file 6: Dataset S5). The designed SSR primer pairs included 1,481 di-, 1,867 tri-, 88 tetra-, 16 penta-, and 47 hexa- nucleotide repeats (Table 7). Forty-eight SSR markers were randomly selected for their utility in petunia (Additional file 7: Dataset S6). Of these, 36 (75 %) were successfully amplified and displayed polymorphisms between at least two of the three species examined. For example, 24 SSRs were polymorphic between P. integrifolia and P. axillaris, 28 were polymorphic between P. axillaris and P. exserta, and 24 were polymorphic between P. integrifolia and P. exserta. This success rate was higher than for our previous results developing SSR markers from P. axillaris ESTs .
Representation and in silico expression analysis of Petunia plastochron-related transcripts
Petunia homologs of plastochron-related genes from Arabidopsis
Mapping plastochron-related homologous transcripts in an interspecific Petunia population
We identified Petunia transcripts homologous to numerous plastochron-related genes from Arabidopsis (Table 9). With the previously identified SNPs, 13 of these transcripts were converted to CAPS markers and used for linkage mapping (Additional file 1: Table S2). Together with the previously published SSR and CAPS markers , we constructed a linkage map with a total length of 289 cM consisting of 90 markers. Of those, eight were CAPS markers developed from plastochron-related transcripts located on five chromosomes (1, 2, 3, 6, and 7). The average linkage group length was 41.3 cM with a range from 32.7 cM (Chr7) to 58.7 cM (Chr5). The average marker density was 3.2 cM.
Summary of QTL for development rate (DRate; in nodes day−1) in a P. integrifolia × P. axillaris F2 population
High quality transcriptome assemblies were generated by RNA-seq data from three Petunia spp. and these were mined for molecular markers and plastochron-related transcripts. We used a SNP discovery pipeline to identify over 80,000 SNPs across three species. As the majority of the SNPs were located in exons and over 80 % of the SNP-containing transcripts can be annotated, these SNPs could possess utility in future gene discovery and breeding applications. Thirteen percent of the identified SNPs were polymorphic among all three species, and can be further used for comparative genomics and population diversity studies. Additionally, over 10,000 SNPs were transformed into CAPS markers and an additional set of SSR markers were developed. Plastochron-related transcripts were identified, converted to CAPS markers, and genetically mapped in a P. integrifolia × P. axillaris F2 population. Overall, these data provide sequence, expression information and a large marker resource for the Petunia community.
Plant material and growth conditions
Seeds of Petunia axillaris (PI667515) and Petunia exserta (OPGC943) were obtained from the Germplasm Resources Information Network. Seeds of Petunia integrifolia were purchased from Diane’s Flower Seeds (Ogdon, UT). Five different tissue types were harvested for RNAseq: 1) whole seedlings, 2) callus, 3) shoot apices, 4) flowers of mixed developmental stages, and 5) trichomes. For whole seedlings, P. axillaris, P. integrifolia and P. exserta seeds were sown on 100 mm diameter Petri plates containing half-strength Murashige and Skoog (MS) media. Plates were left unsealed to allow air exchange and grown under a 16-h photoperiod provided by fluorescent lamps (ca. 75 μmol m−2 s−1) at 22 °C. The date when seeds began to germinate was recorded (generally 2–3 d after sowing) and seedlings where harvested 21 d after this date. At harvest, media was removed from the roots, and plants were placed in 1.5 ml tubes and immediately frozen with liquid nitrogen.
For callus generation, seeds of each species were sterilized by: 1) soaking in sterilized water for 30 min, 2) rinsing with 90 % ethanol, 3) soaking in 50 % bleach for 7 min and 4) rinsing with sterilized water. Twenty seeds per species were sown in a magenta box containing 55 ml media [2.22 g MS, 15 g sucrose and 6 g agar per liter] and germinated under a 16 h photoperiod provided by fluorescent lamps (ca. 30 μmol m−2 s−1) at 22 °C. When seedlings reached a sufficient size, 1 cm2 leaf pieces were excised and placed onto 100 × 15 mm plates containing 25 ml callus-induction media [4.44 g MS, 100 μg 2,4-D, 20 g sucrose, and 8 g agar per liter] (Rao et al., 1973) and incubated in the dark at room temperature. Callus developed after approximately four weeks. Every few weeks until harvest, callus pieces were transferred to new plates containing the same media.
For the collection of shoot apices, flowers and trichomes, seed was sown in 288-cell plug trays in a soilless media and placed under in a greenhouse under intermittent mist at 24 °C. After emergence of the first true leaf, seedlings were removed from mist and placed in a glass-covered greenhouse at 20 °C under a 16 h photoperiod under ambient light, supplemented with 60 μmol m−2 s−1 provided by high-pressure sodium lamps. Shoot apices were collected when plants had unfolded 6–8 true leaves. Small leaves were removed with forceps to minimize contamination with leaf petiole tissue. The shoot apices from ca. 100 plants of each species were excised, placed in 15 ml conical tubes, and immediately frozen in liquid nitrogen. Flowers were harvested when a full range of reproductive stages of development were present, ranging from flower buds ca. 5 mm in length through the beginning fruit development. Harvested flowers were placed in paper envelopes, and frozen in liquid nitrogen. In order to harvest trichomes, stems and petioles from flowering Petunia plants were cut into ca. 5 cm pieces and immediately frozen in liquid nitrogen. The frozen tissue pieces were held using chilled forceps and the trichomes were removed from the stems by gentle scraping with a frozen paint brush and were collected in a mortar containing liquid nitrogen. Trichomes were immediately ground and transferred to 2 ml crew cap tubes for storage at −80 °C.
Total RNA was extracted from each tissue sample using the Trizol® reagent (LifeTechnologies) following the manufacturer’s instructions. Forty micrograms of RNA was treated for DNA contamination using RNase-free DNase set (Qiagen). DNase-treated RNA was purified using the RNeasy® MinElute Cleanup kit (Qiagen). RNA yield and quality were determined using gel electrophoresis, spectroscopy and the Agilent 2100 BioAnalyzer RNA 6000 Pico chip with RIN values ≥ 8.0 achieved (Agilent Technologies).
RNA-seq library construction, sequencing, and analysis
A TruSeq RNA Sample Preparation kit (Illumina) was used to construct the cDNA libraries for sequencing. Illumina barcodes allowed pooling of 7–8 libraries per lane. Sequencing was completed on the Illumina HiSeq 2500 next-generation sequencer to 100 nt paired-end at the Michigan State University Research Technology Support Facility (RTSF; http://rtsf.msu.edu/; East Lansing, MI). Raw read quality was assessed with FastQC (http://www.bioinformatics.babraham.ac.uk/ projects/fastqc/). Sequences were filtered and trimmed based on quality metrics and adapter sequences were removed with TrimmomaticPE . The TrimmomaticPE options employed included ILLUMINACLIP, SLIDINGWINDOW:5:20, LEADING:5, TRAILING:5, and HEADCROP:10–14. TrimmomaticPE outputs sequences into paired-end and single-end files requiring no further mate-specific filtering. Cleaned reads were reassessed with FastQC for quality visualization to ensure no further filtering was required.
De novo assembly and expression analysis
Reads meeting quality standards were de novo assembled using the Velvet (version 1.1.07) and Oases (version 0.2.08) packages [69, 70]. K-mer lengths of 55, 57, 59, 61, 63, 65, 67, and 69 were tested and metrics collected for each assembly. Criterion including the total number of transcripts produced, Velvet N50 length, average transcript length, and average transcript read coverage were considered. A k-mer of 61 was selected representing the best balance of metrics for all assemblies. Two de novo assemblies were completed for each Petunia transcriptome to reduce redundancy and limit the number of raw reads confounding assemblies. Reads from the callus, shoot apex and seedling libraries composed the first de novo assembly. Resulting transcripts were concatenated into a single pseudomolecule or ‘artificial chromosome one’ to serve as our core tissue reference transcriptome. Reads from the trichome and flower libraries were then mapped to artificial chromosome one with Bowtie (version 1.0.0) and TopHat (version 1.4.1) . Reads aligning to the artificial chromosome were discarded to eliminate redundancy of data already assembled. All unmapped reads from the flower and trichome libraries were then combined with the original callus, shoot apex, and seedling reads in a second de novo transcriptome assembly ; Supplemental Fig. 1). Transcripts from the second and final assembly were filtered for low complexity. Representative transcripts, comprised of the longest isoform were extracted from the final assemblies. Completeness of the assemblies was assessed by mapping all reads from each tissue library back to their respective de novo assembly (per species) individually with Bowtie and TopHat . CEGMA or Core Eukaryotic Genes Mapping Approach , a pipeline used to accurately annotate core genes in eukaryotic genomes, was also searched to determine the quality and completeness of the assemblies. Reads from each tissue library were mapped back to the de novo assembly with Bowtie and TopHat. Estimation of transcript abundance was then completed with Cufflinks [71, 72]. Cufflinks reports expression abundance in Fragments Per Kilobase of exon model per Million fragments mapped, or FPKM. Values consider the number of reads supporting each model and counts are normalized.
Functional annotation and identification of orthologous and paralogous gene families
Putative functions were assigned using BLASTX queries against the Arabidopsis thaliana protein database (TAIR10; http://www.arabidopsis.org/), UniRef100 (UniProt Knowledgebase Reference Clusters; release 63; http://www.uniprot.org/), RefSeq protein database (NCBI Reference Sequence Database; www.ncbi.nlm.nih.gov/refseq/), Pfam domains (version 27.0; pfam.xfam.org), and the NCBI non-redundant protein database (nr; www.ncbi.nlm.nih.gov/protein). BLASTX queries against all databases were searched using an e-value cutoff of 1e-5 and allowing a single hit. Hits were then processed accepting only those with > 30 % coverage and > 70 % identity to the subject. A total of 126 P. axillaris, 12 P. integrifolia, and two P. exserta representative transcripts were removed as suspected contamination based on their annotation. Functional annotation of the transcripts was also queried using web-based FastAnnotator . FastAnnotator facilitates the integration of the annotation tools Blast2GO, PRIAM, and RPS BLAST, providing gene ontology classifications (GO), enzyme and domain identification. Transcripts were translated in batch using OrfPredictor  and orthologous and paralogous proteins in the three species were assigned with OrthoMCL v2.0.2 with default parameters and an e-value cut-off of 1e-10 .
Sequence reads from each tissue type were pooled for each species and mapped to the Petunia axillaris transcriptome using the Burrows-Wheeler Aligner (BWA) program  with default values. Duplicate reads were removed after the initial alignment, to eliminate reads that mapped to the same position of the transcriptome. Duplicate removal was performed for the aligned reads using picardTools/1.89 (broadinstitute.github.io/picard). The subsequent local realignment to correct misalignments due to the presence of indels was performed using the Genome Analysis Toolkit (GATK)  with the IndelRealigner function. Bam files for each genotype were compressed using ReduceReads from GATK. The initial variant calling was performed by HaplotypeCaller from GATK with a Phred-scaled confidence threshold of 30. After the initial variant calling, SNPs were filtered based on several criteria. The first round of SNP filtering employed the following steps: exclude three SNPs within 100 bp; filter out variants with zero mapping quality constituting more than 10 % of all reads at that site; exclude SNPs covered only by sequences on the same strand with an FS value (Phred-scaled p-value using Fisher’s exact test to detect strand bias) >60; exclude SNPs with a minor allele frequency <0.01; and exclude SNPs with low QD (quality by depth), low quality (<11), and low total depth of coverage (<11) by the default parameters recommended from GATK best practice with the VariantFiltration tool [53, 76]. The second round of SNP filtering was performed by custom Python (Python 2.7.8)  scripts with the following criteria: SNPs within the beginning and end 30 bp of the reference transcripts were excluded; SNPs were selected with at least two genotypes having a depth of coverage greater than 10; exclude loci with a heterozygous genotype call in Petunia axillaris. The intron-exon boundaries were predicted by aligning the Petunia axillaris transcripts to the draft Petunia axillaris genome scaffold sequences with GMAP . SNPs within 30 bp of the exon-intron boundary region were then filtered out. Indels were not called because alternative splicing may impede reliable indel discovery. A custom Python script was used to extract the 5ʹ-UTR and 3ʹ-UTR and to calculate the distribution of SNPs in these regions. Depth of coverage was calculated by BEDTools . Annotated unigenes containing SNPs were visualized by BGI WEGO (web gene ontology annotation plotting) (http://wego.genomics.org.cn/cgi-bin/wego/index.pl). KEGG pathways were assigned to unigenes containing SNPs using the online KEGG Automatic Annotation Server (KAAS) (http://www.genome.jp/tools/kaas/) . KEGG Orthology (KO) assignment was applied using Bi-directional Best Hit (BBH) method with plant organisms (Arabidopsis thaliana (thale cress), Arabidopsis lyrata (lyrate rockcress), Theobroma cacao (cacao), Glycine max (soybean), Fragaria vesca (woodland strawberry), Vitis vinifera (wine grape), Solanum lycopersicum (tomato), and Oryza sativa L. ssp. japonica (Japanese rice) (RefSeq)) as references.
A total of 55 SNPs were selected for validation. Primers were designed by first choosing SNPs where at least two species exhibited polymorphism. Unigenes with exons greater than 600 bp were then selected, SNPs or alleles in the sequences were masked using IUBI./IUPAC nucleotide acid code. Batch Primer3 (http://probes.pw.usda.gov/cgi-bin/batchprimer3/batchprimer3.cgi) was used for primer design using the “SNP flanking primers design” option. The expected amplicon sizes ranged from 250 to 350 bp, with primer sizes ranging from 18 to 27 bp, and primer Tm ranging from 57 to 63 °C. Genomic DNA was extracted using the DNeasy plant mini kit (QIAGEN, Valencia, CA, USA). PCR amplification of genomic DNA was carried out in a 10 μl reaction containing 1 × PCR buffer, 0.2 μl 10 μM dNTPs, 0.6 U of DNA polymerase, and 5 ng template DNA. The following PCR profile was used: an initial denaturation at 95 °C for 30 s, followed by 11 cycles of 95 °C for 30 s, 65 °C to 55 °C for 30 s, with a decrease of 1 °C per cycle, and 72 °C for 30 s, followed 30 cycles of 95 °C for 30 s, 58 °C for 30 s, and 72 °C for 30 s, followed by a final extension at 72 °C for 5 min. Amplicons were purified by Agencourt AMPure XP beads (Beckman Coulter, Pasadena, CA) and quantified by BioDrop Duo (Isogen Life Science, the Netherlands). Sanger sequencing of the amplicons was performed on an ABI 3730xl platform (Life Technologies, Carlsbad, CA) at the Michigan State University Genomics Core Facility. Sequences amplified from the same primer set were initially aligned with CLC Sequence Viewer 6.9.1 (CLC Bio, Valencia, CA, USA) for preliminary SNP verification. Sequence chromatograms were visualized by Seq Scanner 2 (Life Technologies, Carlsbad, CA) for further SNP confirmation.
CAPS and SSR marker design
SNP markers were transformed into CAPS markers using the following pipeline. First, the SNP output was transformed into VCF (variant call format) format. Then, by using the previous GMAP output as reference, a new VCF file was generated by transforming the current SNP locations and genotypes to corresponding SNP locations and variant calls in the P. axillaris draft genome scaffolds by a custom Python script. Primers were designed with the galaxy-pcr-markers pipeline (https://github.com/cfljam/galaxy-pcr-markers) with the modification of producing only one pair of primers, minimum amplicon size of 200 bp and maximum amplicon size of 400 bp. A total of 20 commonly used restriction enzymes: AluI, ApaI, BamHI, BbrPI, BfrI, ClaI, DdeI, DpnII, DraI, EcoRI, HaeIII, HincII, HinfI, HpaI, PvuII, RsaI, SacI, Sau3AI, SmaI, and TaqI were selected for the pipeline. The unigene dataset were used for SSR identification with MISA (microsatellite identification tool) (http://pgrc.ipk-gatersleben.de/misa/). The SSR identification criteria were 6, 5, 5, 5, 5 repeats for di-, tri-, tetra-, penta-, and hexanucleotides, respectively. Primer pairs were designed from primer3_core (Primer3 2.3.6) (http://primer3.sourceforge.net/releases.php). Primer parameters were: designated amplicon size 100–280 bp, primer annealing temperature 55 to 65 °C, primer length 18 to 28 bp, and GC content from 45 to 55 %. Unigenes with homology to known genes involved in specifying plant development rate [16–18, 20, 25] were searched for possible CAPS markers from the entire CAPS output. The CAPS markers were used to genotype the three species by following the previous PCR protocol. CAPS markers were digested by the above mentioned restriction enzymes (New England Biolabs, Beverly, MA) and separated on 2 % NuSieve GTG agarose (Lonza, Basel, Switzerland) plus 1 % UltraPure agarose (Life Technologies, Carlsbad, CA) gel with 100 V for 2 h at room temperature.
Candidate gene mapping and QTL mapping of development rate in an interspecific hybrid Petunia population
An interspecific hybrid F2 population developed from a cross between P. integrifolia and P. axillaris containing 164 individuals was used to identify QTL for development rate. Population development, the measurement of development rate, and the SSR-based genetic linkage map were reported previously . CAPS markers developed from P. axillaris plastochron-related transcripts based on SNPs between P. integrifolia and P. axillaris were used to screen the F2 population. If no CAPS were previously designed from our SNP set, CAPS markers were manually designed with CAPS Designer (http://solgenomics.net/tools/caps_designer/caps_input.pl) and Primer3Plus (http://www.bioinformatics.nl/cgi-bin/primer3plus/primer3plus.cgi/). The updated genetic linkage map was generated with JoinMap 4.0  with the Kosambi mapping function . The recombination threshold was set to 0.3 and the logarithmic odds (LOD) score was set to 5. The locus order was calculated by the regression module of JoinMap4. Linkage group numbers were assigned according to the previous study . QTL mapping was conducted with MapQTL 6 . A permutation test with 1000 permutations was used to establish the LOD significance threshold. The initial QTL mapping was performed by interval mapping (a single-QTL model). Multiple QTL model mapping (MQM) was then used to reduce the background noise. Significant QTL and regions were graphically visualized using MapChart 2.1 .
Sequences are available under Genbank/EMBL/DDBJ BioProject numbers PRJNA262254, PRJNA262142, and PRJNA261953. These include SRA (Sequence Read Archive) runs for P. axillaris: SRR1585615, SRR1585635, SRR1585830, and SRR1585954-1585955; P. exserta: SRR1586492-1586494, SRR1586500, and SRR1586504; P. integrifolia: SRR1587109, SRR1587150-1587151, and SRR1587153-1587154. Transcriptome Shotgun Assembly projects were also deposited under accessions GBRU00000000 (P. axillaris), GBRT00000000 (P. exserta), and GBRV00000000 (P. integrifolia). The versions described in this paper are the first versions of these assemblies (GBRU01000000, GBRT01000000, GBRV01000000). The SNPs are available at NCBI dbSNP under accession numbers ss1750993386 - ss1751108187.
Funding was provided by USDA Specialty Crops Research Initiative Award number 2011-51181-30666 to R.M.W. and C.S.B. and National Science Foundation award number IOS-1025636 and a Michigan State University (MSU) Foundation Strategic Partnership Grant to C.S.B. R.M.W. and C.S.B. are supported in part by Michigan AgBioResearch and through USDA National Institute of Food and Agriculture, Hatch project numbers MICL02121 and MICL02265. The authors would like to thank the staff at the Michigan State University Research Technology Support Facility for TruSeq library preparation and transcriptome sequencing. We would also like to thank Dr. Jin-Ho Kang for collection of trichomes and Stephanie Rett for her technical support and assistance with SNP validation.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Reck-Kortmann M, Silva-Arias GA, Segatto ALA, Mader C, Bonatto SL, de Freitas LB. Multilocus phylogeny reconstruction: New insights into the evolutionary history of the genus Petunia. Mol Phylogenet Evol. 2014;81:19–28.View ArticlePubMedGoogle Scholar
- Stehmann JR, Lorenz-Lemke AP, Freitas LB, Semir J. The genus Petunia. In: Gerats T, Strommer J, editors. Petunia Evolutionary, Developmental and Physiological Genetics. New York: Springer; 2009. p. 1–28.Google Scholar
- Morita Y, Saito R, Ban Y, Tanikawa N, Kuchitsu K, Ando T, et al. Tandemly arranged chalcone synthase A genes contribute to the spatially regulated expression of siRNA and the natural bicolor floral phenotype in Petunia hybrida. Plant J. 2012;70(5):739–49.View ArticlePubMedGoogle Scholar
- Kelly RO, Deng ZA, Harbaugh BK. Evaluation of 125 petunia cultivars as bedding plants and establishment of class standards. Horttechnology. 2007;17(3):386–96.Google Scholar
- Kroon J, Souer E, Degraaff A, Xue YB, Mol J, Koes R. Cloning and structural analysis of the anthocyanin pigmentation locus Rt of Petunia hybrida: Characterization of insertion sequences in two mutant alleles. Plant J. 1994;5(1):69–80.View ArticlePubMedGoogle Scholar
- Holton TA, Brugliera F, Lester DR, Tanaka Y, Hyland CD, Menting JGT, et al. Cloning and expression of cytochrome P450 genes controlling flower color. Nature. 1993;366(6452):276–9.View ArticlePubMedGoogle Scholar
- Quattrocchio F, Wing J, van der Woude K, Souer E, de Vetten N, Mol J, et al. Molecular analysis of the anthocyanin2 gene of petunia and its role in the evolution of flower color. Plant Cell. 1999;11(8):1433–44.PubMed CentralView ArticlePubMedGoogle Scholar
- Bartok JW. Energy conservation for commercial greenhouses. New York: Cooperative Extension, Ithaca; 2001.Google Scholar
- Blanchard MG, Runkle ES. The influence of day and night temperature fluctuations on growth and flowering of annual bedding plants and greenhouse heating cost predictions. Hortscience. 2011;46(4):599–603.Google Scholar
- Warner RM, Walworth AE. Quantitative inheritance of crop timing traits in interspecific hybrid Petunia populations and interactions with crop quality parameters. J Hered. 2010;101(3):308–16.View ArticlePubMedGoogle Scholar
- Vallejo VA, Tychonievich J, Lin W-K, Wangchu L, Barry CS, Warner RM. Identification of QTL for crop timing and quality traits in an interspecific Petunia population. Molecular Breeding. 2015, 35(2):doi:10.1007/s11032-11015-10218-11034.Google Scholar
- Werner T, Motyka V, Strnad M, Schmulling T. Regulation of plant growth by cytokinin. Proc Natl Acad Sci U S A. 2001;98(18):10487–92.PubMed CentralView ArticlePubMedGoogle Scholar
- Lohmann D, Stacey N, Breuninger H, Jikumaru Y, Muller D, Sicard A, et al. SLOW MOTION is required for within-plant auxin homeostasis and normal timing of lateral organ initiation at the shoot meristem in Arabidopsis. Plant Cell. 2010;22(2):335–48.PubMed CentralView ArticlePubMedGoogle Scholar
- Chen MK, Wilson RL, Palme K, Ditengou FA, Shpak ED. ERECTA family genes regulate auxin transport in the shoot apical meristem and forming leaf primordia. Plant Physiol. 2013;162(4):1978–91.PubMed CentralView ArticlePubMedGoogle Scholar
- Guenot B, Bayer E, Kierzkowski D, Smith RS, Mandel T, Zadnikova P, et al. PIN1-independent leaf initiation in Arabidopsis. Plant Physiol. 2012;159(4):1501–10.PubMed CentralView ArticlePubMedGoogle Scholar
- Miyoshi K, Ahn BO, Kawakatsu T, Ito Y, Itoh JI, Nagato Y, et al. PLASTOCHRON1, a timekeeper of leaf initiation in rice, encodes cytochrome P450. Proc Natl Acad Sci U S A. 2004;101(3):875–80.PubMed CentralView ArticlePubMedGoogle Scholar
- Kawakatsu T, Itoh J, Miyoshi K, Kurata N, Alvarez N, Veit B, et al. PLASTOCHRON2 regulates leaf initiation and maturation in rice. Plant Cell. 2006;18(3):612–25.PubMed CentralView ArticlePubMedGoogle Scholar
- Wang JW, Schwab R, Czech B, Mica E, Weigel D. Dual effects of miR156-targeted SPL genes and CYP78A5/KLUH on plastochron length and organ size in Arabidopsis thaliana. Plant Cell. 2008;20(5):1231–43.PubMed CentralView ArticlePubMedGoogle Scholar
- Grigg SP, Canales C, Hay A, Tsiantis M. SERRATE coordinates shoot meristem function and leaf axial patterning in Arabidopsis. Nature. 2005;437(7061):1022–6.View ArticlePubMedGoogle Scholar
- Prigge MJ, Wagner DR. The Arabidopsis SERRATE gene encodes a zinc-finger protein required for normal shoot development. Plant Cell. 2001;13(6):1263–79.PubMed CentralView ArticlePubMedGoogle Scholar
- Chaudhury AM, Letham S, Craig S, Dennis ES. amp1 - a mutant with high cytokinin levels and altered embryonic pattern, faster vegetative growth, constitutive photomorphogenesis and precocious flowering. Plant J. 1993;4(6):907–16.View ArticleGoogle Scholar
- Helliwell CA, Chin-Atkins AN, Wilson IW, Chapple R, Dennis ES, Chaudhury A. The Arabidopsis AMP1 gene encodes a putative glutamate carboxypeptidase. Plant Cell. 2001;13(9):2115–25.PubMed CentralView ArticlePubMedGoogle Scholar
- Li SB, Liu L, Zhuang XH, Yu Y, Liu XG, Cui X, et al. MicroRNAs inhibit the translation of target mRNAs on the endoplasmic reticulum in Arabidopsis. Cell. 2013;153(3):562–74.PubMed CentralView ArticlePubMedGoogle Scholar
- Evans MMS, Poethig RS. The viviparous8 mutation delays vegetative phase change and accelerates the rate of seedling growth in maize. Plant J. 1997;12(4):769–79.View ArticleGoogle Scholar
- Kawakatsu T, Taramino G, Itoh JI, Allen J, Sato Y, Hong SK, et al. PLASTOCHRON3/GOLIATH encodes a glutamate carboxypeptidase required for proper development in rice. Plant J. 2009;58(6):1028–40.View ArticlePubMedGoogle Scholar
- Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009;10(1):57–63.PubMed CentralView ArticlePubMedGoogle Scholar
- Ozsolak F, Milos PM. RNA sequencing: advances, challenges and opportunities. Nat Rev Genet. 2011;12(2):87–98.PubMed CentralView ArticlePubMedGoogle Scholar
- Schilmiller AL, Miner DP, Larson M, McDowell E, Gang DR, Wilkerson C, et al. Studies of a biochemical factory: tomato trichome deep expressed sequence tag sequencing and proteomics. Plant Physiol. 2010;153(3):1212–23.PubMed CentralView ArticlePubMedGoogle Scholar
- Varshney RK, Nayak SN, May GD, Jackson SA. Next-generation sequencing technologies and their implications for crop genetics and breeding. Trends Biotechnol. 2009;27(9):522–30.View ArticlePubMedGoogle Scholar
- Hamilton JP, Sim SC, Stoffel K, Van Deynze A, Buell CR, Francis DM. Single nucleotide polymorphism discovery in cultivated tomato via sequencing by synthesis. Plant Genome. 2012;5(1):17–29.View ArticleGoogle Scholar
- Barbazuk WB, Emrich SJ, Chen HD, Li L, Schnable PS. SNP discovery via 454 transcriptome sequencing. Plant J. 2007;51(5):910–8.PubMed CentralView ArticlePubMedGoogle Scholar
- Hirsch CN, Foerster JM, Johnson JM, Sekhon RS, Muttoni G, Vaillancourt B, et al. Insights into the maize pan-genome and pan-transcriptome. Plant Cell. 2014;26(1):121–35.PubMed CentralView ArticlePubMedGoogle Scholar
- Bedewitz MA, Gongora-Castillo E, Uebler JB, Gonzales-Vigil E, Wiegert-Rininger KE, Childs KL, et al. A root-expressed L-phenylalanine: 4-hydroxyphenylpyruvate aminotransferase is required for tropane alkaloid biosynthesis in Atropa belladonna. Plant Cell. 2014;26(9):3745–62.PubMed CentralView ArticlePubMedGoogle Scholar
- Galliot C, Hoballah ME, Kuhlemeier C, Stuurman J. Genetics of flower size and nectar volume in Petunia pollination syndromes. Planta. 2006;225(1):203–12.View ArticlePubMedGoogle Scholar
- Klahre U, Gurba A, Hermann K, Saxenhofer M, Bossolini E, Guerin PM, et al. Pollinator choice in petunia depends on two major genetic loci for floral scent production. Curr Biol. 2011;21(9):730–9.View ArticlePubMedGoogle Scholar
- Bossolini E, Klahre U, Brandenburg A, Reinhardt D, Kuhlemeier C. High resolution linkage maps of the model organism Petunia reveal substantial synteny decay with the related genome of tomato. Genome. 2011;54(4):327–40.View ArticlePubMedGoogle Scholar
- Strommer J, Gerats AGM, Sanago M, Molnar SJ. A gene-based RFLP map of Petunia. Theor Appl Genet. 2000;100(6):899–905.View ArticleGoogle Scholar
- Strommer J, Peters J, Zethof J, De Keukeleire P, Gerats T. AFLP maps of Petunia hybrida: building maps when markers cluster. Theor Appl Genet. 2002;105(6–7):1000–9.PubMedGoogle Scholar
- Villarino GH, Bombarely A, Giovannoni JJ, Scanlon MJ, Mattson NS. Transcriptomic analysis of Petunia hybrida in response to salt stress using high throughput RNA sequencing. Plos One. 2014;9(4):e94651.PubMed CentralView ArticlePubMedGoogle Scholar
- Williams JS, Der JP, dePamphilis CW, Kao TH. Transcriptome analysis reveals the same 17 S-locus F-box genes in two haplotypes of the self-incompatibility locus of Petunia inflata. Plant Cell. 2014;26(7):2873–88.PubMed CentralView ArticlePubMedGoogle Scholar
- Zenoni S, D’Agostino N, Tornielli GB, Quattrocchio F, Chiusano ML, Koes R, et al. Revealing impaired pathways in the an11 mutant by high-throughput characterization of Petunia axillaris and Petunia inflata transcriptomes. Plant J. 2011;68(1):11–27.View ArticlePubMedGoogle Scholar
- Gongora-Castillo E, Fedewa G, Yeo Y, Chappell J, DellaPenna D, Buell CR. Genomic approaches for interrogating the biochemistry of medicinal plant species. Method Enzymol. 2012;517:139–59.View ArticleGoogle Scholar
- Parra G, Bradnam K, Korf I. CEGMA: a pipeline to accurately annotate core genes in eukaryotic genomes. Bioinformatics. 2007;23(9):1061–7.View ArticlePubMedGoogle Scholar
- Zerbe P, Hamberger B, Yuen MMS, Chiang A, Sandhu HK, Madilao LL, et al. Gene discovery of modular diterpene metabolism in nonmodel systems. Plant Physiol. 2013;162(2):1073–91.PubMed CentralView ArticlePubMedGoogle Scholar
- Nakasugi K, Crowhurst RN, Bally J, Wood CC, Hellens RP, Waterhouse PM. De novo transcriptome sequence assembly and analysis of RNA silencing genes of Nicotiana benthamiana. Plos One. 2013;8(3):e59534.PubMed CentralView ArticlePubMedGoogle Scholar
- Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009;25(9):1105–11.PubMed CentralView ArticlePubMedGoogle Scholar
- Mather K. Specific differences in Petunia I. Incompatibility. J Genet. 1943;45(3):215–35.View ArticleGoogle Scholar
- Garg R, Patel RK, Tyagi AK, Jain M. De novo assembly of chickpea transcriptome using short reads for gene discovery and marker identification. DNA Res. 2011;18(1):53–63.PubMed CentralView ArticlePubMedGoogle Scholar
- Yates SA, Swain MT, Hegarty MJ, Chernukin I, Lowe M, Allison GG, et al. De novo assembly of red clover transcriptome based on RNA-Seq data provides insight into drought response, gene discovery and marker identification. BMC Genomics. 2014;15:453.PubMed CentralView ArticlePubMedGoogle Scholar
- Li L, Stoeckert CJ, Roos DS. OrthoMCL: Identification of ortholog groups for eukaryotic genomes. Genome Res. 2003;13(9):2178–89.PubMed CentralView ArticlePubMedGoogle Scholar
- McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, et al. The genome analysis toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res. 2010;20(9):1297–303.PubMed CentralView ArticlePubMedGoogle Scholar
- Quinn EM, Cormican P, Kenny EM, Hill M, Anney R, Gill M, et al. Development of strategies for SNP detection in RNA-seq data: application to lymphoblastoid cell lines and evaluation using 1000 genomes data. Plos One. 2013;8(3):e58815.PubMed CentralView ArticlePubMedGoogle Scholar
- DePristo MA, Banks E, Poplin R, Garimella KV, Maguire JR, Hartl C, et al. A framework for variation discovery and genotyping using next-generation DNA sequencing data. Nat Genet. 2011;43(5):491–8.PubMed CentralView ArticlePubMedGoogle Scholar
- Niu S-H, Li Z-X, Yuan H-W, Chen X-Y, Li Y, Li W. Transcriptome characterisation of Pinus tabuliformis and evolution of genes in the Pinus phylogeny. BMC Genomics. 2013;14:263.PubMed CentralView ArticlePubMedGoogle Scholar
- Miguel Blanca J, Canizares J, Ziarsolo P, Esteras C, Mir G, Nuez F, et al. Melon transcriptome characterization: Simple sequence repeats and single nucleotide polymorphisms discovery for high throughput genotyping across the species. Plant Genome. 2011;4(2):118–31.View ArticleGoogle Scholar
- Coulondre C, Miller JH, Farabaugh PJ, Gilbert W. Molecular basis of base substitution hotspots in Escherichia coli. Nature. 1978;274(5673):775–80.View ArticlePubMedGoogle Scholar
- Keller I, Bensasson D, Nichols R. Transition-transversion bias is not universal: a counter example from grasshopper pseudogenes. PLoS Genet. 2007;3(2):e22.PubMed CentralView ArticlePubMedGoogle Scholar
- Brousseau L, Tinaut A, Duret C, Lang T, Garnier-Gere P, Scotti I. High-throughput transcriptome sequencing and preliminary functional analysis in four Neotropical tree species. BMC Genomics. 2014;15:238.PubMed CentralView ArticlePubMedGoogle Scholar
- Chen J, Uebbing S, Gyllenstrand N, Lagercrantz U, Lascoux M, Kallman T. Sequencing of the needle transcriptome from Norway spruce (Picea abies Karst L.) reveals lower substitution rates, but similar selective constraints in gymnosperms and angiosperms. BMC Genomics. 2012;13:589.PubMed CentralView ArticlePubMedGoogle Scholar
- Duarte J, Riviere N, Baranger A, Aubert G, Burstin J, Cornet L, et al. Transcriptome sequencing for high throughput SNP development and genetic mapping in Pea. BMC Genomics. 2014;15:126.PubMed CentralView ArticlePubMedGoogle Scholar
- Howe GT, Yu J, Knaus B, Cronn R, Kolpak S, Dolan P, et al. A SNP resource for Douglas-fir: de novo transcriptome assembly and SNP detection and validation. BMC Genomics. 2013;14:137.PubMed CentralView ArticlePubMedGoogle Scholar
- Hamilton JP, Hansey CN, Whitty BR, Stoffel K, Massa AN, Van Deynze A, et al. Single nucleotide polymorphism discovery in elite north american potato germplasm. BMC Genomics. 2011;12:302.PubMed CentralView ArticlePubMedGoogle Scholar
- Raji AA, Anderson JV, Kolade OA, Ugwu CD. Gene-based microsatellites for cassava (Manihot esculenta): prevalence, polymorphisms, and cross-taxa utility. BMC Plant Biol. 2009;9:118.PubMed CentralView ArticlePubMedGoogle Scholar
- Poncet V, Rondeau M, Tranchant C, Cayrel A, Hamon S. SSR mining in coffee tree EST databases: potential use of EST-SSRs as markers for the Coffea genus. Mol Genet Genomics. 2006;276(5):436–49.View ArticlePubMedGoogle Scholar
- Guo Y, Khanal S, Tang S, Bowers J, Heesacker A. Comparative mapping in intraspecific populations uncovers a high degree of macrosynteny between A-and B-genome diploid species of peanut. BMC Genomics. 2012;13:608.PubMed CentralView ArticlePubMedGoogle Scholar
- Tychonievich J, Wangchu L, Barry C, Warner RM. Utilizing wild species for marker-assisted selection of crop timing and quality traits in Petunia. Acta Hortic. 2013;1000:465–9.Google Scholar
- Zubko E, Adams CJ, Machaekova I, Malbeck J, Scollan C, Meyer P. Activation tagging identifies a gene from Petunia hybrida responsible for the production of active cytokinins in plants. Plant J. 2002;29(6):797–808.View ArticlePubMedGoogle Scholar
- Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.PubMed CentralView ArticlePubMedGoogle Scholar
- Zerbino DR, Birney E. Velvet: algorithms for de novo short read assembly using de Bruijn graphs. Genome Res. 2008;18(5):821–9.PubMed CentralView ArticlePubMedGoogle Scholar
- Schulz MH, Zerbino DR, Vingron M, Birney E. Oases: robust de novo RNA-seq assembly across the dynamic range of expression levels. Bioinformatics. 2012;28(8):1086–92.PubMed CentralView ArticlePubMedGoogle Scholar
- Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, et al. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc. 2012;7(3):562–78.PubMed CentralView ArticlePubMedGoogle Scholar
- Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotech. 2010;28(5):511–5.View ArticleGoogle Scholar
- Chen TW, Gan RCR, Wu TH, Huang PJ, Lee CY, Chen YYM, et al. FastAnnotator- an efficient transcript annotation web tool. BMC Genomics. 2012;13(7):S9.Google Scholar
- Min XJ, Butler G, Storms R, Tsang A. OrfPredictor: predicting protein-coding regions in EST-derived sequences. Nucleic Acids Res. 2005;33:W677–80.PubMed CentralView ArticlePubMedGoogle Scholar
- Li H, Durbin R. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics. 2009;25(14):1754–60.PubMed CentralView ArticlePubMedGoogle Scholar
- Van der Auwera GA, Carneiro MO, Hartl C, Poplin R, del Angel G, Levy-Moonshine A, Jordan T, Shakir K, Roazen D, Thibault J et al.: From FastQ data to high-confidence variant calls: The genome analysis toolkit best practices pipeline. Current Protocols in Bioinformatics. 2013;43:11.10.1-11.10-33. Google Scholar
- van Rossum G, de Boer J. Interactively testing remote servers using the Python programming language. CWI Quarterly. 1991;4(4):283–303.Google Scholar
- Wu TD, Watanabe CK. GMAP: a genomic mapping and alignment program for mRNA and EST sequences. Bioinformatics. 2005;21(9):1859–75.View ArticlePubMedGoogle Scholar
- Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26(6):841–2.PubMed CentralView ArticlePubMedGoogle Scholar
- Moriya Y, Itoh M, Okuda S, Yoshizawa AC, Kanehisa M. KAAS: an automatic genome annotation and pathway reconstruction server. Nucleic Acids Res. 2007;35:W182–5.PubMed CentralView ArticlePubMedGoogle Scholar
- Van Ooijen JW. JoinMap® 4, Software for the calculation of genetic linkage maps in experimental populations. Wageningen, Netherlands: Kyazma BV; 2006.Google Scholar
- Kosambi DD. The estimation of map distances from recombination values. Ann Eugenics. 1944;12:172–5.View ArticleGoogle Scholar
- Van Ooijen JW. MapQTL 6®, Software for the mapping of quantitative trait loci in experimental populations of diploid species. Wageningen, Netherlands: Kyazma BV; 2009.Google Scholar
- Voorrips RE. MapChart: software for the graphical presentation of linkage maps and QTLs. J Hered. 2002;93(1):77–8.View ArticlePubMedGoogle Scholar