Integrative analysis of the shikonin metabolic network identifies new gene connections and reveals evolutionary insight into shikonin biosynthesis

Abstract

Plant specialized 1,4-naphthoquinones present a remarkable case of convergent evolution. Species across multiple discrete orders of vascular plants produce diverse 1,4-naphthoquinones via one of several pathways using different metabolic precursors. Evolution of these pathways was preceded by events of metabolic innovation and many appear to share connections with biosynthesis of photosynthetic or respiratory quinones. Here, we sought to shed light on the metabolic connections linking shikonin biosynthesis with its precursor pathways and on the origins of shikonin metabolic genes. Downregulation of Lithospermum erythrorhizon geranyl diphosphate synthase (LeGPPS), recently shown to have been recruited from a cytoplasmic farnesyl diphosphate synthase (FPPS), resulted in reduced shikonin production and a decrease in expression of mevalonic acid and phenylpropanoid pathway genes. Next, we used LeGPPS and other known shikonin pathway genes to build a coexpression network model for identifying new gene connections to shikonin metabolism. Integrative in silico analyses of network genes revealed candidates for biochemical steps in the shikonin pathway arising from Boraginales-specific gene family expansion. Multiple genes in the shikonin coexpression network were also discovered to have originated from duplication of ubiquinone pathway genes. Taken together, our study provides evidence for transcriptional crosstalk between shikonin biosynthesis and its precursor pathways, identifies several shikonin pathway gene candidates and their evolutionary histories, and establishes additional evolutionary links between shikonin and ubiquinone metabolism. Moreover, we demonstrate that global coexpression analysis using limited transcriptomic data obtained from targeted experiments is effective for identifying gene connections within a defined metabolic network.

References

[1] Papageorgiou VP, Assimopoulou AN, Couladouros EA et al. The chemistry and biology of alkannin, shikonin, and related naphthazarin natural products. Chemistry (Weinheim an der Bergstrasse, Germany). 1999; 38: 270-301.

[2] Skoneczny D, Weston P, Zhu X et al. Metabolic profiling and identification of Shikonins in root periderm of two invasive Echium spp. weeds in Australia. Molecules. 2017; 22: 330.

[3] Zhu X, Skoneczny D, Weidenhamer JD et al. Identification and localization of bioactive naphthoquinones in the roots and rhizosphere of Paterson’s curse (Echium plantagineum), a noxious invader. J Exp Bot. 2016; 67: 3777-88.

[4] Widhalm JR, Rhodes D . Biosynthesis and molecular actions of specialized 1,4-naphthoquinone natural products produced by horticultural plants. Hortic Res. 2016; 3: 16046.

[5] Wang F, Yao X, Zhang Y et al. Synthesis, biological function and evaluation of Shikonin in cancer therapy. Fitoterapia. 2019; 134: 329-39.

[6] Schmid HV, Zenk MH . P-hydroxybenzoic acid and mevalonic acid as precursors of the plant naphthoquinone alkannin. Tetrahedron Lett. 1971; 12: 4151-5.

[7] Loscher R, Heide L . Biosynthesis of p-Hydroxybenzoate from p-Coumarate and p-Coumaroyl-coenzyme a in cell-free extracts of Lithospermum erythrorhizon cell cultures. Plant Physiol. 1994; 106: 271-9.

[8] Block A, Widhalm JR, Fatihi A et al. The origin and biosynthesis of the Benzenoid moiety of ubiquinone (coenzyme Q) in Arabidopsis. Plant Cell. 2014; 26: 1938-48.

[9] Soubeyrand E, Johnson TS, Latimer S et al. The Peroxidative cleavage of Kaempferol contributes to the biosynthesis of the Benzenoid moiety of ubiquinone in plants. Plant Cell. 2018; 30: 2910-21.

[10] Yazaki K, Kataoka M, Honda G et al. cDNA cloning and gene expression of phenylalanine ammonia-Lyase in Lithospermum erythrorhizon. Biosci Biotechnol Biochem. 1997; 61: 1995-2003.

[11] Yamamura Y, Ogihara Y, Mizukami H . Cinnamic acid 4-hydroxylase from Lithospermum erythrorhizon: cDNA cloning and gene expression. Plant Cell Rep. 2001; 20: 655-62.

[12] Singh RS, Gara RK, Bhardwaj PK et al. Expression of 3-hydroxy-3-methylglutaryl-CoA reductase, p-hydroxybenzoate-m-geranyltransferase and genes of phenylpropanoid pathway exhibits positive correlation with shikonins content in arnebia [Arnebia euchroma (Royle) Johnston]. BMC Mol Biol. 2010; 11: 1471-2199.

