De novo assembly and transcriptome analysis of five major tissues of Jatropha curcas L. using GS FLX titanium platform of 454 pyrosequencing

Background Jatropha curcas L. is an important non-edible oilseed crop with promising future in biodiesel production. However, factors like oil yield, oil composition, toxic compounds in oil cake, pests and diseases limit its commercial potential. Well established genetic engineering methods using cloned genes could be used to address these limitations. Earlier, 10,983 unigenes from Sanger sequencing of ESTs, and 3,484 unique assembled transcripts from 454 pyrosequencing of uncloned cDNAs were reported. In order to expedite the process of gene discovery, we have undertaken 454 pyrosequencing of normalized cDNAs prepared from roots, mature leaves, flowers, developing seeds, and embryos of J. curcas. Results From 383,918 raw reads, we obtained 381,957 quality-filtered and trimmed reads that are suitable for the assembly of transcript sequences. De novo contig assembly of these reads generated 17,457 assembled transcripts (contigs) and 54,002 singletons. Average length of the assembled transcripts was 916 bp. About 30% of the transcripts were longer than 1000 bases, and the size of the longest transcript was 7,173 bases. BLASTX analysis revealed that 2,589 of these transcripts are full-length. The assembled transcripts were validated by RT-PCR analysis of 28 transcripts. The results showed that the transcripts were correctly assembled and represent actively expressed genes. KEGG pathway mapping showed that 2,320 transcripts are related to major biochemical pathways including the oil biosynthesis pathway. Overall, the current study reports 14,327 new assembled transcripts which included 2589 full-length transcripts and 27 transcripts that are directly involved in oil biosynthesis. Conclusion The large number of transcripts reported in the current study together with existing ESTs and transcript sequences will serve as an invaluable genetic resource for crop improvement in jatropha. Sequence information of those genes that are involved in oil biosynthesis could be used for metabolic engineering of jatropha to increase oil content, and to modify oil composition.


Background
The genus Jatropha of Euphorbiaceae family contains about 175 known species [1] of which Jatropha curcas L. is the most promising and suitable species for biodiesel production worldwide. Biodiesel is produced by transesterification of the oil to fatty acid methyl esters. Being non-edible crop, use of jatropha oil for biodiesel production does not threaten food security. Its cultivation in degraded soil can control erosion, and help in land reclamation. It is used for manufacturing soap, purgative agents, candles, coloring dyes and astringents [2].
The oil content of Jatropha which is around 30% could be increased to about 50% in order to make it commercially viable. In terms of oil composition, decreasing the content of unsaturated fatty acids (to increase oxidative stability), free fatty acids (to prevent soap formation and to increase the yield of biodiesel), and 18-carbon fatty acids (to lower the viscosity for better atomization of the biodiesel) would help to improve the quality of jatropha biodiesel. Jatropha oil cake -the by product of oil extraction, is not used as animal feed due to the presence of toxic curcin and phorbol esters.
Since it is not economical to remove these compounds, it would be desirable to block their accumulation itself by seed-specific silencing of the relevant genes. This will add economic value to the crop. Large scale cultivation of jatropha also requires the development of varieties that are tolerant to drought, and resistant to pests and diseases.
Genetic engineering methods could play a major role in jatropha crop improvement, because the scope for classical breeding is limited due to longer breeding cycle. Genetic manipulations in jatropha require genomic information and cloning of all the important genes. The genome size of J. curcas is estimated to be 410Mb [7]. The jatropha genome has been fully sequenced by Synthetic Genomics Inc, California, USA but it is not available for public use [8]. Partial genomic sequence of 285.8Mb is available in public databases but it is not assembled and annotated [9]. Natarajan et al., 2010 [10] have reported 12,084 ESTs using a normalized cDNA library from developing seeds, and Costa et al., 2010 [11] have reported 13,249 ESTs using non-normalized cDNA libraries from developing and germinating endosperm. Contig assembly of these ESTs showed the presence of 10,983 unigenes (6,361 from Natarajan et al., 2010 [10] and 4,622 from Costa et al., 2010 [11]). Recently, Sato et al., 2011 [9] have reported 991,050 pyrosequencing reads from leaf and callus transcriptome which were assembled in to 4,751 contigs. Hybrid assembly of the unigenes and the assembled transcripts revealed the presence of 14,467 unique sequences (10,983 unigenes and 3,484 unique assembled transcripts). This is far below the number of sequences required to represent the whole transcriptome of jatropha. Annotation of the currently available unique sequences showed that several important genes are yet to be cloned.
High throughput 454 pyrosequencing of uncloned cDNAs is a powerful method for whole genome transcriptome analysis and gene discovery. 454 pyrosequencing using GS FLX platform was done in many plants including arabidopsis [12], artemisia [13], cucumber [14], medicago [15], maize [16] and barley [17], and an average read length of 100-200 bases was obtained. Shorter read length was a major problem in assembly, especially in case of de novo assembly of the data from novel organisms, which do not have previously assembled and annotated reference sequences [18,19]. This limitation has been quickly overcome by using the newly launched titanium platform of GS FLX which increased the average read length to as high as 422 bases [20]. Concurrent improvements in de novo assembly software made it possible to assemble large number of full-length transcripts from novel organisms using 454 pyrosequencing data. Here, we report 454 pyrosequencing of normalized cDNAs from roots, mature leaves, flowers, developing seeds, and mature embryo of J. curcas using GS FLX titanium platform, and report de novo assembly of 17,457 transcripts.

