BACKGROUND
a microrna  is a small  non-coding rna molecule  that modulates the stability of mrna targets and their rate of translation into proteins  <cit> .
mirnas are present in the genome of vertebrates, plants, algae and even viruses and are involved in diverse and complex biological processes, like development and cell differentiation  <cit> , tumorigenesis  <cit>  and immunity  <cit> .
they can also alter plant gene expression in response to environmental stresses  <cit> .
in animals, maturation of canonical mirnas occurs in two steps: first, the long primary mirna transcript is processed within the nucleus into a ∼60– <dig> nucleotides  stem-loop hairpin precursor  by the enzyme drosha  <cit> .
afterwards, within the cytoplasm, the enzyme dicer cleaves the pre-mirna into a double stranded rna duplex  and into a loop.
the loop is degraded as a by-product  <cit> , whereas the rna duplex is unwound by helicase activity, releasing the mature mirna and the star sequence  <cit> .
the last is typically degraded whereas the mature mirna guides the microribonucleo-protein complex  to target messengers rnas  by partial sequence complementarity  <cit> .
machine learning in mirna recognition
differently from protein-coding genes, ncrna genes do not contain easily detectable signals in the sequence level  <cit> .
therefore, computational pipelines for mirna discovery rely on characteristic features of pre-mirnas.
coupled with comparative genomic computational pipelines, using rnaseq read libraries or completely ab-initio methods, machine learning  algorithms have played an important role in mirna discovery  <cit> .
ml algorithms induce models which are able to predict novel mirnas based on patterns learned from known pre-mirna sequences and from other rna hairpin-like sequences, such as pseudo pre-mirnas, transfer rna  and mrnas.
promir  <cit> , a probabilistic method, searches for pre-mirna in genomic sequences using sequence and structure features.
naïve bayes based probabilistic models were proposed to score mirnas  <cit>  and pre-mirna  <cit>  candidates.
support vector machines   <cit> , random forest   <cit> , relaxed variable kernel density estimator   <cit>  and generalized gaussian components based density estimation   <cit>  were used to induce classifiers for pre-mirna prediction.
a limitation of working directly with the sequence is that the information available may not be sufficient to infer accurate models for novel mirna prediction.
for example, the approximate number of rnas folding into hairpin-like secondary structures in the human genome, without filtering based on phylogenetic conservation, was estimated as  <dig> million  <cit> .
filters derived from other pre-mirna features reduced the number of pre-mirna candidates to around  <dig> .
therefore, the use of features that consider different aspects may reduce false positive rates in mirna detection.
feature sets investigated in the literature
some of the features commonly extracted from rna sequences for pre-mirna recognition may not help to distinguish between positive  and negative  classes.
therefore, the feature sets considered may have an important effect in the learning process.
several feature sets have been proposed for pre-mirna recognition  <cit> .
to comparatively evaluate the effect of these different sets, we investigated seven feature sets proposed in the literature, named here fsi, i∈{ <dig> ..,7}.
they were used to induce six classifiers and they contain most of the features employed in computational pipelines for pre-mirna discovery.
next, we briefly define each of these sets.
the first set, named fs <dig>  <cit> , has  <dig> sequence and structural features.
the second feature set fs <dig>  which corresponds to a subset of  <dig> features of fs <dig>  was used by  <cit>  to induce a classifier, named micropred.
the third feature set fs <dig>  was used to induce a classifier called g 2de  <cit>  and is composed by seven features also present in fs <dig>  fs <dig> is a set of  <dig> sequence-structural features used by the triplet-svm  <cit> .
fs <dig> is a set of  <dig>  sequence-structure motifs used by another classifier, mirident  <cit> .
fs <dig> is the feature set used by the mipred  <cit> .
this feature set merged fs <dig>  the minimum free energy of folding  and a stability measure .
finally, fs <dig> merged fs <dig> with four features of fs <dig> and three other features: percentage of low complexity regions detected in the sequence, maximal length of the amino acid string without stop codons and cumulative size of internal loops.
it was used to induce the huntmi  <cit> .
the classification performances reported by the mentioned tools are among the highest in previous works.
their approximate specificities and sensitivities are: triplet-svm , mipred , micropred , g 2de , mirident  and huntmi .
proposal and key findings
in this study, we investigated the discriminative power of seven rna feature sets, previously used in six tools developed for pre-mirna prediction.
among them are two sets composed of sequence-structure features  and five sets are a miscellany of rna features .
the investigation of a specific feature set, using a particular training data and learning algorithm, may insert learning biases.
as a consequence, the predictive performance for other training sets and learning algorithms could be different.
to deal with this problem, we evaluated each feature set using the learning algorithms, svms, rf and g 2de, which were used in the publications proposing those feature sets.
according to the experimental results, the miscellaneous feature sets produced more accurate predictive models than features sets composed from only sequence-structure patterns.
however, the differences in accuracy among miscellaneous feature sets are small or insignificant, despite their large differences in composition and dimensionality.
inspired by these results, we selected a subset of  <dig> features, of lower computational cost, but with a similar classification performance, when compared with fs1-fs <dig> and fs6-fs <dig> feature sets.
the classes of positive and negative test sets used in the experiments presented in this paper were predicted by the tools that we used as reference.
except for one tool, higher sensitivity was tied to lower specificity and vice versa.
methods
our goal was to investigate the predictive performance of rna features in distinguishing pre-mirnas from pseudo hairpins.
as such, we adopted seven feature sets and three learning algorithms.
the feature sets were used to induce classifiers for pre-mirnas in triplet-svm , mipred , micropred , g 2de , mirident  and huntmi .
svms were used in triplet-svm, micropred and mirident, whereas rf was used in mipred and huntmi.
generalized gaussian density estimator   <cit>  is not a tool for pre-mirna prediction in the sense that the features have to be computed by the user in his/her own pipeline.
nevertheless, we included g 2de because of its predictive performance and class distribution interpretability.
the subsections below provide details of our experiments.
data sets
human pre-mirna sequences were downloaded from mirbase  <dig>  <cit>  as the primary source of positive examples.
in an attempt to avoid overfitting, we removed redundant sequences.
as such, we clustered the available  <dig>  human pre-mirnas sequences using dnaclust  <cit>  such that sequences within a cluster shared 80% similarity.
then, one sequence of each cluster was randomly picked.
this yielded the set of positives composed of  <dig>  non-redundant pre-mirnas sequences.
the negative examples were the  <dig>  pseudo hairpins from human refseq genes, originally obtained by  <cit>  and subsequently used by  <cit> .
these sequences were obtained in order to keep basic features such as length distribution and minimal free energy of folding , similar to those observed in human pre-mirnas.
moreover, this set has no redundant sequences.
but, in order to adopt a uniform criterion for redundancy removal, we applied the same procedure adopted in the positive set.
only singleton clusters were formed.
in practice, it is expected that high similarity between positive and negative examples leads to higher specificity  <cit> .
experiments
the predictive performance of any model is dependent on the training set representativeness, which usually increases with the increase in the training set size.
typically, density based algorithms, such as g 2de, are more sensitive to the course of dimensionality and larger training sets are more likely to provide higher predictive performances.
we determined experimentally the training set size which would be suitable for any algorithm and feature set.
each experiment was repeated  <dig> times, in order to provide standard deviations of each classification performance estimation.
one repetition consisted of a test set, named here as gen and  <dig> training sets.
for a given repetition, gen was composed by  <dig> sequences of each class, which corresponded to 1/ <dig> of the  <dig>  non-redundant pre-mirna sequences.
the remaining positive and negative sequences were used to sample increments of  <dig> sequences of each class to compose the training sets of  <dig> ,…, <dig>  instances.
each feature set was computed from the same training and test sets.
in total, we worked with  <dig> test sets and 13× <dig> training sets.
the classification performances of the three algorithms converged to a threshold for training set sizes equal to  <dig> , for all feature sets.
therefore, we presented the results for the largest training set, which contained  <dig>  sequences.
classification performance measures
classification performance was measured as accuracy , sensitivity , specificity , f-measure  and mathew correlation coefficient ; see below.
this measures can be computed as given below, such that tp, fn, tn and fp are the numbers of true positives, false negatives, true negatives and false positives, respectively.
acc=100×tp+tntp+fn+tn+fpse=100×tptp+fnsp=100×tptn+fpfm=100×2×tp2×tp+fn+fpmcc=100×tp×tn−fp×fn×××)
the first three measures are commonly used whereas fm is prefered when a compromise between sensitivity and precision is desirable.
mcc measures the correlation between real and predicted classes and it is considered less biased towards class imbalance.
we presented the predictive performances by the mean and the standard deviation , over the  <dig> repetitions.
features
features used in this work are presented in table  <dig>  with references for the detailed descriptions.
the prediction of the secondary structure in this work considered the energy model, as implemented in rnafold  <cit>  and unafold  <cit> .
they predict the structure which gives the minimum free energy of folding .
we kept the same parameters used in the original publications.
detailed descriptions can be found in the corresponding references.
the sequence-structure features combined sequence nucleotide information and its predicted state at the secondary structure.
in fs <dig>  each feature represents the relative frequency of three contiguous nucleotides states at the secondary structure, fixing the middle character .
because the triplet-svm script excludes sequences with multiple loops, we implemented a python script to compute fs <dig> in any sequence.
the motifs in fs <dig> give the counts of its occurrence in the sequence-structure string.
this string is obtained by padding the nucleotide sequence with its respective predicted state at the secondary structure .
this set was obtained using the python script provided in the authors’ website.
to compute fs <dig>  we implemented a python script, based on the micropred perl pipeline.
fs <dig> contains the largest diversity of features and depends on several independent scripts.
we used the same scripts and options used in micropred.
however, we used rnafold from viennarna- <dig>   <cit>  and unafold v <dig>   <cit> , instead of the older versions used in micropred.
initially, we obtained implausible values for features based on pair-probabilities.
they were linked to the function get _pr() from the rna perl package of viennarna- <dig> .
we bypassed this problem by restarting perl for each new sequence.
fs <dig> and fs <dig> were obtained from fs <dig>  a customized python script computed the randfold  and the mfe with the randfold  <cit>  package and merged these two features with fs <dig> to obtain fs <dig>  fs <dig> was obtained merging fs <dig> and the seven additional features obtained using the python script obtained from the authors’ website.
algorithms
the algorithms we adopted have different learning biases.
this is important for the present work, since learning biases can play in favor of a feature set over others.
svm and rf are the two most applied algorithms for pre-mirna classification, whereas g 2de offered class distribution interpretability.
similar interpretation would be obtained using rvkde but  <cit>  showed that rvkde produced accuracies similar to g 2de and slightly lower than svm, even though the number of kernels constructed by each algorithm were on average  <dig> ,  <dig>  and six .
support vector machines
svms deal with classification tasks by finding a hyperplane that separates training instances from two different classes with the maximum margin.
the examples used to determine the hyperplane are the support vectors.
because many problems are not linearly separable, for these problems, the original feature space is mapped into a higher-dimensional space, where linear separation becomes feasible.
points from the original space are mapped to the new space by a kernel function.
the rbf  kernel is a reasonable choice as it performs well for a wide range of problems  <cit> .
for the training of svms, we used a python interface for the library libsvm  <dig>   <cit> .
this interface implements the c-svm algorithm using the rbf kernel.
the kernel parameters γ and c were tuned by 5-fold cross validation  over the grid 2− <dig> − <dig> …,215×2− <dig> − <dig> …, <dig>  the pair  that led to the highest cv accuracy was used to train the svms using the whole training set.
the induced model was then applied to the corresponding gen test set.
random forest
rf is an ensemble learning algorithm that induces a set of decision trees based on the concepts of “bagging” and random feature selection.
bagging is an important approach to estimate generalization error, whereas the latter is important to generate tree diversity.
it was shown  <cit>  that the strength of the ensemble depends on the strength of individual trees and the correlation between any two trees in the forest  <cit> .
the number of features affects the strength of individual trees as well as the tree diversity, while the number of trees affects the generalization error.
in order to obtain an ensemble with lower generalization error, a sufficiently large number of trees shall be chosen, taking into consideration two facts: rfs do not overfit, but limit the generalization error.
this means that the number of trees has to be large enough to ensure lower generalization error, but after a certain value it does not have any effect on the generalization error estimate.
for our experiments, we adopted the r package randomforest <cit> .
each ensemble was generated over the grid ×∗d, representing respectively the number of trees and the number of features.
the  is the default number of features tried in each node split and d is the dimension of the feature space.
we chose the ensemble with the lowest generalization error over the grid and applied it to the corresponding gen test set.
generalized gaussian density estimator
g 2de  <cit>  was designed to predict an instance class based on the probability density functions  of both positive and negative classes.
each pdf is fitted as mixture of generalized gaussian components, using a limited user-defined number of components.
one important feature of g 2de is to provide the coefficients and parameters associated with these generalized components  <cit> .
the learning process of g 2de involves the estimation of the pdf parameters of each class, in addition to the weights of each component.
if k is the maximum number of components and the feature space has dimension d, the number of parameters will be k/ <dig>  an evolutionary optimization algorithm finds the solution by maximizing the number of instances correctly classified in the training set plus the likelihood of class distributions.
it requires two user-defined parameters: the number of gaussian components  and the number of individuals for the initial population in the genetic algorithm .
the first was kept six as in  <cit> , and n was set to 100k, instead of 10k.
high values of n implicate in more running time.
on the other hand, it is expected that high n increases chances of findings an optimal solution.
since the solution is not deterministic, we ran g 2de five times and chose the solution which gave the highest cv accuracy.
the number five was determined by us through computational experiments.
feature selection
once the prediction model is induced, the highest computational cost in the evaluation of putative pre-mirnas is the feature extraction from the sequences to be classified.
since this procedure is performed on millions of sequences, we performed feature selection on different data sets, excluding features which depend on shuffled sequences, which have a higher extraction cost.
we also analyzed the importance of the whole set of  <dig> features obtained by combining all features from fs <dig> and fs <dig>  and three features from fs <dig>  the analyses of feature importance was performed using the feature importance estimated by randomforest  <cit> , and the f-score, as described in  <cit> .
briefly, randomforest estimates the importance of a feature x
i
by computing the difference between the number of correctly classified out-of-bag  vectors before and after the permutation of x
i
in those vectors, during the training phase.
for example, in the data used for this study, the number of oob vectors was approximately  <dig>  if an ensemble correctly classifies on average  <dig> and  <dig> oob vectors before and after the permutation of x
i
, the estimated importance is approximately  <dig>  this value indicates that x
i
was crucial for the correct classification.
however, the interpretation must consider that the importance is conditioned to the induced ensemble.
thus, another feature x
j
with importance  <dig> in the ensemble including x
i
could obtain much higher importance in another ensemble excluding x
i
.
nevertheless, this measure provides a criterion to evaluate the relevance of each feature given the whole set.
differently, the f-score estimates the ratio of between and within classes distances and is computed before the learning step.
features with higher f-scores are more likely to be more discriminative, even though there is no objective criterion to decide on a specific score cut-off.
in our experiments, we trained svms eliminating features with f-score below different score thresholds.
RESULTS
effect of feature sets and training algorithm
data dimensionality may affect the learning process, particularly for parametric models.
in our experiments, g 2de only converged to a predictive model for the feature set fs <dig>  which has only seven features.
this algorithm uses a genetic algorithm to estimate the parameters of gaussian components and their corresponding weights.
in order to obtain the individuals for the initial population, the genetic algorithm uses another algorithm which generates random covariance matrices.
in our experiments, this algorithm generated non-positive definite matrices, which caused the non convergence of g 2de for higher dimensions.
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
±
predicted accuracies , sensitivities , specificities , f-measures  and mathew correlation coefficients  of classifiers trained with  <dig>  examples, presented as the mean and standard deviation .
capital letters in columns indicate the performance cluster of each feature set, within algorithm .
lower case letters in columns indicate the cluster of each algorithms, within feature sets.
bold numbers represents the highest performances, which were not significantly different according to the clustering criteria in  <cit> .
nevertheless, it must be stressed that the highest predictive performances by the classifiers were not significantly different when the feature sets fs <dig>  fs <dig>  and fs <dig> were used .
as previously mentioned, fs <dig> and fs <dig> are subsets of fs <dig>  while fs <dig> merged fs <dig> with seven additional features.
these characteristics together with the very similar results obtained when using fs <dig>  fs <dig>  and fs <dig> suggest that the increase in the number of features leads to a limited increase of the predictive performance, even though the additional features were shown to be distinct features of pre-mirnas  <cit> .
feature discrimination and feature selection
initially, we analyzed the importance of each one of the  <dig> features, obtained by merging fs <dig>  fs <dig> and three features from fs <dig>  when they were all used to induce classification ensembles by rf.
in parallel, we also computed the f-scores.
the pearson correlation coefficient between the averages of these two measures was 75%, showing a relatively high correlation between importance and discrimination.
figure  <dig> shows the features whose importance was considered higher than  <dig>  it can be seen that only  <dig> features obtained average importances higher than  <dig>  corresponding to approx.
6% of the  <dig> oob training instances.
this results suggested that most of the  <dig> features may be redundant or irrelevant.
interestingly, the features depeding on shuffled sequences appeared among those with the hightest importance or f-score.
however, these features were not included in our feature selection step, due to their high computational cost and redundancy to the selected features.
moreover, since the features sets fs <dig> and fs <dig> share the  <dig> features of fs <dig> and they all produced classifiers with the highest predictive performances, we assumed that the relevant features for pre-mirna classification were among these  <dig> common features.
therefore, the feature selection was performed using fs <dig>  eliminating zd, which depends on shuffled sequences.
the features selected  from this set, in order of relevance, were: mfei <dig>  mfei <dig>  dg, dq, df, nefe, diff, ds, ds/l, |g-c|/l, |g-u|/l, %/stems and mfei <dig>  interestingly, six features of this set are energy-based measures .
the other relevant features are: entropy , compactness of the tree graph representation , two thermodynamical features , normalized frequencies of g-c and g-u pairing /stems).
the computational cost of each feature set was estimated by the computation time for a data set composed of  <dig> pre-mirnas sequences randomly sampled from mirbase  <dig>  among the feature sets that produced the highest classification performance, fs <dig> had the highest cost , followed by fs <dig> , fs <dig>  and fs <dig> .
because select and fs <dig> do not contain any stability measure, their costs are significantly lower than the cost for fs <dig>  they took 2: <dig> min and  <dig> s to be computed.
among the sequence-structure based feature sets, fs <dig> took 3: <dig> min to be computed, whereas fs <dig> took  <dig> s.
the next comparison evaluates how the predictive performance associated with each feature set is affected by the use of different classifiers.
for fs <dig>  fs <dig>  and fs <dig>  the maximum difference in sensitivity and specificity between svm and rf was  <dig> % and  <dig> %, respectively.
the predictive performances of the three classifiers using fs <dig> were very similar.
thus, the learning biases of the three learning algorithms did not seem to have a significant effect on the predictive performance.
this small effect of the learning biases is explained by the use of a sufficiently large training set, since most learning algorithms present a clear difference in their predictive performance only when small training sets are used.
comparison with tools using the same algorithms and feature sets
in order to compare our results with tools used as references for our experiments, we predicted gen test sets with those tools.
their main characteristics are summarized in table  <dig>  in table  <dig> we show the predictive performance on the gen sets obtained by the triplet-svm, mipred, micropred, g 2de, mirident and huntmi classifiers in our experiments.
according to table  <dig>  except for the g 2de tool, which uses the g 2de algorithm, the predictive performance values obtained were much lower than the published values or the sensitivity was compromised by the specificity, or vice versa.
the comparison in table  <dig> used test sets obtained from mirbase  <dig>  whereas the classifiers in those were induced with sequences of older releases.
as the representativeness of the pre-mirna population increases in newer releases, it is likely that the underlying distribution of the positive class would also changes.
therefore, the low sensitivities obtained by tools trained with old releases of mirbase, such as triplet-svm, micropred, mipred and mirident are not surprising.
however, the low specificity values obtained by micropred and huntmi, when compared to other older tools, were not expected.
the specificity obtained by huntmi in  <cit>  was 72%, while we obtained 94% in our experiments, using the same algorithm and feature set.
in contrast, the corresponding sensitivities were 99% and 88%.
likewise, the specificity of micropred  <cit>  was 68%, while we found 95% using the same algorithm and feature set.
different results are usually obtained for experiments ran by different research groups, but not so different.
bp = number of base pairs on the stem, mfe = minimum free energy of the secondary structure, noml = no multiple loops, rr = removed redundancies, e-value ≤  <dig> = expected value in blastn against mirbase, expval = only experimentally validated precursors and rf = random forest.
results are presented as the mean and the standard deviation .
acc = accuracy; se = sensitivity; sp = specificity; fm = f-measure; mcc = mathew correlation coefficient.
the largest loss in specificity is observed for micropred and huntmi tools, which both correct for class imbalance.
that is, they attempt to correct the bias due to training sets formed of imbalanced number of examples in each class.
this imbalance may cause a bias towards the majority class in the learning algorithm.
ideally, the class imbalance correction would increase the sensitivity without dropping the specificity.
micropred increased sensitivity from 83% to 90%, while the specificity dropped slightly from 99% to 97%.
the imbalance rate, the ratio of positive to negative examples, in the data set was 1: <dig>  likewise, huntmi reported sensitivity and specificity values of 94% and 95%, working under an imbalance rate of 1: <dig>  since the ideal imbalance rate is determined experimentally, it is plausible that the class imbalance correction methods applied by micropred and huntmi caused generalization problems.
a contributing effect, or alternative explanation, might be that huntmi uses negative sequence which differ greatly from positive examples, whereas the negative sequences used in micropred and in our experiments were selected to be similar to positive sequences.
the observed loss of specificity might be countered with modifications to the training procedures.
however, the lack of generalization of micropred was also mentioned in  <cit> .
g+c content effect
g+c content is an important feature during the folding of hairpin-like rna sequences.
because g+c-rich sequences have more alternative high-energy stable binding-pairs, the prediction of the corresponding secondary structure is more complex.
we drew the slopes of sensitivity and specificity correspondents to  <dig> %-g+c content quantiles, to have a picture of the predictive performance of the feature sets and the algorithms in predicting g+c-rich sequences.
figure  <dig> shows that the variation in specificity throughout the intervals is random.
nevertheless, the sensitivities depended on the feature set and on the algorithm.
all feature sets dropped the sensitivities of rf classifiers in g+c-rich pre-mirnas.
however, when fs <dig>  fs <dig>  fs <dig> and select  were used to train svm classifiers, only random variations in sensitivity along the  <dig> %-g+c content quantile intervals were obtained.
these four feature sets have %g+c related features in common, such as mfei <dig> and normalized frequencies of g-c and g-u pairing /stems).
as figure  <dig> shows, except for mfei <dig>  the other features appear with relatively low importance in the induction of ensembles by rf.
on the other hand, the support vectors from the svm model contain all the features used.
these results confirmed the importance of including %g+c related features to detect g+c-rich pre-mirna.
scope of the investigation
the research reported in this paper was carried out using human sequences, one of the species with the highest abundance of positive sequences.
assuming that the larger the amount of positive sequences, the larger the amount of information about the human pre-mirna population, the investigation performed with human sequences allowed a more fair comparison of the features sets and learning algorithms.
nevertheless, as it has been indicated that the rising of novel mirnas is highly correlated with morphological complexity  <cit> , our results may vary for more distantly related species.
CONCLUSIONS
a considerable part of the computional cost involved in pre-mirna prediction is due to the feature extraction from candidate sequences.
aiming to recommend effective and less costly sets of features, we investigated the discriminant power of seven rna feature sets, under controlled sources of variation.
throughout extensive computational experiments, we showed that feature diversity is an important requirement in pre-mirna recognition.
nevertheless, despite the discriminant power of individual features, higher dimensional sets did not produce higher classification performance classifiers.
based on these results, we proposed a smaller and less costly to compute subset of features, which produced classification performances as high as the produced by higher dimensional and more expensive sets.
because we attempted to avoid all possible sources of bias, we believe that the maximum classification performances reported here are the state-of-the-art for pre-mirna prediction.
since these maximum classification performances are below experimentally feasible rates, other approaches to increase classification performance are welcome.
as our tests showed, the tools used as references in our work either obtained low accuracies or the sensitivities or specificities were compromised.
competing interests
the authors declare that they have no competing interests.
authors’ contributions
as and ac conceived and supervised the study.
il assembled the data, implemented the scripts, ran the experiments and summarized the results.
the three authors wrote and approved the final manuscript.

 a microrna is a small non-coding rna molecule that modulates the stability of mrna targets and their rate of translation into proteins. mirnas are present in the genome of vertebrates, plants, algae and even viruses and are involved in diverse and complex biological processes, like development and cell differentiation, tumorigenesis and immunity. they can also alter plant gene expression in response to environmental stresses. in animals, maturation of canonical mirnas occurs in two steps: first, the long primary mirna transcript is processed within the nucleus into a ∼60– nucleotides stem-loop hairpin precursor by the enzyme drosha. afterwards, within the cytoplasm, the enzyme dicer cleaves the pre-mirna into a double stranded rna duplex and into a loop. the loop is degraded as a by-product, whereas the rna duplex is unwound by helicase activity, releasing the mature mirna and the star sequence. the last is typically degraded whereas the mature mirna guides the microribonucleo-protein complex to target messengers rnas by partial sequence complementarity. machine learning in mirna recognition differently from protein-coding genes, ncrna genes do not contain easily detectable signals in the sequence level. therefore, computational pipelines for mirna discovery rely on characteristic features of pre-mirnas. coupled with comparative genomic computational pipelines, using rnaseq read libraries or completely ab-initio methods, machine learning algorithms have played an important role in mirna discovery. ml algorithms induce models which are able to predict novel mirnas based on patterns learned from known pre-mirna sequences and from other rna hairpin-like sequences, such as pseudo pre-mirnas, transfer rna and mrnas. promir, a probabilistic method, searches for pre-mirna in genomic sequences using sequence and structure features. naïve bayes based probabilistic models were proposed to score mirnas and pre-mirna candidates. support vector machines, random forest, relaxed variable kernel density estimator and generalized gaussian components based density estimation were used to induce classifiers for pre-mirna prediction. a limitation of working directly with the sequence is that the information available may not be sufficient to infer accurate models for novel mirna prediction. for example, the approximate number of rnas folding into hairpin-like secondary structures in the human genome, without filtering based on phylogenetic conservation, was estimated as million. filters derived from other pre-mirna features reduced the number of pre-mirna candidates to around. therefore, the use of features that consider different aspects may reduce false positive rates in mirna detection. feature sets investigated in the literature some of the features commonly extracted from rna sequences for pre-mirna recognition may not help to distinguish between positive and negative classes. therefore, the feature sets considered may have an important effect in the learning process. several feature sets have been proposed for pre-mirna recognition. to comparatively evaluate the effect of these different sets, we investigated seven feature sets proposed in the literature, named here fsi, i∈{..,7}. they were used to induce six classifiers and they contain most of the features employed in computational pipelines for pre-mirna discovery. next, we briefly define each of these sets. the first set, named fs, has sequence and structural features. the second feature set fs which corresponds to a subset of features of fs was used by to induce a classifier, named micropred. the third feature set fs was used to induce a classifier called g 2de and is composed by seven features also present in fs fs is a set of sequence-structural features used by the triplet-svm. fs is a set of sequence-structure motifs used by another classifier, mirident. fs is the feature set used by the mipred. this feature set merged fs the minimum free energy of folding and a stability measure. finally, fs merged fs with four features of fs and three other features: percentage of low complexity regions detected in the sequence, maximal length of the amino acid string without stop codons and cumulative size of internal loops. it was used to induce the huntmi. the classification performances reported by the mentioned tools are among the highest in previous works. their approximate specificities and sensitivities are: triplet-svm, mipred, micropred, g 2de, mirident and huntmi. proposal and key findings in this study, we investigated the discriminative power of seven rna feature sets, previously used in six tools developed for pre-mirna prediction. among them are two sets composed of sequence-structure features and five sets are a miscellany of rna features. the investigation of a specific feature set, using a particular training data and learning algorithm, may insert learning biases. as a consequence, the predictive performance for other training sets and learning algorithms could be different. to deal with this problem, we evaluated each feature set using the learning algorithms, svms, rf and g 2de, which were used in the publications proposing those feature sets. according to the experimental results, the miscellaneous feature sets produced more accurate predictive models than features sets composed from only sequence-structure patterns. however, the differences in accuracy among miscellaneous feature sets are small or insignificant, despite their large differences in composition and dimensionality. inspired by these results, we selected a subset of features, of lower computational cost, but with a similar classification performance, when compared with fs1-fs and fs6-fs feature sets. the classes of positive and negative test sets used in the experiments presented in this paper were predicted by the tools that we used as reference. except for one tool, higher sensitivity was tied to lower specificity and vice versa. methods our goal was to investigate the predictive performance of rna features in distinguishing pre-mirnas from pseudo hairpins. as such, we adopted seven feature sets and three learning algorithms. the feature sets were used to induce classifiers for pre-mirnas in triplet-svm, mipred, micropred, g 2de, mirident and huntmi. svms were used in triplet-svm, micropred and mirident, whereas rf was used in mipred and huntmi. generalized gaussian density estimator is not a tool for pre-mirna prediction in the sense that the features have to be computed by the user in his/her own pipeline. nevertheless, we included g 2de because of its predictive performance and class distribution interpretability. the subsections below provide details of our experiments. data sets human pre-mirna sequences were downloaded from mirbase as the primary source of positive examples. in an attempt to avoid overfitting, we removed redundant sequences. as such, we clustered the available human pre-mirnas sequences using dnaclust such that sequences within a cluster shared 80% similarity. then, one sequence of each cluster was randomly picked. this yielded the set of positives composed of non-redundant pre-mirnas sequences. the negative examples were the pseudo hairpins from human refseq genes, originally obtained by and subsequently used by. these sequences were obtained in order to keep basic features such as length distribution and minimal free energy of folding, similar to those observed in human pre-mirnas. moreover, this set has no redundant sequences. but, in order to adopt a uniform criterion for redundancy removal, we applied the same procedure adopted in the positive set. only singleton clusters were formed. in practice, it is expected that high similarity between positive and negative examples leads to higher specificity. experiments the predictive performance of any model is dependent on the training set representativeness, which usually increases with the increase in the training set size. typically, density based algorithms, such as g 2de, are more sensitive to the course of dimensionality and larger training sets are more likely to provide higher predictive performances. we determined experimentally the training set size which would be suitable for any algorithm and feature set. each experiment was repeated times, in order to provide standard deviations of each classification performance estimation. one repetition consisted of a test set, named here as gen and training sets. for a given repetition, gen was composed by sequences of each class, which corresponded to 1/ of the non-redundant pre-mirna sequences. the remaining positive and negative sequences were used to sample increments of sequences of each class to compose the training sets of,…, instances. each feature set was computed from the same training and test sets. in total, we worked with test sets and 13× training sets. the classification performances of the three algorithms converged to a threshold for training set sizes equal to, for all feature sets. therefore, we presented the results for the largest training set, which contained sequences. classification performance measures classification performance was measured as accuracy, sensitivity, specificity, f-measure and mathew correlation coefficient; see below. this measures can be computed as given below, such that tp, fn, tn and fp are the numbers of true positives, false negatives, true negatives and false positives, respectively. acc=100×tp+tntp+fn+tn+fpse=100×tptp+fnsp=100×tptn+fpfm=100×2×tp2×tp+fn+fpmcc=100×tp×tn−fp×fn×××) the first three measures are commonly used whereas fm is prefered when a compromise between sensitivity and precision is desirable. mcc measures the correlation between real and predicted classes and it is considered less biased towards class imbalance. we presented the predictive performances by the mean and the standard deviation, over the repetitions. features features used in this work are presented in table with references for the detailed descriptions. the prediction of the secondary structure in this work considered the energy model, as implemented in rnafold and unafold. they predict the structure which gives the minimum free energy of folding. we kept the same parameters used in the original publications. detailed descriptions can be found in the corresponding references. the sequence-structure features combined sequence nucleotide information and its predicted state at the secondary structure. in fs each feature represents the relative frequency of three contiguous nucleotides states at the secondary structure, fixing the middle character. because the triplet-svm script excludes sequences with multiple loops, we implemented a python script to compute fs in any sequence. the motifs in fs give the counts of its occurrence in the sequence-structure string. this string is obtained by padding the nucleotide sequence with its respective predicted state at the secondary structure. this set was obtained using the python script provided in the authors’ website. to compute fs we implemented a python script, based on the micropred perl pipeline. fs contains the largest diversity of features and depends on several independent scripts. we used the same scripts and options used in micropred. however, we used rnafold from viennarna- and unafold v, instead of the older versions used in micropred. initially, we obtained implausible values for features based on pair-probabilities. they were linked to the function get _pr() from the rna perl package of viennarna-. we bypassed this problem by restarting perl for each new sequence. fs and fs were obtained from fs a customized python script computed the randfold and the mfe with the randfold package and merged these two features with fs to obtain fs fs was obtained merging fs and the seven additional features obtained using the python script obtained from the authors’ website. algorithms the algorithms we adopted have different learning biases. this is important for the present work, since learning biases can play in favor of a feature set over others. svm and rf are the two most applied algorithms for pre-mirna classification, whereas g 2de offered class distribution interpretability. similar interpretation would be obtained using rvkde but showed that rvkde produced accuracies similar to g 2de and slightly lower than svm, even though the number of kernels constructed by each algorithm were on average, and six. support vector machines svms deal with classification tasks by finding a hyperplane that separates training instances from two different classes with the maximum margin. the examples used to determine the hyperplane are the support vectors. because many problems are not linearly separable, for these problems, the original feature space is mapped into a higher-dimensional space, where linear separation becomes feasible. points from the original space are mapped to the new space by a kernel function. the rbf kernel is a reasonable choice as it performs well for a wide range of problems. for the training of svms, we used a python interface for the library libsvm. this interface implements the c-svm algorithm using the rbf kernel. the kernel parameters γ and c were tuned by 5-fold cross validation over the grid 2− − …,215×2− − …, the pair that led to the highest cv accuracy was used to train the svms using the whole training set. the induced model was then applied to the corresponding gen test set. random forest rf is an ensemble learning algorithm that induces a set of decision trees based on the concepts of “bagging” and random feature selection. bagging is an important approach to estimate generalization error, whereas the latter is important to generate tree diversity. it was shown that the strength of the ensemble depends on the strength of individual trees and the correlation between any two trees in the forest. the number of features affects the strength of individual trees as well as the tree diversity, while the number of trees affects the generalization error. in order to obtain an ensemble with lower generalization error, a sufficiently large number of trees shall be chosen, taking into consideration two facts: rfs do not overfit, but limit the generalization error. this means that the number of trees has to be large enough to ensure lower generalization error, but after a certain value it does not have any effect on the generalization error estimate. for our experiments, we adopted the r package randomforest. each ensemble was generated over the grid ×∗d, representing respectively the number of trees and the number of features. the is the default number of features tried in each node split and d is the dimension of the feature space. we chose the ensemble with the lowest generalization error over the grid and applied it to the corresponding gen test set. generalized gaussian density estimator g 2de was designed to predict an instance class based on the probability density functions of both positive and negative classes. each pdf is fitted as mixture of generalized gaussian components, using a limited user-defined number of components. one important feature of g 2de is to provide the coefficients and parameters associated with these generalized components. the learning process of g 2de involves the estimation of the pdf parameters of each class, in addition to the weights of each component. if k is the maximum number of components and the feature space has dimension d, the number of parameters will be k/ an evolutionary optimization algorithm finds the solution by maximizing the number of instances correctly classified in the training set plus the likelihood of class distributions. it requires two user-defined parameters: the number of gaussian components and the number of individuals for the initial population in the genetic algorithm. the first was kept six as in, and n was set to 100k, instead of 10k. high values of n implicate in more running time. on the other hand, it is expected that high n increases chances of findings an optimal solution. since the solution is not deterministic, we ran g 2de five times and chose the solution which gave the highest cv accuracy. the number five was determined by us through computational experiments. feature selection once the prediction model is induced, the highest computational cost in the evaluation of putative pre-mirnas is the feature extraction from the sequences to be classified. since this procedure is performed on millions of sequences, we performed feature selection on different data sets, excluding features which depend on shuffled sequences, which have a higher extraction cost. we also analyzed the importance of the whole set of features obtained by combining all features from fs and fs and three features from fs the analyses of feature importance was performed using the feature importance estimated by randomforest, and the f-score, as described in. briefly, randomforest estimates the importance of a feature x i by computing the difference between the number of correctly classified out-of-bag vectors before and after the permutation of x i in those vectors, during the training phase. for example, in the data used for this study, the number of oob vectors was approximately if an ensemble correctly classifies on average and oob vectors before and after the permutation of x i, the estimated importance is approximately this value indicates that x i was crucial for the correct classification. however, the interpretation must consider that the importance is conditioned to the induced ensemble. thus, another feature x j with importance in the ensemble including x i could obtain much higher importance in another ensemble excluding x i. nevertheless, this measure provides a criterion to evaluate the relevance of each feature given the whole set. differently, the f-score estimates the ratio of between and within classes distances and is computed before the learning step. features with higher f-scores are more likely to be more discriminative, even though there is no objective criterion to decide on a specific score cut-off. in our experiments, we trained svms eliminating features with f-score below different score thresholds.  effect of feature sets and training algorithm data dimensionality may affect the learning process, particularly for parametric models. in our experiments, g 2de only converged to a predictive model for the feature set fs which has only seven features. this algorithm uses a genetic algorithm to estimate the parameters of gaussian components and their corresponding weights. in order to obtain the individuals for the initial population, the genetic algorithm uses another algorithm which generates random covariance matrices. in our experiments, this algorithm generated non-positive definite matrices, which caused the non convergence of g 2de for higher dimensions. ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± predicted accuracies, sensitivities, specificities, f-measures and mathew correlation coefficients of classifiers trained with examples, presented as the mean and standard deviation. capital letters in columns indicate the performance cluster of each feature set, within algorithm. lower case letters in columns indicate the cluster of each algorithms, within feature sets. bold numbers represents the highest performances, which were not significantly different according to the clustering criteria in. nevertheless, it must be stressed that the highest predictive performances by the classifiers were not significantly different when the feature sets fs fs and fs were used. as previously mentioned, fs and fs are subsets of fs while fs merged fs with seven additional features. these characteristics together with the very similar results obtained when using fs fs and fs suggest that the increase in the number of features leads to a limited increase of the predictive performance, even though the additional features were shown to be distinct features of pre-mirnas. feature discrimination and feature selection initially, we analyzed the importance of each one of the features, obtained by merging fs fs and three features from fs when they were all used to induce classification ensembles by rf. in parallel, we also computed the f-scores. the pearson correlation coefficient between the averages of these two measures was 75%, showing a relatively high correlation between importance and discrimination. figure shows the features whose importance was considered higher than it can be seen that only features obtained average importances higher than corresponding to approx. 6% of the oob training instances. this results suggested that most of the features may be redundant or irrelevant. interestingly, the features depeding on shuffled sequences appeared among those with the hightest importance or f-score. however, these features were not included in our feature selection step, due to their high computational cost and redundancy to the selected features. moreover, since the features sets fs and fs share the features of fs and they all produced classifiers with the highest predictive performances, we assumed that the relevant features for pre-mirna classification were among these common features. therefore, the feature selection was performed using fs eliminating zd, which depends on shuffled sequences. the features selected from this set, in order of relevance, were: mfei mfei dg, dq, df, nefe, diff, ds, ds/l, |g-c|/l, |g-u|/l, %/stems and mfei interestingly, six features of this set are energy-based measures. the other relevant features are: entropy, compactness of the tree graph representation, two thermodynamical features, normalized frequencies of g-c and g-u pairing /stems). the computational cost of each feature set was estimated by the computation time for a data set composed of pre-mirnas sequences randomly sampled from mirbase among the feature sets that produced the highest classification performance, fs had the highest cost, followed by fs, fs and fs. because select and fs do not contain any stability measure, their costs are significantly lower than the cost for fs they took 2: min and s to be computed. among the sequence-structure based feature sets, fs took 3: min to be computed, whereas fs took s. the next comparison evaluates how the predictive performance associated with each feature set is affected by the use of different classifiers. for fs fs and fs the maximum difference in sensitivity and specificity between svm and rf was % and %, respectively. the predictive performances of the three classifiers using fs were very similar. thus, the learning biases of the three learning algorithms did not seem to have a significant effect on the predictive performance. this small effect of the learning biases is explained by the use of a sufficiently large training set, since most learning algorithms present a clear difference in their predictive performance only when small training sets are used. comparison with tools using the same algorithms and feature sets in order to compare our results with tools used as references for our experiments, we predicted gen test sets with those tools. their main characteristics are summarized in table in table we show the predictive performance on the gen sets obtained by the triplet-svm, mipred, micropred, g 2de, mirident and huntmi classifiers in our experiments. according to table except for the g 2de tool, which uses the g 2de algorithm, the predictive performance values obtained were much lower than the published values or the sensitivity was compromised by the specificity, or vice versa. the comparison in table used test sets obtained from mirbase whereas the classifiers in those were induced with sequences of older releases. as the representativeness of the pre-mirna population increases in newer releases, it is likely that the underlying distribution of the positive class would also changes. therefore, the low sensitivities obtained by tools trained with old releases of mirbase, such as triplet-svm, micropred, mipred and mirident are not surprising. however, the low specificity values obtained by micropred and huntmi, when compared to other older tools, were not expected. the specificity obtained by huntmi in was 72%, while we obtained 94% in our experiments, using the same algorithm and feature set. in contrast, the corresponding sensitivities were 99% and 88%. likewise, the specificity of micropred was 68%, while we found 95% using the same algorithm and feature set. different results are usually obtained for experiments ran by different research groups, but not so different. bp = number of base pairs on the stem, mfe = minimum free energy of the secondary structure, noml = no multiple loops, rr = removed redundancies, e-value ≤ = expected value in blastn against mirbase, expval = only experimentally validated precursors and rf = random forest. results are presented as the mean and the standard deviation. acc = accuracy; se = sensitivity; sp = specificity; fm = f-measure; mcc = mathew correlation coefficient. the largest loss in specificity is observed for micropred and huntmi tools, which both correct for class imbalance. that is, they attempt to correct the bias due to training sets formed of imbalanced number of examples in each class. this imbalance may cause a bias towards the majority class in the learning algorithm. ideally, the class imbalance correction would increase the sensitivity without dropping the specificity. micropred increased sensitivity from 83% to 90%, while the specificity dropped slightly from 99% to 97%. the imbalance rate, the ratio of positive to negative examples, in the data set was 1: likewise, huntmi reported sensitivity and specificity values of 94% and 95%, working under an imbalance rate of 1: since the ideal imbalance rate is determined experimentally, it is plausible that the class imbalance correction methods applied by micropred and huntmi caused generalization problems. a contributing effect, or alternative explanation, might be that huntmi uses negative sequence which differ greatly from positive examples, whereas the negative sequences used in micropred and in our experiments were selected to be similar to positive sequences. the observed loss of specificity might be countered with modifications to the training procedures. however, the lack of generalization of micropred was also mentioned in. g+c content effect g+c content is an important feature during the folding of hairpin-like rna sequences. because g+c-rich sequences have more alternative high-energy stable binding-pairs, the prediction of the corresponding secondary structure is more complex. we drew the slopes of sensitivity and specificity correspondents to %-g+c content quantiles, to have a picture of the predictive performance of the feature sets and the algorithms in predicting g+c-rich sequences. figure shows that the variation in specificity throughout the intervals is random. nevertheless, the sensitivities depended on the feature set and on the algorithm. all feature sets dropped the sensitivities of rf classifiers in g+c-rich pre-mirnas. however, when fs fs fs and select were used to train svm classifiers, only random variations in sensitivity along the %-g+c content quantile intervals were obtained. these four feature sets have %g+c related features in common, such as mfei and normalized frequencies of g-c and g-u pairing /stems). as figure shows, except for mfei the other features appear with relatively low importance in the induction of ensembles by rf. on the other hand, the support vectors from the svm model contain all the features used. these results confirmed the importance of including %g+c related features to detect g+c-rich pre-mirna. scope of the investigation the research reported in this paper was carried out using human sequences, one of the species with the highest abundance of positive sequences. assuming that the larger the amount of positive sequences, the larger the amount of information about the human pre-mirna population, the investigation performed with human sequences allowed a more fair comparison of the features sets and learning algorithms. nevertheless, as it has been indicated that the rising of novel mirnas is highly correlated with morphological complexity, our results may vary for more distantly related species.  a considerable part of the computional cost involved in pre-mirna prediction is due to the feature extraction from candidate sequences. aiming to recommend effective and less costly sets of features, we investigated the discriminant power of seven rna feature sets, under controlled sources of variation. throughout extensive computational experiments, we showed that feature diversity is an important requirement in pre-mirna recognition. nevertheless, despite the discriminant power of individual features, higher dimensional sets did not produce higher classification performance classifiers. based on these results, we proposed a smaller and less costly to compute subset of features, which produced classification performances as high as the produced by higher dimensional and more expensive sets. because we attempted to avoid all possible sources of bias, we believe that the maximum classification performances reported here are the state-of-the-art for pre-mirna prediction. since these maximum classification performances are below experimentally feasible rates, other approaches to increase classification performance are welcome. as our tests showed, the tools used as references in our work either obtained low accuracies or the sensitivities or specificities were compromised. competing interests the authors declare that they have no competing interests. authors’ contributions as and ac conceived and supervised the study. il assembled the data, implemented the scripts, ran the experiments and summarized the results. the three authors wrote and approved the final manuscript. 
