In silico identification of novel and conserved microRNAs and targets in peppermint (Mentha piperita) using expressed sequence tags (ESTs)

2.1 Computational reorganization of mint miRNA

In total 74 mature miRNA taken from miRNA sequence database miRBASE Version 21 (https://www.mirbase.org). To find conserved and novel miRNAs in peppermint, 43 distinctive miRNAs were predicted and chosen for next studies. On the whole 1216 EST sequences of peppermint (Mentha piperita) are taken in NCBI dbEST database (https://www.ncbi.nlm.nih.gov/nucest/). Blast 2.2.30 + an offline version of BLAST was used for further studies (ftp://ftp.ncbi.nlm.nih.gov/blast). Nearly 1188 contigs were identified using CD HIT server from which 1061 Unigenes were annotated.

2.2 Expressed sequence tag assembly and functional annotation

The ESTs were assembled employing CLC Genomics Workbench 8.0 (https://www.clcbio.com/products/ clc-genomics-workbench/). Gene ontologies and functional annotation were analyzed by BLAST2GO suite21. BLASTn and BLASTx servers were employed in similarity search against the viridiplantae family. The unigenes were annotated based on similarity search using of e-10.

2.3 Mature miRNA identification

Pair-wise alignment of assembled EST sequences against a unique miRNA database with a threshold value of e-10 and window size of 7 was carried out using the BLAST algorithm using mature miRNAs extracted from the miRBase database. miRNA candidates were selected using homology search approach of the threshold of 18nt length with no gaps. Further study was followed by the known miRNAs were selected that are closely matched with the assembled sequence and subject to BLASTx search against non-redundant protein database. Fig. 1 represents the miRNA scheme and EST analysis in peppermint.

Close

Fig. 1

Illustrative representation of the method followed for the identification of miRNAs from mint ESTs.

Export to PPT

2.4 miRNA secondary structure prediction

The non-protein coding predecessor sequence was employed to envisage the secondary structure using Mfold software which uses energy minimization method, protein coding sequences were removed prior to prediction (Kozomara and Griffiths-Jones, 2014; Zuker, 2003), filter criteria was followed to select pre-miRNA with the length of 60 nt properly folded to form stem-loop-hairpin secondary structure with maximum of 6 nt mismatches. There must be no breaks allowed between the miRNA/miRNA* duplex, A + U content must be between 30 and 70 percent, and the projected secondary structure must have a higher minimal folding free energy index (MFEI) and negative minimal folding free energy (MFE), which distinguish miRNAs from short RNAs. The negative free was calculated using the fact that: $$\mathrm{M}\mathrm{F}\mathrm{E}\mathrm{I}=\frac{\mathrm{M}\mathrm{F}\mathrm{E}\times 100}{\mathrm{L}\mathrm{e}\mathrm{n}\mathrm{g}\mathrm{t}\mathrm{h}\mathrm{o}\mathrm{f}\mathrm{p}\mathrm{r}\mathrm{e}\mathrm{m}\mathrm{i}\mathrm{R}\mathrm{N}\mathrm{A}\mathrm{s}\mathrm{e}\mathrm{q}\mathrm{u}\mathrm{e}\mathrm{n}\mathrm{c}\mathrm{e}}$$

2.5 miRNA target genes prediction

The important criteria for miRNA target prediction are seed match, free energy, conservation, target site abundance and site accessibility. The target miRNAs were predicted using homology algorithm (Rhoades et al., 2002), psRNA target server was used to find miRNA targets with default parameters. The prospective miRNAs as query are searched against with mRNAs and assembled ESTs of peppermint. Target accessibility was granted maximum energy to unpair the target site (UPE) by adhering to the specifications of maximum exception being 3.0, length of complementarily scoring (hspsize) being 20, and flanking length approximately about the size of target site for accessibility analysis which is 13 bp downstream and 17 bp upstream. The highest mismatch at the complementary site is 4 without gaps, while the choice of the central mismatch resulting to translational inhibition is 9–11 nucleotides.

2.6 Phylogenetic analysis of mature miRNAs

The standalone blast search tool was used against miRBase database and finally matured miRNAs were predicted. Multiple sequence alignment tool ClustalW, used to align sequence of precursor miRNA and newly identified miRNAs. MEGA version 6.25 was used to conduct phylogenetic and molecular evolutionary analyses and phylogenetic tree was built using maximum likelihood method (Olsen et al., 1994).

2.7 Functional annotation and pathway analysis of miRNA target genes

The target miRNA gene sequences were matched to the non-redundant protein sequence database of NCBI GO keywords using the BLAST2GO 3.0 tool (https://balst2go.com/blast2go-pro/download-b2g). The Gene Ontology (GO) database was used to look into biological processes, cellular components, and molecular functions. The Kyoto Encyclopedia of Genes and Genomes (KEGG) (Kanehisa and Goto, 2000) KASS server was used to annotate target miRNAs, its metabolic pathways, and their networks (Moriya et al., 2007).

2.8 Structure prediction of key metabolites and its binding pockets of M. Piperata

The 3D structure and active sites of identified key metabolites like (s) -limone 3- monooxygenase, (+) menthofuran synthase, (+) pulegone reductase, (+) neomenthol dehyrogenase and vestitone reductase were modeled and predicted by RaptorX (https://raptorx.uchicago.edu/BindingSite/) (Källberg et al., 2012).

3.1 Identification of miRNA in peppermint ESTs

A total of 1188 unigenes were used to identify miRNAs, with the majority of mature plant miRNAs being evolutionarily conserved. MiRBase (https://www.mirbase.org/cgi-bin/browse.pl) was used to obtain 74 mature miRNAs. After eliminating the duplicated sequences, 43 unique miRNAs were used as a reference, and 1188 unigenes were BLASTed against them using settings such as e-10, word, and window size of 7. With 9 unigenes, the search yielded 24 possible mature miRNA showing homologs. The anticipated miRNA length is between 16 and 22 nucleotides, which is consistent with earlier plant research. Mfold software was used to create secondary structures for selected miRNAs, and various filters were employed to pick candidate miRNAs. Nine miRNAs had higher free energy and content as prospective miRNAs, out of a total of 24. The stability of the hairpin stem loop secondary structure of the precursor miRNA was determined using the terms minimum free energy (MFE) and anticipated value (E-value). Therefore, a lower value of MFE implies that the secondary hairpin stem loop structure of the precursor miRNA is more stable. Fig. 3 provides the anticipated mRNA structures.

Close

Fig. 3

Structures of predicted mRNA.

Export to PPT

This study shows the occurrence of nine miRNA from 1216 ESTs of peppermint, according to Zhang, et al 2007 the likelihood of discovering miRNA is one per 1000 ESTs. The nucleotide percentage for recently recognized miRNA in peppermint ESTs are tabulated in Tables 3 and Table 4.

3.2 miRNA target prediction

Peppermint contigs/unigenes were analysed for the miRNAs prediction which regulates gene expression and actively involved in growth development, stress stimulus, defense mechanism. psRNA target server, a homology based tool was employed in prediction of peppermint miRNA from target mRNAs. A total of 19 targets for the new identified mature miRNAs of peppermint were reported. Among 19, nine were of miRNAs which are represented in the Table 3.

3.3 Phylogenetic analysis

The anticipated precursor miRNA sequences were aligned using Clustal W, and phylogenetic analysis was utilized to establish the evolutionary link using MEGA's maximum likelihood technique. MW298 and ML1312 were discovered to be closely related among 19 predicted miRNA sequences that are active in Monoterpene biosynthesis, indicating (+)-pulegone reductase and (4S)-limonene synthase, respectively, as illustrated in Fig. 4.Fig. 5..

Close

Fig. 4

Phylogenetic tree prepared using maximum likelihood method showing relationship between newly identified pre-miRNAs of mint with their closely related families.

Export to PPT

Close

Fig. 5

Functional annotations of miRNAs in M.piperata.

Export to PPT

3.4 Structure prediction of key metabolites and its binding pockets of M. Piperata

The binding pocket residues were found in the enzymes (s) -limone 3- monooxygenase (EC:1.14.14.143), (+) menthofuran synthase (EC:1.3.1.81), (+) pulegone reductase (EC:1.3.1.81), (+) neomenthol dehyrogenase (EC:1.3.1.81), and vestitone reductase (EC:1.3.1.81) (EC: 1.1.1.348). a. Binding pocket locations were utilised to assess the quality of an enzyme's predicted structure. It denotes the number of times the selected binding pocket was discovered in a set of important metabolites. From 19 to 23 residues, the number of binding pockets and their multiplicity were changed. Figs. 6 and 7 shows the binding residues of important metabolite structures as well as their positions. The hydrophobicity of the binding residues is roughly 50%.

Close

Fig. 6

Modelled structure of key metabolites and binding pockets of M. piperata. (A) Binding sites of (s) -limone 3- monooxygenase [G12 N14 K15 G16 I17 R37 L65 D66 V67 N93 A94 G95 V96 V180 S181 K242 G267 A269 T271 D272 273 T274]. (B) Binding sites of (+) menthofuran synthase [T68 A71 G72 T75 S76 T79 M130 P135 I136 I139 P140 R141 P201 F202 G203 R207 I208 C209 P210 G211 F214 G215 V219]. (C) binding sites of (+) pulegone reductase [G15 T17 G18 F19 L20 I39 R40 L69 M88 A89 H90 P91 L92 T130 S131 S132 Y168 K172 P195 T196 W197 I198 S210]. (D) binding sites of (+) neomenthol dehyrogenase [R99 V115 A116 R128 M292 A295 G296 T299 T300 A303 P359 F360 L363 L364 R366 P426 F427 G428 R432 G433 C434 P435 G436 F439 A440 Y444]. (E) Binding sites of vestitone reductase [P52 M135 T139 G163 S164 V165 A184 G185 K189 Y205 N229 V230 C251 G252 V254 S255 Y257 F281 V282 V283 V326 G329].

Export to PPT

Close

Fig. 7

Target motif metabolites gene in M. piperata.

Export to PPT