cDNA synthesis and normalization
Total RNA was isolated from roots, mature leaves, flowers, developing seeds, and embryos of J. curcas. Quality of the RNA as determined by agarose gel electrophoresis (additional file 1) and OD 260 /OD 280 ratio (1.9 ± 0.05) was found to be suitable for cDNA synthesis. The total RNA from the five tissues were pooled and normalized cDNA was synthesised. Normalization of the cDNA greatly reduces the frequency of abundant transcripts, and increases the rate recovery of unique transcripts [10]. The efficiency of normalization was monitored by doing a parallel normalization reaction using chloramphenicol resistance gene as reporter gene. Redundant rate of the reporter gene was reduced from 1% to 0.025% after normalization which indicated 40 fold reduction in abundance. The normalized cDNA was subjected to quality control experiments before using it for 454 pyrosequencing. A test cDNA library was constructed using an aliquote of the normalized cDNA, and clones were randomly selected for further analysis. The average size of the normalized cDNAs as determined by PCR was about 2.3 kb (additional file 1) which is higher than average transcript size. Sequencing and BLASTX analysis revealed that 90% of the clones were full-length at the 5' end. These results showed that the normalized cDNA is highly suitable for 454 pyrosequencing.

pyrosequencing using GS FLX titanium platform
Half-plate 454 pyrosequencing reaction of the normalized cDNA was done using GS FLX titanium platform. It produced 197.7 Mb data from 383,918 reads with an average read length of 515 bases (53 to 1201 bases). First, the sequences were filtered using a minimum quality cut off value of 40. Then, the additional sequences that were added during cDNA synthesis were trimmed off. Poly (A/T) sequences were retained as they were reported to be useful in transcript assembly [21]. Sequences with less than 50 bases were removed before assembly. Finally, 381,957 ready-to-assemble reads were obtained. Size distribution of these reads is shown in figure 1. Length of these reads ranged between 50 and 1169 bases with an average of 343 bases.

