GENE-39674; No. of pages: 8; 4C: Gene xxx (2014) xxx–xxx

Contents lists available at ScienceDirect

Gene journal homepage: www.elsevier.com/locate/gene

F

5 6

a

7

a r t i c l e

8 9 10 11 12

Article history: Received 7 January 2014 Received in revised form 26 April 2014 Accepted 7 May 2014 Available online xxxx

13 14 15 16 17

Keywords: Gene regulation Adventitious root MeJA Panax ginseng

Graduate School of Biotechnology & Ginseng Bank, College of Life Science, Kyung Hee University, Yongin 449-701, South Korea Insilicogen Inc., #909, Venture Valley, 958, Gosaek-dong, Gwonseon-gu, Suwon, Gyeonggi-do 441-813, South Korea

i n f o

a b s t r a c t

Panax ginseng Meyer is one of the major medicinal plants in oriental countries belonging to the Araliaceae family which are the primary source for ginsenosides. However, very few genes were characterized for ginsenoside pathway, due to the limited genome information. Through this study, we obtained a comprehensive transcriptome from adventitious roots, which were treated with methyl jasmonic acids for different time points (control, 2 h, 6 h, 12 h, and 24 h) and sequenced by RNA 454 pyrosequencing technology. Reference transcriptome 39,304,529 (0.04 GB) was obtained from 5,724,987,880 bases (5.7 GB) of 22 libraries by de novo assembly and 35,266 (58.5%) transcripts were annotated with biological schemas (GO and KEGG). The digital gene expression patterns were obtained from in vitro grown adventitious root sequences which mapped to reference, from that, 3813 (6.3%) unique transcripts were involved in ≥2 fold up and downregulations. Finally, candidates for ginsenoside pathway genes were predicted from observed expression patterns. Among them, 30 transcription factors, 20 cytochromes, and 11 glycosyl transferases were predicted as ginsenoside candidates. These data can remarkably expand the existing transcriptome resources of Panax, especially to predict existence of gene networks in P. ginseng. The entity of the data provides a valuable platform to reveal more on secondary metabolism and abiotic stresses from P. ginseng in vitro grown adventitious roots. © 2014 Published by Elsevier B.V.

E

T

C

b

36 34 33

R

35

1. Introduction

38

Panax ginseng Meyer (Korean ginseng) is a medicinal plant known as adaptogen and a major source for triterpenoidal drugs (ginsenosides) from nature. In addition, it has been used worldwide to cure almost all human diseases from ancient times to now, particularly in Asian countries such as Korea (Oriental medicine) and China (Chinese medicine). Ginseng plant materials are available in the form of fresh roots, red ginseng (steamed and dried roots), dried leaves, dried flower buds, green tea, tablets, and capsules (Asaka et al., 1993). From these different materials, 128 ginsenoside isoforms (ginseng saponins) were isolated and modified (Mathiyalagan et al., 2014) to act as therapeutic components

47

N C O

45 46

U

43 44

R

37

41 42

Q4 18 19 20 21 22 23 24 25 26 27 28 29 30 Q5 31

E

32

39 40

O

4

Sathiyamoorthy Subramaniyam a,b, Ramya Mathiyalagan a, Sathishkumar Natarajan a, Yu-Jin Kim a, Moon-gi Jang a, Jun-Hyung Park b, Deok Chun Yang a,⁎

R O

3Q1

P

2

Transcript expression profiling for adventitious roots of Panax ginseng Meyer

D

1

Abbreviations: AACT, acetoacetyl-CoA thiolase; HMGS, 3-hydroxy-3-methylglutarylCoA (HMG-CoA) synthase; HMGR, 3-hydroxy-3-methylglutaryl-CoA reductase; MVA, mevalonate; MK, mevalonate kinase; MVAP, mevalonate-5P; FPS, farnesyl diphosphate (FPP, C15) synthase; GGPS, geranylgeranyl diphosphate (GGPP, C20) synthase; SS, squalene synthase; SE, squalene epoxidase; BAS, beta-amyrin synthase; CAS, cycloartenol synthase; DS, dammarenediol synthase; GT, glycosyltransferase; CYP, cytochrome; PPT, protopanaxatriol; PPD, protopanaxadiol. ⁎ Corresponding author at: Graduate School of Biotechnology & Ginseng Bank, College of Life Science, Kyung Hee University, 1 Seocheon, Giheung-gu Yongin-si, Gyeonggi-do 449-701, South Korea. E-mail address: [email protected] (D.C. Yang).

