ISSN: 1838-7640Theranostics

* Theranostics *
2018; 8(1):277-291.
doi:10.7150/thno.22065 This issue

Research Paper

1. Electrical and Computer Engineering, Purdue University, West Lafayette, IN, USA;

2. Computer Science, Columbia University, New York, NY, USA;

3. Computer Science, Purdue University, West Lafayette, IN, USA.

This is an open access article distributed under the terms of the Creative Commons Attribution (CC BY-NC) license (https://creativecommons.org/licenses/by-nc/4.0/). See http://ivyspring.com/terms for full terms and conditions.

Koo J, Zhang J, Chaterji S. Tiresias: Context-sensitive Approach to Decipher the Presence and Strength of MicroRNA Regulatory Interactions.

MicroRNAs (miRNAs) are short non-coding RNAs that regulate expression of target messenger RNAs (mRNAs) post-transcriptionally. Understanding the precise regulatory role of miRNAs is of great interest since miRNAs have been shown to play an important role in development, diseases, and other biological processes. Early work on miRNA target prediction has focused on static sequence-driven miRNA-mRNA complementarity. However, recent research also utilizes expression-level data to study context-dependent regulation effects in a more dynamic, physiologically-relevant setting.

**Methods**: We propose a novel artificial neural network (ANN) based method, named Tiresias, to predict such targets in a context-dependent manner by combining sequence and expression data. In order to predict the interacting pairs among miRNAs and mRNAs and their regulatory weights, we develop a two-stage ANN and present how to train it appropriately. Tiresias is designed to study various regulation models, ranging from a simple linear model to a complex non-linear model. Tiresias has a single hyper-parameter to control the sparsity of miRNA-mRNA interactions, which we optimize using Bayesian optimization.

**Results**: Tiresias performs better than existing computational methods such as GenMiR++, Elastic Net, and PIMiM, achieving an F1 score of >0.8 for a certain level of regulation strength. For the TCGA breast invasive carcinoma dataset, Tiresias results in the rate of up to 82% in detecting the experimentally-validated interactions between miRNAs and mRNAs, even if we assume that true regulations may result in a low level of regulation strength.

**Conclusion**: Tiresias is a two-stage ANN, computational method that deciphers context-dependent microRNA regulatory interactions. Experiment results demonstrate that Tiresias outperforms existing solutions and can achieve a high F1 score. Source code of Tiresias is available at

**Keywords**: MicroRNA, messenger RNA, machine learning, artificial neural networks, Bayesian optimization, breast cancer, condition-specific regulation.

MicroRNAs (miRNAs) are non-coding RNAs of ~22-nucleotide length that can regulate gene expression by binding to complementary sequences on target messenger RNAs (mRNAs). So far, only part of the hundreds of known human miRNAs have been experimentally linked to their target mRNAs. Accurately determining the entire miRNA targetome is a major bottleneck in miRNA functional characterization, which is important in unraveling disease mechanisms and determining the efficacy of drugs.

The first generation of computational methods developed to predict miRNA targets are primarily based on nucleotide-sequence information [1-4]. Such computational methods mainly utilize the sequence complementarity between the seed region of the miRNAs, mostly situated at nucleotides 2-7 from the miRNA 5'-end, and the 3'-UTR of the mRNAs. However, in the sequence-based genre of miRNA targeting algorithms, it was later realized that sequence complementarity was neither necessary nor sufficient for the miRNA-mRNA interactions to occur, which was further validated by the presence of a vast swath, specifically >90%, of these so-called non-canonical interactions [5]. These interactions were deciphered from the datasets obtained from high-throughput sequencing (CLIP-seq) experiments, such as, PAR-CLIP (Photoactivatable Ribonucleoside-Enhanced Crosslinking and Immunoprecipitation) and HITS-CLIP (High-throughput sequencing of RNA isolated by crosslinking immunoprecipitation). In fact, our recent work where we used both sequence and thermodynamic features at two levels of granularity—seed level (finer granularity) and site level (coarser granularity around the miRNA-mRNA interaction site)—resulted in superior prediction of these non-canonical interaction sites [6]. Our prediction performance improved further when using a non-linear kernel for our classifier, reducing the bias in our simpler linear model, with an improvement of >150% in the true positive rate for non-canonical sites, over the best competitive protocol [7], and validated using cancer-specific data. These methods are simple and useful to build a putative interaction database, but they are far from perfect, leading to many false positives and some false negatives.

Unlike sequence data, expression data is context-specific and, in one sense, dynamic, providing useful clues on regulation effects that may vary depending on conditions, such as time, pathophysiology, and cell-type (spatial) context. Thus, context-dependent regulation can be found by analyzing expression data. Correlation analyses including Pearson correlation coefficient (PCC) were initially attempted with the expression data [8], measuring the strength of interactions between miRNAs and mRNAs. However, correlation-based detection leads to a high error rate, since correlation coefficients of interacting pairs can be small, falling in the measurement error range, unless the regulation relationship is strong. In addition, these methods cannot adequately model the demonstrated joint effects of several miRNAs on a shared mRNA target. To account for the joint effects of multiple miRNAs on an mRNA, multiple linear regression models with regularization have been proposed. Lasso regression [9] and Elastic Net regression [10] fall into this category. However, they can only provide a sparse solution, *i.e.*, a relatively small set of strong interactions.

Recently, utilizing both sequence and expression data is getting much attention. Representative work is a variational inference model considered in GenMiR++ [11, 12]. The aim of this model is to infer the probabilities of putative targets of being real by approximating the posterior with variational inference techniques. A drawback of this approach is that the expectation maximization (EM) algorithm involved in variational inference is computationally expensive. Further, the search space for the parameters is restricted to be in the realm of down-regulation only. A probabilistic regression model in [13], called PIMiM, introduced a group-based regulatory model, based on both sequence and expression data. This method can form groups of miRNAs and mRNAs that may function together, but is hard to optimize in practice due to its large number of hyper-parameters, *e.g.*, the numbers of groups, miRNAs per group, and mRNAs per group. On the other hand, recent research results are reporting that up-regulations also exist [14] and non-linear regulation relationship fits better to real data in some cases [15,16]. However, existing computational methods cannot be easily extended to study non-linear regulation models because the problem formulation will look significantly different. Thus, we need a new method that can reflect the above scenarios.

To address the above challenges, we propose Tiresias in this paper.1 Tiresias is a new explanatory computational method based on artificial neural networks (ANNs) that analyzes a given dataset, and finds miRNA targets under context-specific conditions by incorporating expression-level data into sequence-based prediction results. Tiresias decouples the problem in a manner that makes the learning easier—the first stage estimates miRNA targets and the second stage estimates the regulation weight based on the previous stage outputs. Tiresias considers up-regulation as well as down-regulation, simultaneously. By virtue of the ANN-based approach, Tiresias can easily remove the assumption of linear regulation, extending it to a complex non-linear regulation model. Thus, Tiresias can support multiple variations in modeling, from each of which we can get insights understanding regulation relationship hiding underneath complex biological processes. Tiresias characterizes the density of miRNA-mRNA interactions using one single hyper-parameter unlike the most relevant prior work [13], which uses a large number of hyper-parameters, *e.g.*, to individually control for how many miRNAs belong to a group. This design makes Tiresias a good fit for Bayesian optimization by which we optimize the hyper-parameter in a way that a regression error is minimized.

We compare Tiresias to a variational inference model, GenMiR++ [12], an Elastic Net regression method [10], and a probabilistic regression approach, called PIMiM [13], all of which, like Tiresias, are well suited to high-dimensional biological data. Our experiments using synthetic data show that Tiresias significantly improves upon these prior methods, and can achieve a F1 score > 0.8, that is, a high true positive rate and a low false positive rate, in predicting regulatory interactions, which are stronger than a certain threshold. Using real breast cancer datasets from the Cancer Genome Atlas (TCGA) [17], Tiresias shows a detection rate of up to 82% in recovering the experimentally-validated interactions between miRNAs and mRNAs, even if we assume that regulations may result in small values of PCC for the interacting miRNA-mRNA pairs.

We consider a set of expression-level data that consists of *I* miRNAs and *J* mRNAs. We denote a sample of an *miRNA expression vector* by ** x** = (

,(1)

where *μ _{j}* > 0 is the background mean of

In this paper, we mainly consider the following regulation model2:

,(2)

where *w _{ij}* ∈

Thus, our goal is to find *w _{ij}* and

(3)

subject to *μ _{j}* > 0 for all

Consider a function vector ** f**(

.(4)

The ** f**(

Assuming that the interaction matrix ** s** is given, the regulation weights

(5)

where *λ*_{1} > 0 is a penalty weight. Note in Equation (5) that the first term comes from -ln*p*(** y**|

Since we can find *w _{ij}*'s from

We denote an estimator of ** s** given

(6)

with *a*(*z*) = 1/(1+exp(-*z*)) meaning a sigmoid activation function shown in Figure S1A of Supplementary material, and *u _{ijk}* ∈

(7)

in which *b _{kp}* ∈

All the coefficients in Equation (7) are pre-trained by an autoencoder shown in Figure 2C, in such a way that the autoencoder output gets as close to (** x**,

We train the ANN ** g**(

(8)

where

(9)

and *ρ* < 1 is a sparsity parameter that controls how many of *ŝ _{ij}*'s activate (

Note that, so far, *ŝ _{ij}* does not have much relationship with

We have introduced two ANNs, ** f**(

(10)

which is simply a weighted sum of Equation (5) and Equation (8) with a weight *λ*_{2} > 0.

Note that although the term in Equation (10) is computed at the end of ** f**(

On the other hand, the regulation weights of ** f**(

For training ** h**(

(11)

for a given sample pair of ** x** and

A simple example for the miRNA vector **x**, the mRNA vector **y**, the putative interaction matrix **q**, the interaction matrix **s**, and the matrix **w** of regulation weights when I = 3 and J = 2. Here, '*' denotes a 'don't care' value.

(A) Regulation weight estimation network **f**(**x,s**), where **x** is a sample of the miRNA expression vector, and **s** is the interaction matrix. Training this network without knowing a reasonably good estimate of **s** will not result in any meaningful regulation weights at all. (B) The miRNA target estimation network **g**(**x,y**) that takes **x** and **y** as inputs, and is intended to result in as the estimate of the ground truth interaction matrix **s** given **x** and **y**. Training this network alone just produces a random matrix whose entry is close to 0 or 1, which is far from an appropriate estimate of **s** given **x** and **y**. (C) An autoencoder where the features are extracted by training the network in such a way that and are, respectively, produced as close to x_{i} and y_{j} as possible. (D) Overall architecture of Tiresias. For a given value of ρ, **f**(**x,s**) and **g**(**x,y**) are trained together for the whole data sample pairs of **x** and **y**. The Bayesian optimization tunes the value of ρ based on .

In *h*(*x*,*y*), *ρ* is a hyper-parameter such that *ρ* is a variable that is not learned by training *h*(*x*,*y*), but affects how *h*(*x*,*y*) works.5 A proper choice of *ρ*, *i.e.*, as close as possible to the ground truth value, is the key to achieving good accuracy of our estimation networks.

For this reason, we find the optimal value of *ρ*, denoted by *ρ*^{*}, using the Bayesian optimization framework, as illustrated at the bottom of Figure 2D. We define *ρ*^{*} as the one that minimizes the following on average:

(12)

which is the first term of *c'*(** x,y**) in Equation (10). Roughly speaking,

In effect, we first try a few random guesses for *ρ*, train ** h**(

The decision on the new value of *ρ* to try is made by compromising the trade-off between exploration and exploitation. Exploration means we try new values where is very uncertain, but could lead to improvement. Exploitation means we choose the new value in the vicinity of prior desirable values of *ρ*, *i.e.*, where we have some certainty that is expected to be low. By this approach, we can minimize the number of evaluations of , each of which corresponds to training ** h**(

Since we are analyzing interaction relationships among multiple miRNAs and mRNAs simultaneously, data scaling that standardizes the range of each miRNA or mRNA is important. Without such scaling, the cost function during training is governed by few specific miRNAs or mRNAs that have a broad range or a large magnitude of values, and eventually, we end up with poorly trained ANNs. Thus, all the datasets considered in this paper are scaled appropriately. Details for this are presented in Section S1 of Supplementary material.

The whole procedure to train Tiresias is summarized in Figure 3, assuming that we use the stochastic gradient descent as an optimization method. The **procedure** 1 corresponds to training ** h**(

After we complete training Tiresias through procedures in Figure 3, we use the converged values of unknowns of ** h**(

,(13)

where implies the element-wise product (also known as the Schur product) between two matrices ** a** and

(14)

assuming that *P*(*x*^{(m)},*y*^{(m)}) = 1/*M*.

Note that Tiresias is not predicting an interaction for a new data sample after training, but it figures out the interaction from the dataset by which it is trained. Namely, Tiresias is an explanatory model that analyzes the given dataset and finds out the potential regulatory interactions. Thus, the overfitting to the training data that is an issue of a predictive model is not relevant here [23].

In Tiresias, the regulation model *r _{j}*(

We evaluate Tiresias using both real and synthetic datasets. In real datasets, only part of interacting pairs between miRNAs and mRNAs have been discovered by biological experiments. This makes the real datasets inadequate to compute accuracy metrics for any algorithms, for which perfect ground truth is necessary. Thus, our experiments first start with a synthetic dataset, where we generate a number of ground-truth pairs. In such a controlled environment, we evaluate Tiresias with comparison to existing relevant algorithms in terms of precision and recall, and also examine the effects of internal parameters of Tiresias. After getting a sense of the general trend in the performance of Tiresias, we then apply Tiresias to multiple real datasets. In the experiments with real datasets for several diseases, we will see how well Tiresias performs in deciphering the validated interacting pairs, and what kind of new pairs are newly predicted to be interacting.

We implemented Tiresias using Tensorflow [24]. For the Bayesian optimization, we used an open-source software downloaded from [25]. In order to train, we generated a synthetic dataset of 1000 samples. Each sample of the dataset contains 20 miRNAs and 200 mRNAs. For each mRNA, we randomly selected 10 miRNAs for putative interactions (*i.e.*, *q _{ij}* = 1), and 5 miRNAs among them to regulate the mRNA (

The expressions of mRNAs were generated following Equation (1) and Equation (2) with *x _{i}* ~

A non-linear regulation model consisting of an ANN with two hidden layers. The input layer is , where flatten(**a**) vectorizes a matrix **a**, and (**x···x**) ∈ **R**^{I×J} denotes a matrix repeating the vector **x** J times. In the output layer, denotes the value of r_{j}(**x**,**s**) defined as . The hidden layers are expressed as for l=0,1. Here, for all p,k,l are coefficients to find by a training process.

Figure 5A presents quality measures of recovering interactions between miRNAs and mRNAs for Tiresias, GenMiR++, Elastic Net, and PIMiM.7 For Tiresias, we decided *e _{ij}* as a positive if |

In the figure, it is worth noting that Tiresias achieves a F1 score >0.8 at |*w _{ij}*| ≥ 0.10, which means that in practice, Tiresias can successfully predict regulatory interactions of over a certain strength with a low false positive rate. Meanwhile, other methods show poor performance in the same situations, resulting in F1 scores ranging from 0.2 to 0.6. GenMiR++ is expected to perform better if there exist down-regulations only in the dataset, since it is strictly limited to search down-regulations. The figure should be understood as what GenMiR++ can achieve when applied to the dataset that may contain up-regulations as well. This indicates why we need an algorithm like Tiresias that is designed to consider up-regulations in addition to down-regulations. Elastic Net usually performs poorly unless the regulation effect is significantly strong as can be seen in [26]. PIMiM is reported in [13] to achieve the F1 score of about 0.37. Our experiment confirms the result, showing a similar level of a F1 score. We think that PIMiM may be good for grouping genes, which is its intended purpose, but may not be suitable for predicting miRNA targets.

Figures 5B, 5C and 5D show the F1 score of Tiresias according to the number of samples, miRNAs, and mRNAs. Increase in the number of samples provides Tiresias with more information, and thus it is always helpful for improving the performance, as one may imagine. In the meantime, as we can see from Equation (2) that the number of miRNAs that we are analyzing at a time is the number of terms that are summed in *r _{j}*(

Figure 5E shows how the value of changes according to *K*, the number of features extracted from the autoencoder in Figure 2C. As can be seen, , the average of *c*_{1}(** x**,

(A) Box plots for precision, recall, and F1 score showing performance of Tiresias versus GenMiR++, Elastic Net, and PIMiM. Tiresias significantly outperforms others in all measures except at |*w _{ij}*| = 0.05. (B)-(D) Tiresias's F1 score with increase in the number of samples, miRNAs, and mRNAs (|

Figures 5F and 5G show the variation of according to *λ*_{1} and *λ*_{2}, the weights for the penalty terms in Equation (10). We can see from the figures that *λ*_{1} and *λ*_{2} do not significantly affect . This is because near the end of training ** h**(

In order to test Tiresias on real biological data, we first consider a breast invasive carcinoma dataset (project ID TCGA-BRCA) from the Cancer Genome Atlas (TCGA) portal [17], where 1098 samples of expression data from breast cancer patients are available. Evaluation is done by checking if Tiresias can result in a significant value of *e _{ij}* with a right sign when the

The ground truth interactions between miRNAs and mRNAs for breast cancer were checked using the miRTarBase portal [27], which provides experimentally-validated miRNA-target databases. To increase the reliability of the ground truth, we only counted the pairs of miRNAs and mRNAs that are validated more than 3 times by multiple papers, given the lack of high-confidence validated miRNA-target predictions. This resulted in 24 interacting pairs between 19 miRNAs and 19 mRNAs, *i.e.*, *I* = 19 and *J* = 19. The names of selected miRNAs and mRNAs are shown in Figure 6C. Note here that the absence of an interaction in miRTarBase for other pairs does not necessarily mean that the pairs are not actually interacting. Some of such pairs may be interacting, but are not yet validated by experiments. This implies that in the experiment, we cannot evaluate the precision and recall of Tiresias, because some none-validated interacting pairs may be counted as false positives. For this reason, we focus here on evaluating how well Tiresias can decipher the 24 validated interacting pairs from the expression dataset by a detection rate, which we define as a ratio of the number of detection to the total number of validated interacting pairs.

The ground truth directions (up or down) and strengths of regulations were determined based on PCC, a measure of the linear correlation, between an miRNA and an mRNA for the known interacting pairs, as in other literature [28]. However, it should be noted that although the PCC can decide the regulation strength and direction for the validated interacting pairs, the PCC value itself cannot tell us which pairs are interacting, since non-interacting pairs may also show a significant value of PCC. As shown in Figure 6B, most putative interactions indeed have a non-negligible value of PCC. We will see that Tiresias, which finds interaction relationships based on ** e**, not PCC, can filter out such putative interactions that have significant PCCs.

The putative interaction matrix ** q** was made by sequence-based prediction results, downloaded from TargetScan [2] version 7.1, which includes non-canonical sites. The optimal value of

The regulatory network edge matrix ** e** in Equation (13), estimated by Tiresias, is expressed as a heat map in Figure 6C, where a positive indicates the up-regulation and a negative

Figure 6E that shows Tiresias's detection rate for the validated pairs, after filtering by PCC values in a way that the ground-truth interaction is assumed if PCC ≥ 0.1 or PCC ≥ 0.001. We can see that Tiresias's detection performance is naturally dependent on the regulation strength, *i.e.*, easier to detect if a regulation is strong. However, we can also see that the detection rate of Tiresias is still high enough (>0.8) even if we assume that true regulation relationship may result in a small value of PCC (PCC ≥ 0.001). Figure 6F also shows the detection rate of Tiresias when we increment the numbers of miRNAs and mRNAs under consideration. For this experiment, the candidates of ground-truth interaction pairs before filtering by PCC values are selected if only one paper is reported interacting in miRTarBase, which leads to 119 interacting pairs between 60 miRNAs and 99 mRNAs. As we can see in Figure 5C, the increase in the number of miRNAs means a higher noise level and thus lower performance, which is naturally unavoidable. Thus, we can see from the figure that Tiresias's detection rate is lower than in Figure 6E. Although Tiresias can still achieve ~0.6 detection rate in this case, this suggests us to use Tiresias with a practical limit in the number of miRNAs at a time to build high confidence. We may consider dividing the dataset into multiple ones, each of which contains partial lists of miRNAs, and applying Tiresias to each smaller dataset. By changing the way of separating the dataset, and iterating detection on each partial dataset, we can still see the joint effect among miRNAs. From Figure 6C, we can also see that since PCCs are usually small even for the ground truth pairs (non-black triangle markers in Figure 6C), corresponding values of |*e _{ij}*| are mostly small and thus the detection threshold of Tiresias must be set small enough (<0.01).

(A) Bayesian optimization process, where the utility means the value of the acquisition function by which we can decide the next value of *ρ* to try, based on the history so far. (B) PCCs for TCGA-BRCA dataset. (C) Heat map of the regulatory network edge matrix ** e** estimated by Tiresias for TCGA-BRCA dataset. A putative interacting pair is marked as ○. Down-regulation (of a negative PCC) and up-regulation (of a positive PCC) are marked as ▽ and △, respectively, for the ground truth interactions. For the pairs whose PCC magnitude is >0.1, the triangle markers get larger, proportional to their PCC magnitude, and they are colored black. (D) Regulation strength and direction. Here, the strength in the vertical axis is the PCC values for validated interacting pairs in case of the ground truth, and the estimated values of

We test the normality of mRNAs expressions for the TCGA breast cancer dataset, since we began with assuming that *y _{j}* for any

Figure 6I shows regression error reduction for each mRNA by switching from the linear regulation model in Equation (2) to the non-linear regulation model in Figure 4. Here, the regression error for *y _{j}* is defined as:

,(15)

in which the superscript (*m*) over a value means that the value comes from the *m*-th data sample, as before. We can see that *f _{j}*(

We also applied Tiresias to other real datasets from TCGA, whose project IDs in TCGA site are TCGA-LUAD (Lung Adenocarcinoma), and TCGA-UCEC (Uterine Corpus Endometrial Carcinoma). The TCGA-LUAD dataset has 585 samples, and the TCGA-UCEC dataset consists of 560 samples. The validated interacting pairs were again checked from the miRTarBase portal, but this time we took all pairs into account that were reported interacting at least once in literature, due to their lower popularity. From Figure 5B, performance of Tiresias is expected to be lower than in TCGA-BRCA, since the number of samples available is much lower. Other factors of experimental setup are the same as before, unless otherwise stated.

Figures 7A and 7B shows PCCs for the TCGA-LUAD dataset and the regulatory network edge matrix ** e** estimated by Tiresias, where the value of

Figures 7C and 7D shows PCCs for the TCGA-UCEC dataset and the regulatory network edge matrix ** e** estimated by Tiresias with

In this paper, we proposed Tiresias, a new computational method based on ANNs for predicting which mRNAs are targeted by which miRNAs. Tiresias separately considers the miRNA target estimation and the regulation weight estimation by building an ANN for each problem. However, in the training phase, Tiresias combines them together as a single two-stage ANN, which enables us to associate the two problems and solve them simultaneously.

Our prediction is context-specific, *i.e.*, takes into account what kind of diseased state the cell is in and incorporates expression-level data in addition to sequence-based prediction data. Experimental results showed that Tiresias performs better than existing computational methods such as GenMiR++, Elastic Net, and PIMiM, achieving an F1 score of >0.8 for a certain level of regulation strength. For the TCGA breast cancer dataset, Tiresias showed a true positive rate of up to 82% in recovering the ground truth regulatory interactions between miRNAs and mRNAs, even when we assume that the regulations may lead to small PCCs.

Tiresias predicted new interacting pairs in Figures 6C, 7B, and 7D (non-zero *e _{ij}* for putative interactions). Since Tiresias is expected to have a low false positive rate from Figure 5A, those may be considered suitable candidates for future validation via biological experiments. The full list of miRNA-mRNA pairs that are predicted to be interacting by Tiresias can be found in the Supplementary material.

In order to further improve the F1 score, we can also consider incorporating other regulatory factors such as transcription factors into our model as an input to the miRNA target estimation network ** g**(

(A) PCCs for TCGA-LUAD data set. (B) Regulatory network edge matrix ** e** for TCGA-LUAD, estimated by Tiresias. Note that (hsa-let-7e-5p, HMGA2) and (hsa-miR-7-5p, PK3R3) are detected with significant strength, although their PCCs are negligibly small. (C) PCCs for TCGA-UCEC data set. (D) Regulatory network edge matrix

Note that although Tiresias is an explanatory model that analyzes the given dataset, we can also make use of the trained Tiresias as an mRNA expression level predictor for given miRNA expressions. With the capability to extend into a non-linear regulation model, Tiresias can achieve a highly accurate prediction to the level of mRNA expression.

We have currently used TargetScan for determining the necessary conditions for miRNA-mRNA interactions. This can be replaced with other sequence-based methods that are being actively developed for predicting such putative interactions, including our recent work [7], which is able to look for non-canonical sequence matches as well. A more accurate determination of such priors will further improve the performance of Tiresias.

We have used stochastic gradient descent for training, iterating over the TCGA samples 100 times, at which point, convergence is reached. On a Dell PowerEdge R320 server of 12 Intel Xeon processors, it takes Tiresias ~10 min to complete training for a given value of *ρ*, including feature extraction by an autoencoder. We then use Bayesian optimization to find an optimal value of *ρ*, which in all our experiments, happens within 10 iterations.

miRNA: microRNA; mRNA: messenger RNA; ANN: artificial neural network; PAR-CLIP: photoactivatable ribonucleoside-enhanced crosslinking and immunoprecipitation; HITS-CLIP: high-throughput sequencing of RNA isolated by crosslinking immunoprecipitation; PCC: Pearson correlation coefficient; EM: expectation maximization; TCGA: the Cancer Genome Atlas; pdf: probability density function; MSE: mean squared error; BRCA: breast invasive carcinoma; LUAD: lung adenocarcinoma; UCEC: uterine corpus endometrial carcinoma.

Supplementary figures and tables.

The authors have declared that no competing interest exists.

1.
Griffiths-Jones S, Grocock RJ, Dongen S, Bateman A, Enright AJ. **mirbase: microrna sequences, targets and gene nomenclature**. *Nucleic Acids Res. * 2006;**34**:D140-D144

2.
Agarwal V, Bell GW, Nam J-W, Bartel DP. **Predicting effective microrna target sites in mammalian mrnas**. *eLife. * 2015;**4**:e05005

3.
Betel D, Koppal A, Agius P, Sander C, Leslie C. **Comprehensive modeling of microrna targets predicts functional non-conserved and non-canonical sites**. *Genome Biol. * 2010;**11**(8):R90

4.
Khorshid M, Hausser J, Zavolan M, Nimwegen E. **A biophysical mirna-mrna interaction model infers canonical and noncanonical targets**. *Nat methods. * 2013;**10**(3):253-255

5.
Asirvatham AJ, Gregorie CJ, Hu Z, Magner WJ, Tomasi TB. **Microrna targets in immune genes and the dicer/argonaute and are machinery components**. *Mol Immunol. * 2008;**45**(7):1995-2006

6.
Ghoshal A, Shankar R, Bagchi S, Grama A, Chaterji S. **Microrna target prediction using thermodynamic and sequence curves**. *BMC Genomics. * 2015;**16**:999

7.
Ghoshal A, Grama A, Bagchi S, Chaterji S. *An ensemble svm model for the accurate prediction of non-canonical microrna targets. *
ACM BCB. 2015:403-412

8.
Fulci V, Colombo T, Chiaretti S, Messina M, Citarella F. *et al*. **Characterization of b- and t-lineage acute lymphoblastic leukemia by integrated analysis of microrna and mrna expression profiles**. *Genes Chromosomes Cancer. * 2009;**48**(12):1069-1082

9.
Lu Y, Zhou Y, Qu W, Deng M, Zhang C. **A lasso regression model for the construction of microrna-target regulatory networks**. *Bioinformatics. * 2011;**27**(17):2406-2413

10.
Beck D, Ayers S, Wen J, Brandl MB, Pham TD. *et al*. **Integrative analysis of next generation sequencing for small non-coding rnas and transcriptional regulation in myelodysplastic syndromes**. *BMC Med Genomics. * 2011;**4**(1):19

11.
Huang JC, Morris QD, Frey BJ. **Detecting MicroRNA Targets by Linking Sequence, MicroRNA and Gene Expression Data, RECOMB**.
*2006; 114-129. *

12.
Huang JC, Babak T, Corson TW, Chua G, Khan S. *et al*. **Using expression profiling data to identify human microrna targets**. *Nat Methods. * 2007;**4**(12):1045-1049

13.
Le H-S, Bar-Joseph Z. **Integrating sequence, expression and interaction data to determine condition-specific mirna regulation**. *Bioinformatics. * 2013;**29**(13):i89-i97

14.
Vasudevan S. **Posttranscriptional upregulation by micrornas**. *Wiley Interdiscip Rev RNA. * 2012;**3**(3):311-330

15.
Shu J, Xia Z, Li L, Liang ET, Slipek N. *et al*. **Dose-dependent differential mrna target selection and regulation by let-7a-7f and mir-17-92 cluster micrornas**. *RNA Biology. * 2012;**9**(10):1275-1287

16.
Lai X, Wolkenhauer O, Vera J. **Understanding microrna-mediated gene regulatory networks through mathematical modelling**. *Nucleic Acids Res. * 2016;**44**(13):6019-6035

17.
*The Cancer Genome Atlas. *
http://cancergenome.nih.gov/

18.
*Law of large numbers. *
https://en.wikipedia.org/wiki/Law_of_large_numbers

19.
Dayan P, Hinton GE, Neal RN, Zemel RS. **The Helmholtz machine**. *Neural Computation. * 1995;**7**:889-904

20.
*Ng A. Sparse autoencoder. *
https://web.stanford.edu/class/cs294a/sparseAutoencoder.pdf

21.
Schmidhuber J. **Deep learning in neural networks: An overview**. *Neural Networks. * 2015;**61**:85-117

22.
Brochu E, Cora VM, Freitas N. *A tutorial on bayesian optimization of expensive cost functions, with application to active user modeling and hierarchical reinforcement learning. *
CoRR. 2010 arXiv:1012.2599

23.
Shmueli G. **To explain or to predict?**. *Statist Sci. * 2010;**25**(3):289-310

24.
*TensorFlow: A system for large-scale machine learning. *
https://www.tensorflow.org

25.
*Bayesian Optimization. *
https://github.com/fmfn/BayesianOptimization

26.
Sass S, Pitea A, Unger K, Hess J, Mueller NS, Theis FJ. **Microrna-target network inference and local network enrichment analysis identify two microrna clusters with distinct functions in head and neck squamous cell carcinoma**. *Int J Mol Sci. * 2015;**16**(12):30204-30222

27.
*miRTarBase. The experimentally validated miRNA-target interactions database. *
http://mirtarbase.mbc.nctu.edu.tw

28.
Wang YP, Li KB. **Correlation of expression profiles between micrornas and mrna targets using nci-60 data**. *BMC Genomics. * 2009;**10**:218

29.
*Q-Q plot. *
https://en.wikipedia.org/wiki/Q-Q_plot

30.
Yan X, Su XG. *Linear Regression Analysis: Theory and Computing. * River Edge, NJ, USA:
World Scientific Publishing Co Inc. 2009

31.
Liang H, Li WH. **Microrna regulation of human protein-protein interaction network**. *RNA. * 2007;**13**(9):1402-1408

32.
Hsu CW, Juan HF, Huang HC. **Characterization of microrna-regulated protein-protein interaction network**. *Proteomics. * 2008;**8**(10):1975-1979

Tiresias is a Greek character who was famous for his ability to see precisely into the future. He achieved this without needing to see into all the details of the present. Similarly, our model is unaware of the exact interaction priors but is able to infer them from the data and generate new predictions for miRNA-mRNA interactions.

Later, we will discuss how to extend this to non-linear regulation models.

Although this network does not currently contain any activation function like a sigmoid function, we stick to calling it an ANN in case we replace *r _{j}*(

The features can also be extracted by using multi-layer ANNs, *i.e.*, deeper ANNs in order to reflect more complex non-linearity.

Technically, *λ*_{1} and *λ*_{2} are also hyper-parameters, but we have found by experiments that they are not critical compared to *ρ* (see Figures 5F and 5G). Thus, we can simply set them all to 1, and focus on tuning *ρ*.

Regarding the rareness of up-regulations, the probabilities of generating up-regulation and down-regulation are set to 0.1 and 0.9, respectively.

Note that GenMiR++ and PIMiM were originally designed to find down-regulations only. For this experiment, we modified PIMiM to detect up-regulations as well. However, GenMiR++ was used as is, since this method is inherently restricted to down-regulations by their modeling. Protein interaction matrix required for PIMiM was assumed to be an all-zero matrix, since it is usually very sparse.

A formal definition of the Q-Q plot can be found at [29].

Corresponding author: schaterjiorg

Received 2017-7-24

Accepted 2017-10-9

Published 2018-1-1