De novo assembly
Reference assembly with the partial genomic sequence of jatropha [9] showed mapping of 95.87% of the reads, and the consensus accuracy was 99.15%. The details of reference assembly are given in additional file 2.
However, it could not be used for transcript assembly because the genome is not assembled and annotated. Therefore, de novo assembly was done using GS de novo assembler. This assembler can assemble the data under genomic or cDNA option. The genomic option assembles overlapping reads and reports consensus sequences as contigs. The cDNA option continues with linking the contigs together to form isotigs (splice variants) and isogroups (all the splice variants of individual transcripts). Since the objective of this study was only to assemble the transcripts, the genomic option was used. GS de novo assembler version 2.5p1 was reported to be better than the version 2.3 for de novo assembly using cDNA option [21]. Therefore, we have evaluated the currently available three versions of the assembler (v.2.3, v.2.5p1 and v.2.5.2) using genomic option. Contig assembly was done under the stringent assembly parameters of 40 bases overlap and 95% overlap identity. The results did not show significant differences among the versions in terms of the number of reads assembled into contigs, number of contigs assembled, average contig length, and the length of the longest contig (Table 1). Therefore, we have used the output from the latest version of the assembler (v.2.5.2) for further analysis.
The number of ESTs that are available from jatropha based on Sanger sequencing is 25,333 [10,11]. From the first 454 pyrosequencing study in jatropha, 991,050 reads with an average read length of 407 bases were reported by Sato et al., 2011 [9]. Analysis of these reads showed the presence of 822,387 repeat reads (83%), thereby, reducing the effective number of reads to 168,663. We have obtained 381,957 reads with an average read length of 343 bases, and there were only 36 repeat reads (0.01%). This may be because we have used cDNAs that were normalized against abundant transcripts. In addition, we have used a more diverse pool of cDNAs from five different tissues. These two factors have also played a role in the number of contigs assembled and the contig length. Contig assembly showed 4,751 contigs with an average contig length of 749 bases, and 17,457 contigs with an average contig length of 916 bases from the effective reads of Sato et al., 2011 [9] and the current study, respectively. The sequences of the assembled contigs are given in additional file 3 and 4. The size distribution of the contigs obtained from the current study is shown in figure 2. Hybrid assembly of the existing ESTs [10,11] and assembled contigs [9] revealed the presence of 10,983 unique ESTs from Sanger sequencing, and 3,484 unique transcripts from 454 pyrosequencing. Similar analysis of the 17,457 contigs that were assembled in the present study uncovered 14,327 unique transcripts. Together, at least partial sequences of 28,794 expressed genes are currently known from jatropha.

Annotation
Similarity search for all the contigs was done against non-redundant protein sequences database (nr) using BLASTX, and the 15,925 contigs that showed significant similarity (e-value cut off of <10 -10 ) were  annotated. It was found that 2,589 annotated contigs were full-length, and their size ranged from 247 to 5,465 bp. These sequence information could be used for direct amplification of the desired genes. The remaining 13,336 contigs were partial at 3' end (3,889 contigs) or 5' end (3,980 contigs) or both ends (5,467 contigs). These sequence information will be useful in cloning full-length genes using methods like rapid amplification of cDNA ends (RACE) PCR. However, these methods are cumbersome, and are not amenable for high throughput applications. Therefore, it would be easier to increase the number of full-length transcripts by increasing the number of reads from 454 pyrosequencing.

Pathway mapping of contigs by KEGG
Ortholog assignment and mapping of the contigs to the biological pathways were performed using KEGG automatic annotation server (KAAS). All the contigs were compared against the KEGG database using BLASTX with threshold bit-score value of 60 (default). It assigned EC numbers for 2,320 contigs, and they were mapped to respective pathways. The mapped contigs represented metabolic pathways of major biomolecules such as carbohydrates, lipids, nucleotides, amino acids, glycans, cofactors, vitamins, terpenoids, polyketides, etc. The KEGG pathway analysis also showed that 549 and 241 contigs represent the major metabolic pathways and biosynthetic pathways of secondary metabolites, respectively. The mapped contigs also represented the genes involved in genetic information processing (transcription, translation, folding, sorting and degradation, replication and repair, RNA family), environmental information processing (membrane transport, signal transduction, signaling molecules and interaction) and cellular processes (transport and catabolism, cell motility, cell growth and death, cell communication).