and these are only present in Panax plant materials. Ginsenosides were classified into protopanaxadiol, protopanaxatriol, and oleanane types and further sub-classified into major and minor ginsenosides with respect to type/position of the sugar molecules (Choi, 2008). Every single ginsenoside has its own pharmacological activities for treating cancer, diabetes, obesity, osteoporosis, Alzheimer's, Parkinson's disease, HIV, and sexual dysfunctions. Moreover, ginseng also contains other antioxidant elements, including polysaccharides, phenolic compounds, flavonoids, peptides, polyacetylenic alcohols, fatty acids, and gintonins (Hwang et al., 2012; Shi et al., 2007). (See Table 1.) Different types of ginseng roots i.e., wild ginseng (grown more than 100 years in deep forest), field ginseng (4–15 year cultivation in field), and in vitro plant biomaterials (adventitious roots and embryogenic calli) which are produced from in vitro tissue culture systems using large bioreactors (Sathiyamoorthy et al., 2011; Wu and Zhong, 1999), and ginsenoside contents also differ based on root types. Among those, in vitro plant materials were used in our study, because the plant signaling molecules were used to enhance the ginsenoside productivity in such tissue culture systems. Elicitors such as ethylene, chitosan (Jeong and Park, 2005), salicylic acid (Jeong et al., 2005), jasmonic acid derivatives (Ali et al., 2006b; Bae et al., 2006; Hu and Zhong, 2007, 2008; Lu et al., 2001; Qian et al., 2004), copper (Ali et al., 2006a), germanium (Yu et al., 2005), NiSO4, sodium chloride, yeast, ethephon, and nitric oxide (Hu

http://dx.doi.org/10.1016/j.gene.2014.05.024 0378-1119/© 2014 Published by Elsevier B.V.

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

48 49 50 51 52 53 54 55 56 57 Q6 58 59 60 61 62 63 64 65 66 67 68 69 70

Reads Assembled contigs (CLC Assembly Cell) Reference transcriptome (CAP3)

t1:9 t1:10 t1:11 t1:12 t1:13 t1:14 t1:15

Sequence with blast hits Sequence annotated with (BLAST, GO and KEGG) Sequence Annotated with (BLAST, GO) Sequence without GO mapping Sequence without blast hits

t1:16 t1:17 t1:18

Number of gene ontologies Number of KEGG maps Number of unique UniprotKB sequences

534,324 71,256 60,236

B. Annotation 35,266 (58.5) 3391 (5.6) 22,432 (37.2) 12,834 (21.3) 24,970 (41.5)

123

Individual samples were taken for total RNA extraction using TRIzol® Reagent kit according to the manufacturer's instructions. mRNA was isolated from total RNA using Dynabeads oligo(dT) 25 kit according to the manufacturer's instruction. cDNA was prepared from total mRNA, using cDNA rapid library preparation method according to the manufacturer's instructions (http://www.roche-applied-science. com/) (Supplementary Table 1). Finally, total cDNA was used for sequencing on a GS FLX sequencer (454/Roche) and obtained sequences were added to de novo assembly. De novo assembly was performed to produce reference transcriptome for P. ginseng from available high throughput sequences in NCBI SRA as public and sequence from Ginseng Genetic Resource Bank (Supplementary Table 2). De novo hybrid assemblies were performed by CLC Assembly Cell v.4.0 (http://www. clcbio.com) with customized parameters. First the whole Illumina sequences were subjected to preprocessing steps with default parameters. Preprocessing steps which are quality trim (Q ≥ 20), and adapter trim (universal adapters) were performed with CLC Assembly Cell v.4.0, and 454 pyro-sequences were processed with SFF-extract in MIRA assembler (http://www.chevreux.org/projects_mira.html). Total, preprocessed sequences were assembled with CLC assembler using different word sizes (24–51), and bubble sizes (300–500) to optimize assembly (Supplementary Table 3). Among the assembled contig sets, the best one was selected as a reference transcriptome of P. ginseng and it was reassembled with CAP3 using the default parameters to remove the isoforms (Huang and Madan, 1999). Resulted singletons and contigs were merged together and renamed with unique transcript name. Sequences obtained for this study were submitted to NCBI SRA (SRP039367).

124

2.3. Functional annotations

152

Reference transcriptome was subjected to functional annotation using BLASTX mapping (e-value cut-off 1e− 10) against UniprotKB (viridiplanta) database. Descriptions were obtained from the highest blast score among the mapped sequences for each P. ginseng transcript. Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway map ids were transferred to P. ginseng transcripts from Uniprot database with respect to BLAST sequences. These analyses were performed by python scripts.

153

2.4. Digital gene expression (DGE) profiling

161

Final assembled and renamed unique transcripts were used as a reference for reference mapping. The CLC mapper module from CLC Assembly Cell v.4.0 was used for reference mapping. Total preprocessed reads from adventitious root treated with MeJA (control, 2 h, 6 h, 12 h, and 24 h) were mapped to reference individually with 90% similarity and 80% length fraction respectively. Individual read counts for the reference unique transcripts were obtained from CLC-mapper output files by python scripts and transcripts with read counts N5 were selected for further normalization steps. Total samples were grouped in to three as control, first (2 h, 6 h) and second (12 h, 24 h) to obtain the comparable results. Later, read counts were normalized with Audic statistical significant method (Audic and Claverie, 1997) to obtain the significance expression and normalized values were taken to Gene Spring 12.5 GX for fold change analysis with user created library options with

162

F

C. Data base

2.2. Pyro-sequencing and de novo assembly

4204 968 26,917

71

U

N

C

O

R

R

E

C

T

et al., 2003; Jeong and Park, 2006) were optimized. Among these, methyl 72 jasmonate (MeJA) resulted in high ginsenoside contents. Moreover, MeJA 73 is a well-known elicitor to regulate the secondary metabolism in plants 74 (Pauwels et al., 2009) and reported in enhanced secondary metabolites 75 such as terpenoids, flavonoids, alkaloids, and phenylpropenoids. In recent 76 years, ginseng plant genetic engineering studies have turned towards to 77 modulate the secondary metabolites, particularly the triterpenoid 78 (ginsenoside) pathway (Han et al., 2010), defense-related genes 79 (Ok et al., 2011), and biotic, and abiotic stress (Sathiyaraj et al., 2011). 80 Notably, adventitious roots were mainly used as a tool for gene character81 ization and genetic engineering studies due to their physiology and long 82 growth cycle of P. ginseng (Shim et al., 2010). Based on these evidences, 83 adventitious roots with MeJA were used for this study. 84 Moreover, to understand the genetic networks and dynamic present 85 in P. ginseng, and to reveal the answer for the numerous biological ques86 tions such as, how it maintain the Panax specific secondary metabolites 87 (ginsenosides) and genes which are responsible for different types of 88 ginsenosides are more puzzling. To get a deep understanding in the ex89 istence of genetic networks of Panax, a comprehensive investigation is 90 Q7 needed. To date, very few sequencing studies for P. ginseng have been 91 conducted, i.e. Sanger sequencing for different organs, different age 92 group roots, embryogenic callus, and adventitious roots (Choi et al., 93 2005; Sathiyamoorthy et al., 2010a,b, 2011) and next generation se94 quencing technique (NGS) RNA-seq for 11-year-old ginseng root 95 (Chen et al., 2011; Li et al., 2013), microRNA profiling for plant organs 96 (Li et al., 2013; Mathiyalagan et al., 2013; Wu et al., 2012). The plant or97 gans of other species, such as Panax quinquefolius (Sun et al., 2010; Wu 98 Q8 et al., 2013), and Panax notoginseng (Luo et al., 2011), were sequenced. 99 Quantity of ginseng transcripts was not sufficient to support further 100 functional genome studies with respect to the dynamics of genetic net101 Q9 works, particularly for adventitious roots cultured from ginseng roots. 102 The goal of our present study is to investigate gene expression profiles 103 from MeJA treated adventitious root cultures of ginseng using NGS 104 RNA sequencing technology. Consequently, the transcripts were 105 generated by Roche's 454 GS FLX to obtain continuous data sets for 106 time-dependent samples. In silico bioinformatics analyzed was used to 107 ascertain the dynamic of genetic molecular networks and to annotate 108 Q10 transcripts with standard annotation pipeline. Finally, the candidate 109 genes were predicted for ginsenoside pathway. To our knowledge, this 110 is the first study used to examine gene expression patterns from highly 111 expressed transcripts of P. ginseng using RNA-seq.