[13] Yazaki K, Kunihisa M, Fujisaki T et al. Geranyl Diphosphate:4-Hydroxybenzoate Geranyltransferase from Lithospermum erythrorhizon: cloning and characterization of a key enzyme in shikonin biosynthesis. J Biol Chem. 2002; 277: 6240-6.

[14] Yazaki K, Fukui H, Tabata M . Isolation of the intermediates and related metabolites of shikonin biosynthesis from Lithospermum erythrorhizon cell cultures. Chem Pharm Bull (Tokyo). 1986; 34: 2290-3.

[15] Vranová E, Coman D, Gruissem W . Network analysis of the MVA and MEP pathways for isoprenoid synthesis. Annu Rev Plant Biol. 2013; 64: 665-700.

[16] Henry LK, Gutensohn M, Thomas ST et al. Orthologs of the archaeal isopentenyl phosphate kinase regulate terpenoid production in plants. Proc Natl Acad Sci. 2015; 112: 10050-5.

[17] Henry LK, Thomas ST, Widhalm JR et al. Contribution of isopentenyl phosphate to plant terpenoid metabolism. Nat Plants. 2018; 4: 721-9.

[18] Gaisser S, Heide L . Inhibition and regulation of shikonin biosynthesis in suspension cultures of Lithospermum. Phytochemistry. 1996; 41: 1065-72.

[19] Ueoka H, Sasaki K, Miyawaki T et al. A cytosol-localized Geranyl Diphosphate synthase from Lithospermum erythrorhizon and its molecular evolution. Plant Physiol. 2020; 182: 1933-45.

[20] Auber RP, Suttiyut T, McCoy RM et al. Hybrid de novo genome assembly of red gromwell (Lithospermum erythrorhizon) reveals evolutionary insight into shikonin biosynthesis. Hortic Res. 2020; 7: 82.

[21] Tang CY, Li S, Wang YT et al. Comparative genome/transcriptome analysis probes Boraginales’ phylogenetic position, WGDs in Boraginales, and key enzyme genes in the alkannin/shikonin core pathway. Mol Ecol Resour. 2020; 20: 228-41.

[22] Meyer GW, Bahamon Naranjo MA, Widhalm JR . Convergent evolution of plant specialized 1,4-naphthoquinones: metabolism, trafficking, and resistance to their allelopathic effects. J Exp Bot. 2021; 72: 167-76.

[23] McCoy RM, Utturkar SM, Crook JW et al. The origin and biosynthesis of the naphthalenoid moiety of juglone in black walnut. Hortic Res. 2018; 5: 67.

[24] Tholl D . Biosynthesis and Biological Functions of Terpenoids in Plants. Adv Biochem Eng Biotechnol. 2015; 148: 63-106.

[25] Usadel B, Obayashi T, Mutwil M et al. Co-expression tools for plant biology: opportunities for hypothesis generation and caveats. Plant Cell Environ. 2009; 32: 1633-51.

[26] Wisecaver JH, Borowsky AT, Tzin V et al. A global Coexpression network approach for connecting genes to specialized metabolic pathways in plants. Plant Cell. 2017; 29: 944-59.

[27] Song W, Zhuang Y, Liu T . CYP82AR subfamily proteins catalyze C-1′ Hydroxylations of Deoxyshikonin in the biosynthesis of Shikonin and Alkannin. Org Lett. 2021; 23: 2455-9.

[28] Yamamura Y, Sahin FP, Nagatsu A et al. Molecular cloning and characterization of a cDNA encoding a novel apoplastic protein preferentially expressed in a shikonin-producing callus strain of Lithospermum erythrorhizon. Plant Cell Physiol. 2003; 44: 437-46.

[29] Zhao H, Baloch SK, Kong LR et al. Molecular cloning, characterization, and expression analysis of LeMYB1 from Lithospermum erythrorhizon. Biol Plant. 2014; 58: 436-44.

[30] Obayashi T, Kinoshita K . Rank of correlation coefficient as a comparable measure for biological significance of gene coexpression. DNA Res. 2009; 16: 249-60.

[31] Liesecke F, Daudu D, De Bernonville RD et al. Ranking genome-wide correlation measurements improves microarray and RNA-seq based global and targeted co-expression networks. Sci Rep. 2018; 8: Article number 10885. https://doi.org/10.1038/s41598-018-29077-3.

[32] Wu H, Gao L, Dong J, Yang X . Detecting overlapping protein complexes by rough-fuzzy clustering in protein-protein interaction networks. PLoS One. 2014; 9: 471-2.