Genes involved in oil biosynthesis
Jatropha is popularly grown for biodiesel production. However, making it an ideal biodiesel crop requires genetic manipulations for increased oil yield and modified oil composition using the genes that are involved in oil biosynthesis pathway. Therefore, the contigs that were mapped to oil biosynthesis pathway were analyzed in detail. Biosynthesis of oil (triacyl glycerol) in plants takes place in three major steps, (a) biosynthesis of fatty acids in plastids, (b) activation and transport of fatty acids to endoplasmic reticulum and (c) synthesis of triacyl glycerol or oil. In Arabidopsis, 36 enzymes and 2 proteins are involved in oil biosynthesis [22]. Twentyeight of these enzymes are encoded by single genes, and the remaining 8 enzymes and the 2 proteins are encoded by small gene families. In Jatropha, 29 nonredundant sequences coding for 20 enzymes and 2 proteins that are involved in oil biosynthesis were already reported [9][10][11] but all of them were partial sequences. From this study, full-length transcripts were obtained for 17 of these partial sequences and significantly longer transcripts were obtained for another 5 partial sequences ( Figure 3). In addition, 11 full-length and 16 partial transcripts coding for 15 enzymes and 1 protein were also identified. In total, 56 sequences related to oil biosynthesis were found from this study, and 27 of them are reported for the first time in jatropha ( Table 2).

Biosynthesis of fatty acids
Oxidative decarboxylation of pyruvate to acety-CoA by pyruvate dehydrogenase (PDH) is the first committing step in fatty acid biosynthesis. PDH is a four subunit enzyme complex having pyruvate decarboxylase (E1-a and E1-b subunits), dihydrolipoyl acetyltransferase (E2 subunit) and dihydrolipoamide dehydrogenase (E3 subunit) enzymes [23]. Contigs representing all the four subunits were identified from the current study. Fatty acid biosynthesis is initiated by the condensation of acetyl-CoA to malonyl-CoA which is catalyzed by heteromeric acetyl CoA carboxylase (ACCase). It consists of three nuclear encoded subunits (biotin carboxyl carrier protein, biotin carboxylase and carboxyl transferase α-subunit) and a plastid encoded subunit (carboxyl transferase β-subunit) [24]. Accumulation of the mRNAs for these four subunits is co-ordinately regulated so that a constant molar stoichiometric ratio is maintained [25]. However, the plastid encoded subunit seems to be crucial for the accumulation of heteromeric ACCase, and for increasing the oil content [26,27]. Plants also have a homomeric ACCase localized in cytosol that does not play any major role in fatty acid biosynthesis that takes place almost exclusively in plastids. However, by targeting this enzyme to plastids, oil content was increased five-fold in potato tubers [28] and 5% in rapeseed [29]. We have found full-length contigs for all the four subunits of heteromeric ACCase and a 5' partial contig (3,711 bases) for homomeric ACCase. Condensation of malonyl-CoA and acyl carrier protein (ACP) to form 3-carbon malonyl-ACP is catalysed by malony CoA acyl transferase. Polymerization of the fatty acid carbon chains starts with the formation of a 4-carbon butyryl-ACP by the addition of a carbon from acetyl-CoA to the 3-carbon malonyl-ACP by the action of keto acyl ACP synthase III (KASIII). Further chain elongation up to 16-cabon palmitoyl-ACP is carried out by keto acyl ACP synthase I (KASI) in six consecutive steps by adding two carbons in each step. Finally, keto acyl ACP synthase II (KASII) converts the 16-cabon palmitoyl-ACP to 18-carbon stearoyl-ACP. While KASI is essential for seed oil content [30], KASII is useful for manipulating the chain length of the fatty acids. Longer the chain length of fatty acids more is the viscosity of the oil as well as the biodiesel derived from it [31]. Viscosity affects atomization of the biodiesel and lowers the performance of the diesel engines due to engine deposits [32]. Viscosity of jatropha oil can be lowered by reducing its 18-carbon fatty acid content which is reported to be as high as 84.7% [33]. Previous reports suggest that this can be achieved by silencing the activity of KASII that converts the palmitoyl-ACP to stearoyl-ACP [34,35]. Alternatively, the activity of palmitoyl-ACP thioesterase (encoded by FATA1) could be enhanced to accelerate the cleavage of palmitic acid from the palmitoyl-ACP before it is converted to stearoyl-ACP [36]. The contigs for all these genes were found in the current study. and omega-3-fatty acid desaturase (FAD3), respectively. Increased degree of unsaturation affects oxidative stability and ignition quality of the biodiesel [32]. Therefore, the high content of unsaturated fatty acids (78.4%) in jatropha oil is not desirable for biodiesel production. Knutzon et al., 1992 [37] have significantly reduced the

