PgCYP76B93 docks on phenylurea herbicides and its expression enhances chlorotoluron tolerance in Arabidopsis

The phenylurea herbicides are used to control annual and perennial weeds on crop cultivating fields. The excessive usage of these agrochemicals increase many environmental problems. Thus, engineering transgenic plant for herbicide metabolism can provide efficient and eco-friendly means for enhanced phytoremediation capacity. Cytochrome P450 enzymes comprise one of the major plant enzyme families that mediate the oxidative degradation of xenobiotic chemicals, including herbicides. Considering these notions, phytoremediation properties of transgenic ginseng-derived PgCYP76B93 in Arabidopsis to phenylurea herbicides were assessed. Phylogenetic tree of PgCYP76B93 clustered in between close to the herbicide metabolism-related enzyme families and terpenoid biosynthesis-related. The expression of PgCYP76B93 was considerably upregulated upon treatment with phenylurea herbicide, chlorotoluron. Simulated docking using Autodoc program predicted possible interaction with chlorotoluron. Transgenic Arabidopsis plants overexpressing PgCYP76B93 were resulted in slightly reduced plant height with relatively small leaves. The lower plant height in the PgCYP76B93-overexpressing line than in the control revealed that it was linked to the expression of gibberellin oxidases (GAox). The bioassay of transgenic plants growing on herbicide-containing media revealed enhanced resistance against chlorotoluron.


Introduction
Plant cytochrome P450 enzymes (CYPs) are involved in catalyzing various types of monooxygenation/hydroxylation reactions in primary and secondary metabolism [1]. One third of all known CYP enzymes are from the plant kingdom, and the number of enzymes is estimated to be up to 1% of the total number of annotated plant species. This implies that diversification within CYPs has led to the emergence of new metabolic pathways throughout land plant evolution.
A newly diversified species-specific secondary metabolic pathway in Panax ginseng Meyer involves ginseng saponins, also called ginsenosides [2,3]. Considering the pharmaceutical efficacies of ginsenosides, studies on ginseng plants largely focused on ginsenoside biosynthesis [3,4]. Among the 116 reported P450 genes [5], two CYP genes have been reported to be involved in the dammarane-type of ginsenoside biosynthesis [6,7] and one CYP, called β-amyrin 28-oxidase, is involved in the oleananetype of ginsenoside biosynthesis [8]. However, besides their physiological roles in secondary metabolism, other important features of CYPs are their detoxification roles against herbicide [9]. The mechanism of acquiring increased herbicide metabolism has not yet been adequately demonstrated. A perennial plant, such as ginseng, has been grown commercially for up to 6 years or even more, and it can be exposed to herbicides over several years. Continuous exposure to herbicides might be involved in the diversification of more CYP family genes. Thus, taking the concept one step further, the functional characterization of other CYP family members from ginseng might shed light on herbicide metabolism.
Modern herbicides contributed a lot to sustainable global food production by minimizing time to remove weeds and substituting destructive soil cultivation [10]. Thus understanding herbicides resistance mechanism is a major challenge for modern agriculture [10]. In higher plants, the herbicide detoxification genes have accumulated through the evolution of multigene families, like glutathione S-transferase [11] or glycosyl transferase [12] and cytochrome P450 enzyme families [9]. The CYPs are by far the largest family of enzymes that mediate the oxidative degradation of xenobiotic chemicals, including herbicides. The CYPs are heme protein that gains electrons from NADPH to catalyze the formation of a molecule of water and an oxygenated product. Xenobiotics such as aromatics, pesticides, and hydrocarbons are usually synthesized for industrial and agricultural purposes. Since some of them are harmful to living organisms, organisms have evolved efficient system such as CYPs to eliminate xenobiotic absorption. Genetically engineered plants for herbicide metabolism provide efficient and eco-friendly means for enhancement of detoxification of harmful substances [13].
The first evidence of the involvement of plant P450 in the metabolism of the herbicide monuron was reported in the microsomal fraction of cotton seedlings by Frear in 1969 [14]. To date, several CYP family members have been characterized as being involved in herbicide metabolism. CYP76B1 was isolated from Jerusalem artichoke (Helianthus tuberosus) and was found to metabolize herbicides belonging to the phenylurea class [15][16][17]. CYP71A10 was isolated from soybean and determined to metabolize the phenylurea herbicide, chlortoluron, when expressed in yeast, although with an efficiency 10 times lower than that of CYP76B1 [18]. The CYP76C1, CYP76C2, and CYP76C4 genes also provided suitable evidence for herbicide resistance in Arabidopsis [17]. When human-derived CYP1A2 was overexpressed in Arabidopsis, it also increased tolerance to the phenylurea herbicide linuron [19]. Very recently ginseng-derived CYP736A12 was reported to be involved in chlortoluron and isoproturon tolerance when overexpressed in Arabidopsis [20].
In this study, we report another ginseng-derived CYP gene annotated as PgCYP76B93 that are involved in herbicide tolerance. Homology based molecular docking of PgCYP76B9 shows the more possible interaction with chlortoluron rather than isoproturon. Transgenic plants overexpressing PgCYP76B93 displayed relatively reduced plant height and herbicide chlortoluron tolerance with respect to control lines. The observed functions of PgCYP76B93 provide insight into its use as a selectable marker for genetic engineering of multifunctional crops having the ability to detoxify agrochemicals and environmental contaminants.