[33] van Dongen S, Abreu-Goodger C . Using MCL to Extract Clusters from Networks. In: Bacterial Molecular Networks.Springer: New York, NY, 2012, 281-95.

[34] Song W, Zhuang Y, Liu T . Potential role of two cytochrome P450s obtained from Lithospermum erythrorhizon in catalyzing the oxidation of geranylhydroquinone during Shikonin biosynthesis. Phytochemistry. 2020; 175: 112375.

[35] Wang S, Wang R, Liu T et al. CYP76B74 catalyzes the 3′′-hydroxylation of geranylhydroquinone in shikonin biosynthesis. Plant Physiol. 2019; 179: 402-14.

[36] Zhu Y, Chu S-J, Luo Y-L et al. Involvement of LeMRP, an ATP-binding cassette transporter, in shikonin transport and biosynthesis in Lithospermum erythrorhizon. Plant Biol. 2017; 17: 1-9.

[37] Rai A, Nakaya T, Shimizu Y et al. De novo Transcriptome assembly and characterization of Lithospermum officinale to discover putative genes involved in specialized metabolites. Planta Med. 2018; 84: 920-34.

[38] Wu FY, Tang CY, Guo YM et al. Transcriptome analysis explores genes related to shikonin biosynthesis in Lithospermeae plants and provides insights into Boraginales’ evolutionary history. Sci Rep. 2017; 7: Article number 4477. https://doi.org/10.1038/s41598-017-04750-1.

[39] Fang R, Wu F, Zou A et al. Transgenic analysis reveals LeACS-1 as a positive regulator of ethylene-induced shikonin biosynthesis in Lithospermum erythrorhizon hairy roots. Plant Mol Biol. 2016; 90: 345-58.

[40] Yazaki K, Matsuoka H, Shimomura K et al. A novel dark-inducible protein, LeDI-2, and its involvement in root-specific secondary metabolism in Lithospermum erythrorhizon. Plant Physiol. 2001; 125: 1831-41.

[41] Tatsumi K, Yano M, Kaminade K et al. Characterization of Shikonin derivative secretion in Lithospermum erythrorhizon hairy roots as a model of lipid-soluble metabolite secretion from plants. Front Plant Sci. 2016; 7: 1-11.

[42] Ma S, Bachan S, Porto M et al. Discovery of stress responsive DNA regulatory motifs in arabidopsis. PLoS One. 2012; 7: e43198. https://doi.org/10.1371/journal.pone.0043198.

[43] Rushton PJ, Torres JT, Parniske M et al. Interaction of elicitor-induced DNA-binding proteins with elicitor response elements in the promoters of parsley PR1 genes. EMBO J. 1996; 15: 5690-700.

[44] Van Bel M, Diels T, Vancaester E et al. PLAZA 4.0: an integrative resource for functional, evolutionary and comparative plant genomics. Nucleic Acids Res. 2018; 46: D1190-6.

[45] Leebens-Mack JH, Barker MS, Carpenter EJ et al. One thousand plant transcriptomes and the phylogenomics of green plants. Nature. 2019; 574: 679-85.

[46] Moore BM, Wang P, Fan P et al. Robust predictions of specialized metabolism genes through machine learning. Proc Natl Acad Sci U S A. 2019; 116: 2344-53.

[47] Maeda H, Dudareva N . The shikimate pathway and aromatic amino acid biosynthesis in plants. Annu Rev Plant Biol. 2012; 63: 73-105.

[48] Papageorgiou V, Assimopoulou A, Samanidou V et al. Recent advances in chemistry, biology and biotechnology of Alkannins and Shikonins. Curr Org Chem. 2006; 10: 2123-42.

[49] Oshikiri H, Watanabe B, Yamamoto H et al. Two BAHD acyltransferases catalyze the last step in the shikonin/alkannin biosynthetic pathway. Plant Physiol. 2020; 184: 753-61.

[50] Stefely JA, Pagliarini DJ . Biochemistry of mitochondrial coenzyme Q biosynthesis. Trends Biochem Sci. 2017; 42: 824-43.

[51] Oubrie A, Rozeboom HJ, Kalk KH et al. Structure and mechanism of soluble quinoprotein glucose dehydrogenase. EMBO J. 1999; 18: 5187-94.

[52] Fatihi A, Latimer S, Schmollinger S et al. A dedicated type II NADPH dehydrogenase performs the penultimate step in the biosynthesis of vitamin K1 in Synechocystis and Arabidopsis. Plant Cell. 2015; 27: 1730-41.