Activation and transport of fatty acids
Free fatty acids that are released from acyl-ACPs are converted to respective acyl-CoAs (activation) by long chain acyl CoA synthetase (LCACS). The activated fatty acids are bound by acyl CoA binding proteins (ACBPs) to protect them from acyl-CoA hydrolases, and to transport them to endoplasmic reticulum. Arabidopsis contains nine LCACS genes [38] and six ACBP genes [39]. From the current study, we have identified four contigs for LCACS, and three of them were full-length. We also found five contigs for ACBPs, and two of them were full-length.

Synthesis of triacylglycerol or oil
Triacyl glycerol or oil is synthesised in ER by serial incorporation of three acyl groups to the glycerol backbone. The first acyl group is added to the sn-1 position by glycerol-3-phosphate acyl transferase (GPAT). The resulting lysophosphatidic acid is further acylated at sn-2 position by lysophosphatidic acid acyl transferase (LPAT) to form phosphatidic acid. Phosphatidic acid is dephosphorylated by phosphatidic acid phosphatase (PAP) to form diacyl glycerol. Final acylation of diacyl glycerol in sn-3 position to complete the synthesis of triacylglycerol is carried out by diacyl glycerol acyl transferase (DGAT). These steps are highly amenable for genetic manipulations to increase the oil content. For example, seed oil content was increased by 21% in arabidopsis and 3-7% in canola by overexpressing GPAT [40] and DGAT [41], respectively. We have identified 12 contigs for the four enzymes, GPAT (4 contigs). LPAT (5 contigs), PAP (1 contig) and DGAT (2 contigs) that are required for triacyl glycerol synthesis. These sequences will be useful for increasing the oil content in jatropha by genetic engineering.

Functional classification
Gene ontology (GO) annotation and classification of the sequences was done based on three gene ontology terms such as biological processes, cellular components and molecular functions. Distribution of the contigs in different GO categories is given in figure 4. While the contigs covered all the major biological processes, the cellular and metabolic processes alone accounted for 52.7% of the contigs. Regulation of biological processes and response to stimulus accounted for 11.6% and 8.9% of the contigs, respectively. In cellular components, cell and organelles accounted for the maximum number of contigs (60%). In molecular functions, highest numbers of contigs were classified under catalytic activity (36.77%) which was followed by transferase activity, protein binding, hydrolase activity and nucleic acid binding activity (12.33 to 14.27%).

Validation of assembled transcripts
The 3' untranslated region alone was RT-PCR amplified for 11 transcripts to validate the expression of assembled transcripts. All the transcripts were found to be expressed in all the five tissues studied ( Figure 5). Whole transcripts were RT-PCR amplified for 17 transcripts to validate the expression and contig assembly. As shown in figure 6, all the transcripts were amplified successfully and were of expected size (500 to 5,000 bp approximately). In addition, in silico validation was done by comparing the assembled contig sequences with 7,009 high quality and manually edited ESTs and full-length cDNA sequences (100 and 2,372 bases) obtained from Sanger sequencing. We found 550 contigs to be fully mapping with the nucleotide sequence from Sanger sequencing, and the identity was 99.71% over a length of 100 to 1,932 bases. When 2,581 partially mapped contigs were analysed, the identity was 99.90% over a length of 37 to 1,900 bases. Overall, 1,253,655 bases were aligned with 99.87% identity. These results not only validate the accuracy of 454 pyrosequencing vis-à-vis Sanger sequencing but also the accuracy of contig assembly.

Conclusions
We have analyzed the transcriptome from five major tissues of J. curcas using GS FLX titanium platform of 454 pyrosequencing for the purpose of large scale gene discovery. We assembled 17,457 contigs and identified 14,327 new transcripts to add to the existing jatropha genetic resource. Currently, 28,794 non-redundant transcripts sequences are available from jatropha. This includes 56 transcripts that are directly involved in oil biosynthesis. These sequences will be useful in genetic engineering of jatropha for increased oil content, modified oil composition, toxin-free oil cake, resistance to insect pests and diseases, tolerance to drought stress etc.