Plant materials and growth conditions
Panax ginseng Meyer 'Chun-Poong' was used as the plant material in this study. Arabidopsis thaliana (genotype Col-0) was used as a heterologous system. Mature Arabidopsis seeds were surface-sterilized for 5 min in 70% (v/v) ethanol and then rinsed three times with sterile water. The Arabidopsis seeds were sown on 1/2 MS medium (Duchefa Biochemie, Haarlem, Netherlands) containing 0.5 g/L MES (2-[N-morpholino] ethanesulfonic acid, 1% sucrose, and 0.8% phytoagar); the pH was adjusted to 5.7 using KOH. The seed-sown Petri plates were cold-treated for 2 days and then incubated under long-day photoperiod conditions-16 h light and 8 h dark-at 23 °C in a growth chamber.

Sequence analysis
The ginseng EST clones were analyzed by using the BLAST (Basic Local Alignment Search Tool) program against a specialized cytochrome P450 database (https :// drnel son.utmem .edu/Cytoc hrome P450.html). The amino acid sequences were analyzed using the online programs (https ://www.ncbi.nlm.nih.gov/Struc ture/cdd/wrpsb .cgi) and the ProtParam tool (https ://web.expas y.org/protp aram/). The amino acid sequences were aligned using the BioEdit program (version 7.1.9). A phylogenetic tree was generated according to the neighbor-joining method of the MEGA6 (version 6.06) program.

Chemical treatments in ginseng
To record gene expression levels with the treatment of herbicides, 10 μM of chlorotoluron (CAS Registry No. 15545-48-9, Cayman, USA) and isoproturon (CAS Registry No. 34123-59-6, Tokyo Chemical Industry, Japan) were treated to 4 weeks cultivated adventitious ginseng roots. The treated samples and their corresponding control were collected at different time intervals (0 h, 1 h, 4 h, 8 h, 12 h, 24 h, and 48 h) for qPCR.

RNA isolation and quantitative real-time PCR (qPCR)
Total RNA was isolated using the RNeasy Plant Mini Kit (Qiagen, Valencia, CA, USA) according to the manufacturer's instructions with modifications. The DNase I (Takara, Japan) treatment was additionally performed before the washing step in a 100 μL reaction volume containing 1 μL of RNase inhibitor for 1 h. The concentration of the extracted RNA was determined using the Nano-MD UV-Vis spectrophotometer (Scinco, Seoul, Korea). To synthesize the first-strand of cDNA, total RNA (up to 5 μg maximum) was reverse transcribed using the Rever-tAid reverse transcriptase (Thermo, USA). A qPCR was performed using the Thermal Cycle Dice real-time PCR system (Takara, Shiga, Japan) according to the manufacturer's instructions. The qPCR was performed in a 20 μL reaction volume using the following thermal cycling conditions: initial denaturation at 95 °C for 30 s, followed by 40 cycles of 95 °C for 5 s and 60 °C for 30 s, an additional cycle of 95 °C for 15 s and 60 °C for 30 s, and a final dissociation step at 95 °C for 15 s. At the end of the qPCR, a dissociation curve was generated to evaluate the generation of by-products. To determine the absolute and/or relative fold-differences in template abundance for each sample, the threshold cycle (Ct) value of each sample was normalized to that of β-actin, and the formula 2 −△△Ct (relative quantification) was used to calculate the gene expression levels. At least three independent experiments were performed in order to validate the ΔΔCt method. Dissociation analysis at the end of each run was performed to confirm the specificity of the reaction. The gene-specific primers for PgCYP76B93 were 5ʹ-CAC AAG TCC CTG TCC ACC-3′ (forward) and 5′-GTG GTT GTG GAC GTG GAG-3′ (reverse). The control primers for ginseng β-actin (DC03005B05) were 5′-AGA GAT TCC GCT GTC CAG AA-3′ (forward) and 5′-ATC AGC GAT ACC AGG GAA CA-3′ (reverse). Other primer sequences that were used are listed in Additional file 1: Table S1.

Molecular docking study of PgCYP76B93 with isoproturon and chlorotoluron
The structure of CYP76B93 was predicted by ExPASy-SWISS MODEL using Zebra Fish CYP-450 17A1 mutant Abiraterone complex (PDB code: 6B82) [21] as major template which exhibited 48% sequence similarity, because the crystal structure of CYP76B93 have not been defined yet. Molecules that are used to identify interaction with CYP76B93 are isoproturon (CID: 36679) and chlorotoluron (CID: 27375) retrieved from the PubChem database. These molecules were exported as the standard 3-dimentional XML format using the Open Babel 2.4.1. Docking of the two small molecules with CYP76B93 was processed by the AutoDock vina (The Scripps Research Institute, La Jolla, CA, USA) and set as 20 docking poses in the order of the highest negative ΔG. The flexibility of the ligands was kept as possible and placed at the interface region. The binding possibility of isoproturon and chlorotoluron with CYP76B93 were examined according to their binding energy, and then specific residues involved in interaction with CYP76B93 were determined by LigPlot + [22] (The European Bioinformatics Institute).

Vector map construction and in-planta transformation
The full-length cDNA (1551 bp) of PgCYP76B93 was amplified to functionally characterize in the Arabidopsis heterologous system. The in-planta gene transfer was carried out by using the constitutive 35S promoterdriven pCAMBIA1390 [23]. PgCYP76B93 was amplified using the following primers containing added restriction enzyme sites (SalI and MfeI): 5′-AA GTC GAC ATG GAT ATC TTA ACC ATG-3′ and 5′-CG CAA TTG AAT CAT AAC TGG AGT TG-3′. The purified PCR product was inserted into pCAMBIA1390 by performing overnight ligation at 16 °C. The PgCYP76C9-ligated construct was confirmed by nucleotide sequencing and transformed into Arabidopsis using Agrobacterium tumefaciens strain C58C1 (pMP90) [24]. The transgenic transformants were selected on hygromycin-containing plates (50 μg/mL), and more than 15 T1 independent lines were selected. The homozygous transgenic lines carrying one copy of the insert and following a Mendelian segregation ratio were further characterized. For further data analysis, the Col-0 and empty vector lines were used as controls for the PgCYP76B93-overexpressing lines.

Phenylurea herbicide resistance test
Surface-sterilized seeds were sown on the 1/2 MS media containing two phenylurea herbicide, chlorotoluron. After 2 days of cold stratification at 4 °C, the seedlings were grown at 23 °C under long-day photoperiod conditions of 16 h light/8 h dark. Chlorotoluron was added to the medium at different concentrations for preliminary tests, and a physiological concentration was chosen to show a clear-cut difference in tolerance.

Statistical analysis
Significance of the statistical analysis was measured based of the determination of Student's t-test using MS excel software.

Isolation and identification of ginseng PgCYP76B93 genes
All cytochrome P450 (CYP) genes homologous with known CYP family genes were identified from the expressed sequence tag (EST) libraries constructed from the embryogenic calli, 4-and 14-year-old roots, leaves, bud, and methyl jasmonic acid (MeJA)-treated adventitious roots [25]. One EST clone showing similarity with CYP76B family was further chosen, and the recovery of its full-length cDNA was confirmed by rapid amplification of cDNA ends (RACE) polymerase chain reaction (PCR) [26]. The full-length cDNA of CYP76B was annotated again as PgCYP76B93 (CYP76B93 from P. ginseng) after consulting with the CYP nomenclature committee [27]. The original EST clone coding for PgCYP76B93 was obtained from MeJA-treated adventitious roots, indicating that it can be induced by an elicitor and is expressed weakly in normal conditions. PgCYP76B93 is a 1551 bp long gene encoding 517 amino acids. Other amino acid sequences of the characterized CYPs in plants were collected from various gene banks and were used to generate a phylogenic tree (Fig. 1a). PgCYP76B93 was clustered in between the CYP76B family enzymes involved in terpenoid biosynthesis [28,29] and several CYP76B and CYP76C family enzymes [15,16,18], which were previously reported to be involved in phenylurea herbicide resistance (Fig. 1a). The enzymes belonging to the A-type subfamily of the CYP76 family enzymes were clustered relatively closer with PgCYP76B93, and three ginseng CYP enzymes involved in ginsenoside biosynthesis [6][7][8] were clustered far more distantly (Fig. 1a). This indicates that PgCYP76B93 likely plays a role in herbicide resistance and/or terpenoid biosynthesis. The PgCYP76B93 enzyme was predicted to have a molecular weight of 58,390 kDa and a pI of 7.74 using the ProtParam program [30]. The CYP family enzymes were reported to share low similarity in their amino acid sequences, except for several conserved domains, which are necessary for their tertiary structure and enzymatic function [31]. The most well-conserved motif is the heme-binding region FxxGxRxCxG (also known as the CxG motif ), in which C (cysteine) binds to the heme group. The ExxR and PER motifs form the E-R-R triad and play important roles in locking the structure of the heme pocket. The most conserved amino acids in the CYP family are glutamic acid (E) and arginine (R) in the ExxR motif and cysteine in the CxG motif, and these amino acids are also found in PgCYP76B93 (Fig. 1b). The least conserved motif, (A/G)Gx(D/E)T(T/S) (also known as AGxDTT), contributes to oxygen binding and activation. The AGxDTT motif is also well conserved in PgCYP76B93, but the PER motif is replaced by the PKR motif (Fig. 1b). Understanding the significance of the substitution of K (lysine) for E requires further investigation.

Transcript levels of PgCYP76B93 are gradually increased by chlorotoluron
In one step to understand the function of PgCYP76B93, modulation of transcripts against two phenylurea herbicides were evaluated according to the recent report [20]. Differential mRNA levels of PgCYP76B93 was quantified in time-dependent manner by treatment of 10 μM of chlorotoluron and isoproturon (Fig. 2) Relative gene expression   PgCYP76B93 were increased by chlorotoluron gradually from 1 to 24 h after treatments, and were initially upregulated by isoproturon upto 4 h (Fig. 2). It indicates that PgCYP76B93 might play roles in herbicide detoxification, especially chlorotoluron.

Molecular docking study shows putative phenylurea herbicide interaction with PgCYP76B93
The major secondary structure of PgCYP76B93 is pretty identical to Zebra Fish CYP-450 17A1 mutant Abiraterone complex [21]. The heme site of simulated PgCYP76B93 was not appeared but it could be predicted by superimposing the simulated PgCYP76B93 on the Zebra Fish CYP-450 17A1 Abiraterone complex (Additional file 2: Fig. S1). When isoproturon and chlorotoluron were added to the AutoDock program to test the binding possibility with PgCYP76B93, possible docking simulation for the interaction was observed. Total 5 of 20 docking poses for isoproturon and 6 of 20 docking poses for chlorotoluron are resided around heme in active site (Fig. 3a, b). These positions are interestingly located on the similar site of abiraterone which is a substrate of the template protein, Zebra Fish CYP-450 17A1. Residues interacting with isoproturon and chlorotoluron around heme of CYP76B93 was identified by Ligplot+ program [22] (Fig. 3). Isoproturon dominantly forms hydrophobic interaction with Arg102, Val119, Trp120, Met214, Asp315, Val318, Ala319, Thr323, Pro383, Ala384, Leu387, Ile388, Arg456, Pro459, Gly460 and Ile500 in 5 poses, and minor hydrogen bond with Ile388 and Arg390 in 1 poses near heme (Fig. 3c). Chlorotoluron majorly contacts with Val119, Trp120, Asp315, Val318, Ala319, Asp322, Thr323, Pro383, Ala384, Leu387, Ile388, Arg390, Arg456, Ile500 and Thr501 in 6 poses by hydrophobic interaction and minorly interacts with Ile388, Arg390 and Thr501 by hydrogen bond in 3 poses (Fig. 3d). Val119, Trp120, Val318, Ala319, Thr323, Thr384, Leu387, Ile388, Ile500 involved in interaction with isoproturon or chlorotoluron are frequently observed in each poses, so it can be further considered as a target for mutagenesis study. The closest distance from heme to isoproturon is 9.9 Å and to chlorotoluron is 4.3 Å (Fig. 3e). This result is believed to support the data that transgenic Arabidopsis possessing overexpression of CYP76B93 gene exhibits resistance against chlorotoluron other than isoproturon. Closely bound chlorotoluron to heme of PgCYP76B93 might be used as a substrate for undergoing detoxification mechanism.

Overexpression of PgCYP76B93 reduced plant height and leaf size
In order to understand the functional characteristics of PgCYP76B93 in plant growth and development, its fulllength coding sequence was overexpressed in heterologous Arabidopsis plants. Of the several T1 transgenic lines selected on the media containing an antibiotic, two were characterized further. These two transgenic lines showed on average 3.6 times more upregulation of PgCYP76B93 (Fig. 4a) and displayed reduced plant height (Fig. 4b).
In seedling stage, one strongly expressing line No.9 displayed reduced hypocotyl length (Additional (25 seedlings per plate, total 6 plates were analyzed). Scale bars = 1 cm. Data represent the mean ± standard error (SE) of three independent replicates at P < 0.05(*) and P < 0.01(**) by Student's t-test file 2: Fig. S2). It could be resulted in slightly reduced plant height in overexpression lines. Rosette leaves were smaller than those in the control lines (Fig. 4c). This implies that PgCYP76B93 might be involved in plant growth and development via phytohormone regulation which is involved in stem and leaf development. However, the overexpression of PgCYP76B93 did not alter the root growth phenotype (Fig. 4d, e) suggesting that the phytohormone-regulated PgCYP76B93 gene is more functional in aerial parts than in underground parts.

Transcripts of gibberellin oxidases are altered by the overexpression of PgCYP76B93
Gibberellins (GAs) are well-known plant hormones that promote stem elongation as well as seed germination [32]. The PgCYP76B93-mediated reduced plant height (Fig. 4b) led us to analyze whether the transcript levels of gibberellin biosynthesis-related genes were changed. Gibberellins are biosynthesized from geranylgeranyl diphosphate (GGDP), which is a common C20 precursor for diterpenoids, using three different classes of enzymes [33]: (1) terpene synthases (2) cytochrome P450 monooxygenases, and (3) 2-oxoglutarate-dependent dioxygenases. To uncover a possible link with bioactive GA biosynthesis, four GA monooxidase genes were analyzed. Two GA oxidase enzymes, GA20ox1 and GA3ox1, are involved in catalyzing the reaction of bioactive GA biosynthesis, whereas GA2ox1 and GA2ox2 are involved in the conversion of bioactive gibberellins into an inactive form. In two different developmental stages, all transcripts of GA20ox1 and GA3ox1 were significantly downregulated, whereas GA2ox1 and GA2ox2 were more upregulated (Fig. 5a, b). Though the mRNA levels of GA2ox2 were higher in the 28-day-old seedlings than in the 14-day-old seedlings, the general tendency was a decrease in the expression of the more bioactive form of GA oxidase gene, suggesting that the reduced height (Fig. 4b) of the PgCYP76B93 overexpression lines is caused by reduction in the levels of the more bioactive forms of gibberellins.

Constitutive overexpression of PgCYP76B93 confers herbicide resistance
PgCYP76B93 is closely clustered with the CYP76B family of enzymes, which have been reported to be involved in herbicide metabolism [15][16][17]. Also upregulation of PgCYP76B93 transcripts against chlorotoluron supported this notion. Thus, the PgCYP76B93 overexpression lines were evaluated for their ability to display resistance to phenylurea herbicide chlorotoluron. Several concentrations (0.5, 1, 3, 7, 10 μM) of chlorotoluron were analyzed from the germination, and the seedlings exposed to 7 and 10 μM of chlorotoluron showed distinct tolerant phenotypes compared to the control plants (Fig. 6a). Statistical analysis also suggests that the survival rate (%) on 1/2 MS media containing herbicide is significant (Fig. 6b). To evaluate the resistance of already germinated seedlings, 4-day-old seedlings from Col-0 and PgCYP76B93 overexpression lines (OE lines) were analyzed in the presence and absence of two different concentrations of chlorotoluron, and found OE lines show apparent tolerant growth phenotype (Additional file 2: Fig. S3).