O

A. Sequence & assembly

t1:4 t1:5 t1:6 t1:7 t1:8

116 117

R O

t1:3

cv. “Yunpoong” root explants (provided by Ginseng Genetic Resource Bank). Roots were cultured on 1/2 MS medium with 2% sucrose, 0.25% gelrite, and 10 μM indolebutyric acid (IBA). After induction adventitious roots were excised from the maternal explants prior to sub-culturing, which was replaced every 5 weeks with the same medium. Adventitious roots at 4 weeks of culture were elicited by adding MeJA (50 μM) for different time periods (control, 2, 6, 12, and 24 h).

P

Q3

Table 1 Sequence and annotation summary for in vitro cultured adventitious roots.

D

t1:1 t1:2

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

E

2

112

2. Materials and methods

113

2.1. Adventitious root culture and MeJA treatment

114

In this study, we used in vitro cultured 4-week-old adventitious roots for all conditions, which are induced from 6 year old Panax ginseng

115

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

118 119 120 121 122

125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151

154 155 156 157 158 159 160

163 164 165 166 167 168 169 170 171 172 173 174 175

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

198 199 200 201 202 203 204

F

O

207 208 209

213 214 215 216 217 218

3.1. Sequencing and transcriptome assembly

220

cDNAs prepared from MeJA treated and control samples of adventitious roots which cultured from P. ginseng roots were sequenced using Roche 454 titanium FLX. As a result, 534,324 raw reads were obtained totally from all samples. Sequences were preprocessed along with other collected libraries (which are available in public biological database) to prepare reference transcriptome for P. ginseng. Totally, 22 libraries (Supplementary Table 2) with 59,828,590 reads include 5,724,987,880 (5.7 GB) bases were taken for preprocessing, which also includes methyl jasmonic acid treated samples. Initially, total sequences were subjected to preprocess as mentioned in the Materials and

221

C

196 197

219

E

194 195

3. Results and discussion

R

192 193

R

190 191

211 212

N C O

188 189

The individual samples were subjected to saponin extraction with nbutanol saturated with H2O, and evaporated in a vacuum. The residues were dissolved in methanol and analyzed by HPLC using a C18 (250 × 4.6 mm, particle size 5 μm) column with acetonitrile (solvent A) and distilled water (solvent B) as mobile phases at 85% B for 5 min, 79% B for 20 min, 42% B for 55 min, 10% B for 12 min, and 85% B for 18 min all at 1.6 mL min−1. Detection of ginsenosides was at 203 nm (Quan et al., 2012).

U

186 187

210

R O

Total RNA was extracted from frozen samples with the RNeasy plant mini kit (Qiagen, USA), including the DNase I digestion step. Next, 2 μg of the total RNA was reverse transcribed with RevertAid™ H Minus M-MuLV Reverse Transcriptase (Fermentas, USA). RT-PCR was performed in a 25 μl reaction volume consisting of 0.2 μl of the cDNA product and 5 pmol of each primer using Super-Taq DNA polymerase (Super Bio, Korea) by a Bio-Rad PCR machine (3 min at 95 °C, followed by 28 cycles at 94 °C for 20 s, 60 °C for 20 s, and 72 °C for 30, with a final extension at 72 °C for 5 min). The products were analyzed on 1.2% agarose gels. Real-time quantitative PCR was performed using 100 ng of cDNA in a 10 μL reaction volume using SYBR® 28 Green Sensimix Plus Master Mix (Quantace, Watford, England). The following thermal cycler conditions recommended by the manufacturer were used: 10 min at 95 °C, followed by 40 cycles at 95 °C for 10 s, 60 °C for 10 s, and 72 °C for 20 s. The fluorescent product was detected during the final step of each cycle. Amplification, detection, and data analysis were carried out on a Rotor-Gene 6000 real-time rotary analyzer (Corbett Life Science, Sydney, Australia). To determine the relative fold-differences in template abundance for each sample, the Ct value for each of the gene-specific genes was normalized to the Ct value for β-actin, and was calculated relative to a calibrator using the formula

2.6. Ginsenoside analysis using HPLC

P

184 185

205 206

D

2.5. RNA isolation and qRT-PCR

180 181

2−ΔCt or 2−ΔΔCt. Three independent experiments were performed, and the primer efficiencies were determined according to the method of Livak and Schmittgen (2001) in order to validate the ΔΔCt method. The observed slopes were close to zero, indicating that the efficiencies of the gene and the internal control β-actin were equal.

E

183

178 179

T

182

default parameters. Fold changes were calculated for control vs first and control vs second respectively. Finally ≥2 fold up and downregulated were grouped into four clusters (UpUp, DownDown, UpDown, DownUp) to obtain detail expression patterns using python scripts. Finally, co-expressions were calculated with hierarchical clustering with the default parameters using Cluster3 software and trees were drawn by Java Treeview (de Hoon et al., 2004).

176 177

3

Fig. 1. The Gene Ontology (GO) classifications for ≥2 fold regulated transcripts with respect to biological process (BP), molecular function (MF), and cellular components (CC).

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

222 223 224 Q11 225 226 227 228 229 230