Collection of tissues for RNA extraction
We have used Jatropha curcas L. for our study. Tissue from roots was collected from two months old plants grown in a green house at 27 ± 2°C and 90% RH. Tissues from mature leaves, flowers, developing seeds, and embryo from mature dry seeds were collected from a jatropha plantation in Tamil Nadu, India. The tissues were immediately frozen in liquid nitrogen and stored at -86°C.

RNA extraction
Total RNA from developing seeds was prepared using Trizol reagent as described before [10]. Total RNA from other tissues was isolated by following a protocol described by Singh et al., 2009 [42] except that the pellet obtained after lithium chloride precipitation was washed with 75% ethanol twice, and purified using RNeasy Mini Kit following the manufacturer's protocol (Qiagen, Hilden, Germany). Agarose gel electrophoresis and D 260 / OD 280 ratio were used for assessing quality of total RNA.

Synthesis of normalized cDNA
Pooled total RNA of 100 μg was prepared by mixing 5 μg, 25 μg, 15 μg, 30 μg, and 25 μg of total RNA from roots, mature leaves, flowers, developing seeds, and embryo, respectively. About 30 μg of the pooled total RNA was used for cDNAs synthesis and normalization was carried out as described by Patanjali et al., 1991 [43] and Soares et al., 1994 [44] with slight modifications. Tester3 primer (CAGTGGTATCAACGCA-GAGTGGCCGAGGCGGCCT 15 ) and driver3 primer (GGGATAACAGGGTAATGGCCGAGGCGGCCGA-CATGT 15 ) were used to prime first strand tester cDNA and driver cDNA synthesis, respectively. Tester adaptor (GTAACTAGGCCGTAATGGCCACTCTGCGTTGA-TACCACTG) and driver adaptor (GGCCGTAATGG CCTCGCTACCTTAGGA) were ligated to the 3' end of the first strand tester cDNA and driver cDNA,  respectively. These adaptors were 5' phosphorylated and 3' blocked to prevent their ligation to the 5'end. Double strand tester cDNA was synthesised using tester3 primer and tester5 primer (CAGTGGTATCAACGCA-GAGTGGCCATTACGGCCTAGTTACGGG). The tes-ter5 primer is complementary to the 5' tester adaptor and it is phosphorylated. The sense strand of the double strand tester cDNA which is phosphorylated at the 5' end was destroyed by treating it with lambda exonuclease. As a result, only the antisense strands are retained. Double strand driver cDNA was synthesised using driver3 primer which is phosphorylated and dri-ver5 primer (TCCTAAGGTAGCGAGGCCAT-TACGGCCGGG) which is complementary to the 5' driver adaptor. The anti-sense strand of the double strand driver cDNA which is phosphorylated at the 5' end was destroyed by treating it with lambda exonuclease. As a result, only the sense strands are retained.
Hybridization of the anti-sense strands from tester cDNA and sense strands from driver cDNA was carried out in 1× hybridization buffer (50 mM Tris-HCl pH 8.0, 0.5 M NaCl, 0.2 mM EDTA) at 68°C for 6 h. Formation of double strand hybrids depends on the reassociation kinetics (second order kinetics), in which more abundant sequences anneal faster than rare sequences. Therefore, the double strand hybrids were removed by hydroxyapatite chromatography to achieve normalization [44]. The normalized single strand cDNAs were converted to double strand cDNAs and were amplified by Failsafe™ PCR system (Epicentre Biotechnologies, USA) using single amplification primer (CAGTGGTAT-CAACGCAGAGT) for which the binding sites are present only in the tester cDNA (underlined in tester3 and tester5 primers used for cDNA synthesis). Normalization efficiency was monitored by performing a parallel normalization in which chloramphenicol resistance gene  was added to the cDNA at 1.0% redundant rate as internal control. Normalized cDNAs were cloned in modified pBluescript II SKand normalization efficiency was determined by plating the cDNA library on LB plates containing chloramphenicol.