[53] Dresler S, Szymczak G, Wójcik M . Comparison of some secondary metabolite content in the seventeen species of the Boraginaceae family. Pharm Biol. 2017; 55: 691-5.

[54] Langfelder P, Horvath S . Eigengene networks for studying the relationships between co-expression modules. BMC Syst Biol. 2007; 1: 54.

[55] Takanashi K, Nakagawa Y, Aburaya S et al. Comparative proteomic analysis of Lithospermum erythrorhizon reveals regulation of a variety of metabolic enzymes leading to comprehensive understanding of the Shikonin biosynthetic pathway. Plant Cell Physiol. 2019; 60: 19-28.

[56] Guo J, Zhou YJ, Hillwig ML et al. CYP76AH1 catalyzes turnover of miltiradiene in tanshinones biosynthesis and enables heterologous production of ferruginol in yeasts. Proc Natl Acad Sci U S A. 2013; 110: 12108-13.

[57] Miettinen K, Dong L, Navrot N et al. The seco-iridoid pathway from Catharanthus roseus. Nat Commun. 2014; 5: 3606. https://doi.org/10.1038/ncomms4606.

[58] Bakshi M, Oelmüller R . Wrky transcription factors jack of many trades in plants. Plant Signal Behav. 2014; 9: 1-18.

[59] Karimi M, Inzé D, Depicker A . GATEWAY™ vectors for agrobacterium-mediated plant transformation. Trends Plant Sci. 2002; 7: 193-5.

[60] Cui W, Liu W, Wu G . A simple method for the transformation of agrobacterium tumefaciens by foreign DNA. Chin J Biotechnol. 1995; 11: 267-74.

[61] Livak KJ, Schmittgen TD . Analysis of relative gene expression data using real-time quantitative PCR and the 2-ΔΔCT method. Methods. 2001; 25: 402-8.

[62] Ye J, Coulouris G, Zaretskaya I et al. Primer-BLAST: a tool to design target-specific primers for polymerase chain reaction. BMC Bioinformatics. 2012; 13: 134.

[63] Zhao H, Chang QS, Zhang DX et al. Overexpression of LeMYB1 enhances shikonin formation by up-regulating key shikonin biosynthesis-related genes in Lithospermum erythrorhizon. Biol Plant. 2015; 59: 429-35.

[64] McAdam SAM, Brodribb TJ . Mesophyll cells are the main site of abscisic acid biosynthesis in water-stressed leaves. Plant Physiol. 2018; 177: 911-7.

[65] Kim D, Paggi JM, Park C et al. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019; 37: 907-15.

[66] Robinson MD, Oshlack A . A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol. 2010; 11: R25. https://doi.org/10.1186/gb-2010-11-3-r25.

[67] Love MI, Huber W, Anders S . Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology. 2014; 15: 1-21.

[68] Benjamini Y, Hochberg Y . Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Stat Soc Ser B. 1995; 57: 289-300.

[69] Yu G, Wang L-G, Han Y et al. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012; 16: 284-7.

[70] Bray NL, Pimentel H, Melsted P et al. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016; 34: 525-7.

[71] Robinson MD, McCarthy DJ, Smyth GK . edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010; 26: 139-40.

[72] Langfelder P, Horvath S . WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008; 9: 559. https://doi.org/10.1186/1471-2105-9-559.

[73] Zhu Y, Lu GH, Bian ZW et al. Involvement of LeMDR, an ATP-binding cassette protein gene, in shikonin transport and biosynthesis in Lithospermum erythrorhizon. BMC Plant Biol. 2017; 17: 1-10.

[74] Emms DM, Kelly S . OrthoFinder: solving fundamental biases in whole genome comparisons dramatically improves orthogroup inference accuracy. Genome Biol. 2015; 16: 1-14.

[75] Katoh K, Standley DM . MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013; 30: 772-80.

[76] Nguyen LT, Schmidt HA, Von Haeseler A et al. IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Mol Biol Evol. 2015; 32: 268-74.

[77] Kalyaanamoorthy S, Minh BQ, Wong TKF et al. ModelFinder: fast model selection for accurate phylogenetic estimates. Nat Methods. 2017; 14: 587-9.

[78] Price MN, Dehal PS, Arkin AP . FastTree 2 - approximately maximum-likelihood trees for large alignments. PLoS One. 2010; 5: e9490.

[79] Marçais G, Delcher AL, Phillippy AM et al. MUMmer4: a fast and versatile genome alignment system. PLoS Comput Biol. 2018; 14: e1005944. https://doi.org/10.1371/journal.pcbi.1005944.