237 238 239

243

259

Total sequenced reads from five adventitious roots were mapped to reference transcriptome individually to obtain the expression patterns and further its combination into three groups such as control, first (2 h, 6 h) and second (12 h, 24 h). From now sequences were mentioned with group names. Mapped results are 98,363 (79.5%), 139,801 (79.8%) and 188,119 (80%) reads for control, first and second respectively. DGE pattern for each transcript was calculated from the observed read counts as in the Materials and methods section. Totally 2264 (3.7%) sequences from the first and 2564 (4.2%) sequences from the second groups were in ≥2 fold up and downregulations. These regulated sequences were further grouped into four clusters which resulted in 1069 (1.8%), 2220 (3.7%), 308 (0.5%), and 215 (0.4%) for UpUp, DownDown, UpDown, and DownUp respectively (Supplementary Fig. 1). The results showed that most of the genes were downregulated and few genes were upregulated upon elicitations. These differentially regulated genes were further grouped into biological process, molecular functions and cellular components. In general, plant suspension cultures exposed to external elicitation, genes in plasma membrane were mostly regulated. The detail mechanisms were reviewed (Zhao et al., 2005). From our results, mostly transcripts grouped in cellular components are in plasma membrane (GO:0005886), which might activate

260

U

N

C

O

R

R

E

C

T

E

Reference transcripts were subjected to BLASTX against plant 244 UniprotKB database to obtain the biological descriptions. The 35,266 245 (58.5%) sequences were annotated with minimum one biological term 246 from Gene Ontology (GO) or KEGG pathway information. In that, 3391 247 (5.6%) sequences were completely annotated with GO and KEGG, and 248 22,432 (37.2%) sequences were annotated with only GO term, and 249 12,834 sequences were not annotated with both and 24,970 (41.5%) se250 Q12 quences did not obtain any blast similarity within the given blast cut off. 251 So these 24,970 sequences were considered as P. ginseng genome specif252 ic transcripts (Table 1B). The 26,117 unique genes from UniprotKB data253 base, and 4204 functional terms from GOs, and 968 KEGG pathways 254 (Table 1C) were transferred to P. ginseng transcriptome. Sharing knowl255 edge of known functions of these 26,117 genes from other plant

3.3. Digital gene expression

F

3.2. Functional annotation

235 236

O

242

233 234

genomes aids for comparative transcriptome and to characterize 256 P. ginseng. Moreover, annotations shared from GO and KEGG will sup- 257 port for more detailed functional characterizations. 258

R O

240 241

methods section, resulting with 56,186,312 (94%) reads with 4,982,295,603 (87%) bases. Preprocessed (4.9 GB) sequences were taken for de novo assembly using CLC Assembly Cell v.4.0. Total 71,256 contigs (44,087,537 bases) were selected as the best assembly (bubble size 500 and word size 24) and isoforms were removed using CAP3. Finally 60,236 (39,304,529 bases) unique sequences were obtained ranging from 200 to 13,903 with an average length of 652.5 bases considered as a unique reference draft transcriptome (0.04 GB) for P. ginseng (Table 1A) (Supplementary Table 4). To our knowledge, this is the first large scale (4.9 GB) de novo assembly performed for P. ginseng from in house and public data sets.

P

231 232

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

D

4

Q2

Fig. 2. Illustrations of external methyl jasmonic acid treatment in plant cells, and which probably activate the putative ginsenoside pathway genes. Initially, MeJA binds to membrane receptors and activates the G-proteins to trigger the phospholipase A (PLA), later activate α-linolenic acid and endo-methyl jasmonate from PLA. Endogenesis MeJA regulates the HMGR (Pollier et al., 2013) and regulating the down-stream genes to produce the mono-, di-, sesqi-, and tri-terpenoids genes, although ginsenoside pathway (Zhao et al., 2013) and another way oxylipins activate the transcription factors to regulate the defense and secondary metabolite genes (De Geyter et al., 2012; Zhao et al., 2005).

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

261 262 263 Q13 264 265 266 267 Q14 268 269 270 271 272 273 274 275 276 277 278 279 280

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

294

3.4. Expressions of transcription factors

295

327

Particularly, plant materials elicited by MeJA derivatives, the intracellular metabolism activated through G-proteins present in inner membrane to phospholipase and further, the alpha-linoleic pathway will trigger the inner methyl jasmonate pathway genes. From the continuation, signals transferred via cyclic oxylipins to nucleus and activate transcription factors (Fig. 2) and transcription factors are the key switches for secondary metabolite gene regulation (De Geyter et al., 2012). As observed from DGE results, gene modulation takes part mainly in the plasma membrane (GO:0005886) and then the nucleus (GO:0005634) which mimics the general elicitation mechanism. Totally, 181 (0.3%) sequences were putatively annotated for transcription factors which are involved in ≥2 fold up and downregulations which excludes non-informative sequences like putative and uncharacterized proteins (Fig. 3). In result, transcription factor families like WRKY, MYB, AP2, NAP, G/HBF-1 and ERF are mostly upregulated in treated groups than control and GRAS, DELLA, and AP2 are highly expressed in first group samples. Among these, AP2/ERF (PG027192) in Catharanthus roseus (Sawai and Saito, 2011), WRKY22 (PG004611) in C .roseus (Suttipanta et al., 2011) and WRKY1 (PG003175) in P. quinquefolius (Sun et al., 2013) are reported for terpenoids and indole alkaloids. MYB (PG025883) characterized for proanthocyanidin in Diospyros kaki (Akagi et al., 2009) and DELLA for gibberellin in Malus domestica (Foster et al., 2007) and ERF3 (PG002516), WRKY4 (PG015962) and G/HBF-1 (PG016838) are responsible for stress/ defense inducible transcription factors. Moreover, microRNAs were also predicted and validated from P. ginseng, for transcription factors such as GRAS, AP2, and NAC (Wu et al., 2012). These putative transcription factors in P. ginseng transcriptome are highly homologous to above characterized transcription factors and noted as DGE ≥ 2 fold up and downregulations. So, these are considered a regulator for ginsenoside biosynthetic genes. Based on our prediction, these are the candidates for detail characterization to reveal the relationship between transcription factors and ginsenoside biosynthesis.