Quality Control of cDNA
An aliquot of normalized cDNA was used to construct a cDNA library in modified pBluescript II SKfor quality control purpose. Ten colonies from the cDNA library were randomly selected and insert size was determined by colony PCR using M13 forward (GTAAAACGACGGC-CAGT) and M13 reverse primer (CAGGAAACAGCTAT-GAC). Another 20 clones were randomly selected and sequenced from the 5' end of the cDNA using M13 reverse primer and BigDye™ Terminator v3.1 Cycle Sequencing Kit in 3130xl Genetic Analyzer (Applied Biosystems, CA, USA). These sequences were annotated by using BLASTX algorithm and non-redundant database at NCBI [45]. Full-length nature of the clones was determined from the BLASTX result as described before [10].

pyrosequencing
The normalized cDNA was quantified by PicoGreen assay and 3.8 μg of cDNA was used for making fragmentation library using GS FLX Titanium General Library Preparation Kit (Roche 454 Company, CT, USA) as described in the manufacturer's manual. Half-plate reaction of 454 pyrosequencing was done using 454 Genome Sequencer FLX System (Roche, CT, USA).

De novo assembly
The 454 Genome Sequencer FLX System collects the data and generates standard flow gram file (.sff) which contains raw data for all the reads. The raw data was quality-filtered using a quality cut-off value of 40. The primer and adapter sequences that were incorporated during cDNA synthesis and normalization were removed. Sequences with less than 50bp were removed before contig assembly. De Novo contig assembly of the reads was performed using GS De Novo Assembler software versions, v2.3, v2.5p1 and v2.5.2 which were provided by 454 Life Sciences Corp, CT, USA. The assembly parameters used were minimum overlap length of 40 bp and minimum overlap identity of 95%.

Pathway mapping using KEGG
Gene ortholog assignment and pathway mapping of the contigs was done using KEGG (Kyoto Encyclopaedia of Genes and Genomes) automatic annotation server [46]. The contigs were assigned with the unique enzyme commission (EC) numbers based on the similarity hit against KEGG database using BLASTX (default threshold bit-score value of 60) . Distribution of contigs under the respective EC numbers was used to map them to the KEGG biochemical pathways.

Annotation and functional classification
Annotation of the assembled transcript sequences was performed using BLASTX algorithm and non-redundant Table 4 Name of the genes, primers used and expected size of the RT-PCR products for amplification of full-length transcripts protein database at NCBI [45]. The BLASTX results were also used to assess the full-length nature of the contigs. The automated BLASTX analysis was done using BLAST2GO to assign GO terms for the contigs [47]. The transcripts were classified under three GO terms such as molecular function, cellular process and biological process.

Semi-quantitative RT-PCR
For RT-PCR, total RNA from roots, mature leaves, flowers, developing seeds, and embryo were treated with DNase I and purified using RNeasy Mini Kit following the manufacturer's protocol (Qiagen, Hilden, Germany). About 3.0 μg of purified total RNA from each sample was used for first strand cDNA synthesis using oligo-dT (18) primer and PrimeScript™ reverse transcriptase (Takara Bio Inc, Shiga, Japan). Equal quantity of first strand cDNA (from 25 ng total RNA) was used for PCR. Primers were designed to amplify the 3' UTR of the transcripts. Actin gene was used as an internal control. Primer sequences and the expected size of the amplified fragments are given in table 3. Semi-quantitative analysis of the RT-PCR amplified fragments was done by agarose gel electrophoresis.

RT-PCR amplification of full-length transcripts
First strand cDNA from 3.0 μg of pooled total RNA from the five tissues was prepared as described above, and entire length of 17 transcripts was amplified by PCR. Contig number, transcript (gene) name and primer sequences for RT-PCR are given in table 4. Size of the amplified fragments was determined by agarose gel electrophoresis using DNA markers.

Accession number
All the reads generated and used for this study were submitted to the sequence read archive (SRA) at NCBI with the accession number SRP004898.