some feature sets are composed of sequence-structure patterns commonly found in pre-mirnas, while others are a combination of more sophisticated rna features. small or non-significant differences were found among the estimated classification performances of classifiers induced using sets with diversification of features, despite the wide differences in their dimension. analysis of putative pre-mirnas using a computationally expensive feature set would be wasteful or even unfeasible for large genomes. inspired in these results, we obtained a lower-dimensional feature set, which achieved a sensitivity of 90% and a specificity of 95%. the results obtained by the tools used as references in the experiments carried out showed that five out of these six tools have lower sensitivity or specificity. the material to reproduce the main results from this paper can be downloaded from http://bioinformatics.rutgers.edu/static/software/discriminant.tar.gz. in this work, we analyze the discriminant power of seven feature sets, which are used in six pre-mirna prediction tools. the analysis is based on the classification performance achieved with these feature sets for the training algorithms used in these tools. in mirna discovery the number of putative mirna loci is in the order of millions. in this work, we propose a relatively inexpensive feature set and explore most of the learning aspects implemented in current ab-initio pre-mirna prediction tools, which may lead to the development of efficient ab-initio pre-mirna discovery tools.
 computational discovery of micrornas is based on pre-determined sets of features from mirna precursors. some feature sets are composed of sequence-structure patterns commonly found in pre-mirnas, while others are a combination of more sophisticated rna features. in this work, we analyze the discriminant power of seven feature sets, which are used in six pre-mirna prediction tools. the analysis is based on the classification performance achieved with these feature sets for the training algorithms used in these tools. we also evaluate feature discrimination through the f-score and feature importance in the induction of random forests.  small or non-significant differences were found among the estimated classification performances of classifiers induced using sets with diversification of features, despite the wide differences in their dimension. inspired in these results, we obtained a lower-dimensional feature set, which achieved a sensitivity of 90% and a specificity of 95%. these estimates are within % of the maximal values obtained with any feature set while it is times faster to compute. even compared to another feature set, which is the computationally least expensive feature set of those from the literature which perform within % of the maximal values, it is times faster to compute. the results obtained by the tools used as references in the experiments carried out showed that five out of these six tools have lower sensitivity or specificity.  in mirna discovery the number of putative mirna loci is in the order of millions. analysis of putative pre-mirnas using a computationally expensive feature set would be wasteful or even unfeasible for large genomes. in this work, we propose a relatively inexpensive feature set and explore most of the learning aspects implemented in current ab-initio pre-mirna prediction tools, which may lead to the development of efficient ab-initio pre-mirna discovery tools. the material to reproduce the main results from this paper can be downloaded from http://bioinformatics.rutgers.edu/static/software/discriminant.tar.gz. 