328

3.5. Expression of secondary metabolite genes

303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326

329

O R O P

301 302

D

300

T

298 299

C

296 297

E

290 291

R

288 289

R

287

N C O

285 286

U

283 284

F

292 293

multifunctional role in energy metabolism, transcription factors, and DNA binding proteins in nucleus (GO:0005634) (Fig. 1). Further, more detail GO analyses for the four specific clusters (UpUp, DownDown, UpDown, DownUp) were grouped with respect to cellular components, molecular functions, and biological process (Supplementary Fig. 2). In continuous up-regulations (UpUp) most transcripts involved in molecular functions such as electron carrier activity (GO:0009055), heme binding (GO:0020037), and iron ion binding (GO:0005506), and biological process such as lipid metabolic process (GO:0006629), and metabolic process (GO:0008152) were higher than in continuous downregulations (DownDown). Whereas other functions in molecular functions and biological process (as denoted in GO) are dominated in continuous downregulation patterns (DownDown).

E

281 282

5

Primary and secondary metabolites from P. ginseng are ginsenosides 330 Q15 which belong to tri-terpenoid (C30) components. Terpenes such as 331 mono-terpenes (C10), sesqui-terpenes (C15), and di-terpenes (C20) 332 and tri-terpenoids are derived from terpenoid backbone pathway 333 (referred from KEGG pathway MAP0900) (Sawai and Saito, 2011). 334 Ginsenoside biosynthesis was initiated from oxidosqualene to putative 335 derived ginsenoside pathway (Fig. 1). Based on this, genes responsible 336 for the secondary metabolites were clustered. Totally 152 (0.25%) tran337 scripts were found in ≥2 fold up and downregulations (Supplementary 338 Fig. 3) and to assess and correlate the transcript expression and 7 major 339 ginsenoside contents, ginsenosides were quantified by HPLC from same 340 Q16 root samples, which are used for sequencing. Ginsenoside contents re341 sulted with increased in total content from 3.27 ± 0.11 mg/mL 342 to 3.836 ± 0.07 mg/mL (Supplementary Table 5). In specific, PPD

Fig. 3. Hierarchical cluster and co-expression patterns for transcription factors which are involved in ≥2 fold regulations.

(0.07 ± 0.01) was mostly synthesized than PPT (0.02 ± 0.02) (Supplementary Fig. 6). Ginsenosides differ according to growth stages and in different parts of plants (Kim et al., 2014a; Shi et al., 2007). Further, the gene expression patterns which observed from in silico analysis were compared with literature evidence and validated with qRT-PCR. Totally, 7 transcription factors and 5 ginsenoside responsible genes were validated with qRT-PCR, which shows the similar expression patterns with in silico expression and literature (Supplementary Fig. 4). In in silico result showed that tri-terpenoid backbone and putative ginsenoside pathway responsible transcripts like HMGR, HMGS, FPP, HDR, SE1, SE2, and SS, were also continuously up-regulated in MeJA treated samples than the control (Kim et al., 2010, 2013, 2014a,b; Lee et al., 2004; Niu et al., 2014) (Fig. 2). Whereas GPPS initially downregulated in the first group and later expressed more in the second group. Among these most of the genes were shown similar expression pattern from literature, in silico and qRT-PCR. Further, the downstream of ginsenoside pathways is dominated by cytochromes (CYP) and glycosyltransferases (GT) to produce protopanaxadiol and protopanaxatriol type ginsenosides. Scientists working for P. ginseng were more interested to characterize the downstream pathway, particularly on cytochromes (Devi et al., 2011; Han et al., 2011, 2012) and

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363

6

N

C

O

R

R

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

U

E

C

T

E

D P

a. Terpenoid backbone/ginsenoside responsible transcripts

b. Putative cytochrome responsible transcripts

R O

O

c. Putative glucosyl transferase responsible transcripts

Fig. 4. Hierarchical cluster and co-expression patterns for triterpenoid backbone pathway/putative ginsenoside pathway. a) Terpenoid backbone pathway. b) Cytochrome. c) Glycosyltransferase.

F

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

Acknowledgments

407 408

N C O

406

R

R

E

C

This research was supported by the MKE (Ministry of Knowledge Economy), Korea, under the ITRC (Information Technology Research 409 Center) support program supervised by the “NIPA (National IT Industry 410 Q20 Promotion Agency)” (NIPA-2011-(C1090-1121-0003)). Appendix A. Supplementary data

412 413

Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.gene.2014.05.024.

414

References

415 416 417 418 419 420 421 422 423 424

Akagi, T., Ikegami, A., Tsujimoto, T., Kobayashi, S., Sato, A., Kono, A., Yonemori, K., 2009. DkMyb4 is a Myb transcription factor involved in proanthocyanidin biosynthesis in persimmon fruit. Plant Physiology 151, 2028–2045. Ali, M., Hahn, E.-J., Paek, K.-Y., 2006a. Copper-induced changes in the growth, oxidative metabolism, and saponin production in suspension culture roots of Panax ginseng in bioreactors. Plant Cell Reports 25, 1122–1132. Ali, M., Yu, K.-W., Hahn, E.-J., Paek, K.-Y., 2006b. Methyl jasmonate and salicylic acid elicitation induces ginsenosides accumulation, enzymatic and non-enzymatic antioxidant in suspension culture Panax ginseng roots in bioreactors. Plant Cell Reports 25, 613–620.

U

411

Asaka, I., Ii, I., Hirotani, M., Asada, Y., Furuya, T., 1993. Production of ginsenoside saponins by culturing ginseng (Panax ginseng) embryogenic tissues in bioreactors. Biotechnology Letters 15, 1259–1264. Audic, S., Claverie, J.M., 1997. The significance of digital gene expression profiles. Genome Research 7, 986–995. Bae, K.-H., Choi, Y.-E., Shin, C.-G., Kim, Y.-Y., Kim, Y.-S., 2006. Enhanced ginsenoside productivity by combination of ethephon and methyl jasmonate in ginseng (Panax ginseng C.A. Meyer) adventitious root cultures. Biotechnology Letters 28, 1163–1166. Chen, S., Luo, H., Li, Y., Sun, Y., Wu, Q., Niu, Y., Song, J., Lv, A., Zhu, Y., Sun, C., Steinmetz, A., Qian, Z., 2011. 454 EST analysis detects genes putatively involved in ginsenoside biosynthesis in Panax ginseng. Plant Cell Reports 1–9. Choi, K.T., 2008. Botanical characteristics, pharmacological effects and medicinal components of Korean Panax ginseng C A Meyer. Acta Pharmacologica Sinica 29, 1109–1118. Choi, D.-W., Jung, J., Ha, Y.I., Park, H.-W., In, D.S., Chung, H.-J., Liu, J.R., 2005. Analysis of transcripts in methyl jasmonate-treated ginseng hairy roots to identify genes involved in the biosynthesis of ginsenosides and other secondary metabolites. Plant Cell Reports 23, 557–566. De Geyter, N., Gholami, A., Goormachtig, S., Goossens, A., 2012. Transcriptional machineries in jasmonate-elicited plant secondary metabolism. Trends in Plant Science 17, 349–359. de Hoon, M.J., Imoto, S., Nolan, J., Miyano, S., 2004. Open source clustering software. Bioinformatics 20, 1453–1454. Devi, B.S., Kim, Y.J., Sathiyamoorthy, S., Khorolragchaa, A., Gayathri, S., Parvin, S., Yang, D.U., Selvi, S.K., Lee, O.R., Lee, S., Yang, D.C., 2011. Classification and characterization of putative cytochrome P450 genes from Panax ginseng C. A. Meyer. Biochemistry (Mosc) 76, 1347–1359. Foster, T., Kirk, C., Jones, W., Allan, A., Espley, R., Karunairetnam, S., Rakonjac, J., 2007. Characterisation of the DELLA subfamily in apple (Malus × domestica Borkh.). Tree Genetics & Genomes 3, 187–197. Han, J.-Y., In, J.-G., Kwon, Y.-S., Choi, Y.-E., 2010. Regulation of ginsenoside and phytosterol biosynthesis by RNA interferences of squalene epoxidase gene in Panax ginseng. Phytochemistry 71, 36–46. Han, J.-Y., Kim, H.-J., Kwon, Y.-S., Choi, Y.-E., 2011. The cytochrome P450 enzyme CYP716A47 catalyzes the formation of protopanaxadiol from dammarenediol-II during ginsenoside biosynthesis in Panax ginseng. Plant and Cell Physiology. Han, J.Y., Hwang, H.S., Choi, S.W., Kim, H.J., Choi, Y.E., 2012. Cytochrome P450 CYP716A53v2 catalyzes the formation of protopanaxatriol from protopanaxadiol during ginsenoside biosynthesis in Panax ginseng. Plant & Cell Physiology 53, 1535–1545. Hu, F.-X., Zhong, J.-J., 2007. Role of jasmonic acid in alteration of ginsenoside heterogeneity in elicited cell cultures of Panax notoginseng. Journal of Bioscience and Bioengineering 104, 513–516. Hu, F.-X., Zhong, J.-J., 2008. Jasmonic acid mediates gene transcription of ginsenoside biosynthesis in cell cultures of Panax notoginseng treated with chemically synthesized 2-hydroxyethyl jasmonate. Process Biochemistry 43, 113–118. Hu, X., Neill, S.J., Cai, W., Tang, Z., 2003. Nitric oxide mediates elicitor-induced saponin synthesis in cell cultures of Panax ginseng. Functional Plant Biology 30, 901–907. Huang, X., Madan, A., 1999. CAP3: a DNA sequence assembly program. Genome Research 9, 868–877. Hwang, S.H., Shin, T.J., Choi, S.H., Cho, H.J., Lee, B.H., Pyo, M.K., Lee, J.H., Kang, J., Kim, H.J., Park, C.W., Shin, H.C., Nah, S.Y., 2012. Gintonin, newly identified compounds from ginseng, is novel lysophosphatidic acids–protein complexes and activates G protein-coupled lysophosphatidic acid receptors with high affinity. Molecules and Cells 33, 151–162. Jeong, G.-T., Park, D.-H., 2005. Enhancement of growth and secondary metabolite biosynthesis: effect of elicitors derived from plants and insects. Biotechnology and Bioprocess Engineering 10, 73–77. Jeong, G.-T., Park, D.-H., 2006. Enhanced secondary metabolite biosynthesis by elicitation in transformed plant root system. Applied Biochemistry and Biotechnology 130, 436–446. Jeong, G.T., Park, D.H., Ryu, H.W., Hwang, B., Woo, J.C., Kim, D., Kim, S.W., 2005. Production of antioxidant compounds by culture of Panax ginseng C.A. Meyer hairy roots: I. Enhanced production of secondary metabolite in hairy root cultures by elicitation. Applied Biochemistry and Biotechnology 121–124, 1147–1157. Jung, J.D., Park, H.W., Hahn, Y., Hur, C.G., In, D.S., Chung, H.J., Liu, J.R., Choi, D.W., 2003. Discovery of genes for ginsenoside biosynthesis by analysis of ginseng expressed sequence tags. Plant Cell Reports 22, 224–230. Khorolragchaa, A., Kim, Y.J., Rahimi, S., Sukweenadhi, J., Jang, M.G., Yang, D.C., 2014. Grouping and characterization of putative glycosyltransferase genes from Panax ginseng Meyer. Gene 536, 186–192. Kim, O., Kim, S., Ohyama, K., Muranaka, T., Choi, Y., Lee, H., Kim, M., Hwang, B., 2010. Upregulation of phytosterol and triterpene biosynthesis in Centella asiatica hairy roots overexpressed ginseng farnesyl diphosphate synthase. Plant Cell Reports 29, 403–411. Kim, Y.-K., Kim, J.K., Kim, Y.B., Lee, S., Kim, S.-U., Park, S.U., 2013. Enhanced accumulation of phytosterol and triterpene in hairy root cultures of Platycodon grandiflorum by overexpression of Panax ginseng 3-hydroxy-3-methylglutaryl-coenzyme A reductase. Journal of Agricultural and Food Chemistry 61, 1928–1934. Kim, Y.J., Jeon, J.N., Jang, M.G., Oh, J.Y., Kwon, W.S., Jung, S.K., Yang, D.C., 2014a. Ginsenoside profiles and related gene expression during foliation in Panax ginseng Meyer. Journal of Ginseng Research 38, 66–72. Kim, Y.J., Lee, O.R., Oh, J.Y., Jang, M.G., Yang, D.C., 2014b. Functional analysis of HMGR encoding genes in triterpene saponin-producing Panax ginseng Meyer. Plant Physiology.

O

R O

P

D

E

T

glycosyltransferase (Khorolragchaa et al., 2014). Totally, 533 (0.9%) transcripts were annotated with the term “cytochrome” in our study 366 and 33 (0.05%) transcripts were involved in ≥ 2 fold regulations 367 (Fig. 3). The CYP73A100 (PG019176), CYP71D312 (PG000067) and 368 CYP716A47 (PG000594) are found with similar patterns as in literature 369 reports, particularly CYP73A100 was continuously up-regulated (Han 370 et al., 2011, 2012). In our result most of CYP was annotated from 371 P. ginseng and P. notoginseng and two GTs were from P. notoginseng 372 which are found in ≥ 2 fold up and downregulations. Most annotated 373 sequences were from P. notoginseng, the study reported full length 374 CYP and GT from RNA sequencing (Khorolragchaa et al., 2014; Luo 375 et al., 2011). These literature evidences mimic our RNA sequencing 376 CYP and GT gene expression patterns. With these conformations, we 377 propose that all the listed putative genes i.e. 30 transcription factors, 378 20 cytochromes, 11 glycosyltransferases (Figs. 3 and 4) were the candi379 dates for ginsenoside pathway. In addition, the co-expression patterns 380 of CYP and GT transcripts along with transcription factors were obtained 381 (Supplementary Fig. 5). Result showed that, putative uncharacterized 382 CYPs (PG027814, PG024387, PG024073, PG019557, PG002087, 383 PG005498, PG000598, PG001995) and GT2 (PG010742, PG002718, 384 PG025219, PG002650, PG000627) were highly co-expressed with 385 other ginsenoside pathway related transcripts and transcription factors. 386 From these analyzes, our predicted transcripts are the probable candi387 dates for ginsenoside pathway directly or indirectly. Moreover, ≥ 2 388 fold up and downregulated transcripts were grouped according to 389 KEGG secondary metabolite pathways. The results exhibited that the 390 transcripts are highly correlated with phenylpropenoid pathway than 391 others (Supplementary Table 6). According to this result mostly MeJA 392 modulates the phenylpropenoid genes that are similar to other plants 393 Q17 (Pauwels et al., 2009; Zhao et al., 2005). In conclusion, molecular net394 work dynamics present in P. ginseng was not observed or predicted till 395 now by partial or complete data set using high throughput experiments. 396 This study aimed to uncover expression patterns of highly expressed 397 transcripts for P. ginseng by in silico digital gene expression analyses 398 through RNA-seq. A study based on Sanger sequencing, but there was 399 no comparison between samples (Jung et al., 2003) and another study 400 compared two sets of sequences (Sathiyamoorthy et al., 2010a), but 401 Q18 samples were prepared in different conditions. Therefore, this study 402 has given an expression data set from the same type of adventitious 403 root samples with time dependent conditions. Moreover this data set 404 will support for molecular network predictions from P. ginseng and 405 Q19 ginsenoside pathway.

F

364 365

7

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 Q22 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 Q23

D

P

R O

O

F

Sathiyamoorthy, S., In, J.G., Gayathri, S., Kim, Y.J., Yang, D.C., 2010b. Generation and gene ontology based analysis of expressed sequence tags (EST) from a Panax ginseng C. A. Meyer roots. Molecular Biology Reports 37, 3465–3472. Sathiyamoorthy, S., In, J.G., Lee, B.S., Kwon, W.S., Yang, D.U., Kim, J.H., Yang, D.C., 2011. In silico analysis for expressed sequence tags from embryogenic callus and flower buds of Panax ginseng C. A. Meyer. Journal of Ginseng Research 35, 21–30. Sathiyaraj, G., Lee, O., Parvin, S., Khorolragchaa, A., Kim, Y.-J., Yang, D., 2011. Transcript profiling of antioxidant genes during biotic and abiotic stresses in Panax ginseng C. A. Meyer. Molecular Biology Reports 38, 2761–2769. Sawai, S., Saito, K., 2011. Triterpenoid biosynthesis and engineering in plants. Frontiers in Plant Science 2. Shi, W., Wang, Y., Li, J., Zhang, H., Ding, L., 2007. Investigation of ginsenosides in different parts and ages of Panax ginseng. Food Chemistry 102, 664–668. Shim, J.S., Lee, O.R., Kim, Y.J., Lee, J.H., Kim, J.H., Jung, D.Y., In, J.G., Lee, B.-S., Yang, D.C., 2010. Overexpression of PgSQS1 increases ginsenoside production and negatively affects ginseng growth rate in Panax ginseng. Journal of Ginseng Research 34, 99–103. Sun, C., Li, Y., Wu, Q., Luo, H., Sun, Y., Song, J., Lui, E.M., Chen, S., 2010. De novo sequencing and analysis of the American ginseng root transcriptome using a GS FLX titanium platform to discover putative genes involved in ginsenoside biosynthesis. BMC Genomics 11, 262. Sun, Y., Niu, Y., Xu, J., Li, Y., Luo, H., Zhu, Y., Liu, M., Wu, Q., Song, J., Sun, C., Chen, S., 2013. Discovery of WRKY transcription factors through transcriptome analysis and characterization of a novel methyl jasmonate-inducible PqWRKY1 gene from Panax quinquefolius. Plant Cell, Tissue and Organ Culture (PCTOC) 114, 269–277. Suttipanta, N., Pattanaik, S., Kulshrestha, M., Patra, B., Singh, S.K., Yuan, L., 2011. The transcription factor CrWRKY1 positively regulates the terpenoid indole alkaloid biosynthesis in Catharanthus roseus. Plant Physiology 157, 2081–2093. Wu, J., Zhong, J.-J., 1999. Production of ginseng and its bioactive components in plant cell culture: current technological and applied aspects. Journal of Biotechnology 68, 89–99. Wu, B., Wang, M., Ma, Y., Yuan, L., Lu, S., 2012. High-throughput sequencing and characterization of the small RNA transcriptome reveal features of novel and conserved microRNAs in Panax ginseng. PLoS One 7, e44385. Wu, D., Austin, R.S., Zhou, S., Brown, D., 2013. The root transcriptome for North American ginseng assembled and profiled across seasonal development. BMC Genomics 14, 564. Yu, K.W., Murthy, H.N., Jeong, C.S., Hahn, E.J., Paek, K.Y., 2005. Organic germanium stimulates the growth of ginseng adventitious roots and ginsenoside production. Process Biochemistry 40, 2959–2961. Zhao, J., Davis, L.C., Verpoorte, R., 2005. Elicitor signal transduction leading to production of plant secondary metabolites. Biotechnology Advances 23, 283–333.

N

C

O

R

R

E

C

T

Lee, M.-H., Jeong, J.-H., Seo, J.-W., Shin, C.-G., Kim, Y.-S., In, J.-G., Yang, D.-C., Yi, J.-S., Choi, Y.-E., 2004. Enhanced triterpene and phytosterol biosynthesis in Panax ginseng overexpressing squalene synthase gene. Plant and Cell Physiology 45, 976–984. Li, C., Zhu, Y., Guo, X., Sun, C., Luo, H., Song, J., Li, Y., Wang, L., Qian, J., Chen, S., 2013. Transcriptome analysis reveals ginsenosides biosynthetic genes, microRNAs and simple sequence repeats in Panax ginseng C. A. Meyer. BMC Genomics 14, 245. Livak, K.J., Schmittgen, T.D., 2001. Analysis of relative gene expression data using realtime quantitative PCR and the 2(−delta delta C(T)) method. Methods 25, 402–408. Lu, Wong, Teng, 2001. Effects of elicitation on the production of saponin in cell culture of Panax ginseng. Plant Cell Reports 20, 674–677. Luo, H., Sun, C., Sun, Y., Wu, Q., Li, Y., Song, J., Niu, Y., Cheng, X., Xu, H., Li, C., Liu, J., Steinmetz, A., Chen, S., 2011. Analysis of the transcriptome of Panax notoginseng root uncovers putative triterpene saponin-biosynthetic genes and genetic markers. BMC Genomics 12 (Suppl. 5), S5. Mathiyalagan, R., Subramaniyam, S., Natarajan, S., Kim, Y.J., Sun, M.S., Kim, S.Y., Kim, Y.J., Yang, D.C., 2013. In silico profiling of microRNAs in Korean ginseng (Panax ginseng Meyer). Journal of Ginseng Research 37, 227–247. Mathiyalagan, R., Subramaniyam, S., Kim, Y.J., Natarajan, S., Min, J.W., Kim, S.Y., Yang, D.C., 2014. Synthesis and pharmacokinetic characterization of a pH-sensitive polyethylene glycol ginsenoside CK (PEG-CK) conjugate. Bioscience, Biotechnology, and Biochemistry 1–3. Niu, Y., Luo, H., Sun, C., Yang, T.J., Dong, L., Huang, L., Chen, S., 2014. Expression profiling of the triterpene saponin biosynthesis genes FPS, SS, SE, and DS in the medicinal plant Panax notoginseng. Gene 533, 295–303. Ok, R.L., Gayathri, S., Yu, j.k, In, J.G., Kwon, W.S., Kim, J.H., Yang, D.C., 2011. Defense genes induced by pathogens and abiotic stresses in Panax ginseng C.A. Meyer. Journal of Ginseng Research 35, 1–11. Pauwels, L., Inze, D., Goossens, A., 2009. Jasmonate-inducible gene: what does it mean? Trends in Plant Science 14, 87–91. Pollier, J., Moses, T., Gonzalez-Guzman, M., De Geyter, N., Lippens, S., Vanden Bossche, R., Marhavy, P., Kremer, A., Morreel, K., Guerin, C.J., Tava, A., Oleszek, W., Thevelein, J.M., Campos, N., Goormachtig, S., Goossens, A., 2013. The protein quality control system manages plant defence compound synthesis. Nature 504, 148–152. Qian, Z.G., Zhao, Z.J., Xu, Y.F., Qian, X.H., Zhong, J.J., 2004. Novel chemically synthesized hydroxyl-containing jasmonates as powerful inducing signals for plant secondary metabolism. Biotechnology and Bioengineering 86, 809–816. Quan, L.H., Min, J.W., Sathiyamoorthy, S., Yang, D.U., Kim, Y.J., Yang, D.C., 2012. Biotransformation of ginsenosides Re and Rg1 into ginsenosides Rg2 and Rh1 by recombinant beta-glucosidase. Biotechnology Letters 34, 913–917. Sathiyamoorthy, S., In, J.G., Gayathri, S., Kim, Y.J., Yang, D., 2010a. Gene ontology study of methyl jasmonate-treated and non-treated hairy roots of Panax ginseng to identify genes involved in secondary metabolic pathway. Genetika 46, 932–939.

U

511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552

S. Subramaniyam et al. / Gene xxx (2014) xxx–xxx

E

8

Please cite this article as: Subramaniyam, S., et al., Transcript expression profiling for adventitious roots of Panax ginseng Meyer, Gene (2014), http://dx.doi.org/10.1016/j.gene.2014.05.024

553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593