Methods of biological network inference for reverse engineering cancer chemoresistance mechanisms.

Drug Discovery Today Volume 19, Number 2 February 2014

REVIEWS

Reviews INFORMATICS

Methods of biological network inference for reverse engineering cancer chemoresistance mechanisms Paola Lecca Centre for Integrative Biology – University of Trento, Laboratory of Computational Oncology, via Sommarive 14, 38123 Trento, Italy

We review recent Bayesian network inference methodologies we developed to infer genetic and metabolic pathways associated to oncological drug chemoresistance. Bayesian inference is supported by a rigorous and widely accepted mathematical formalization of predictive analytics. It is an inherently integrative approach allowing the incorporation of prior knowledge and constraints. Moreover, it is recommended to treat noisy data, and large amount of data whose dynamics laws are mostly unknown. We focus on variational Bayesian methods for the inference of stochastic reaction processes and we present a compendium of the recent results of inference of gene and metabolic networks presiding at the development of pancreas cancer resistance to gemcitabine. Introduction Development of resistance in chemotherapeutic agents is a big concern in cancer patients. Chemoresistance is due to several mechanisms that decrease drug cytotoxicity. The majority of these mechanisms are still poorly understood. The current attempts to reduce chemotherapy resistance are based on assumptions about the various candidate mechanisms. The inability to identify and understand the molecular mechanisms associated to the pathways entailing with the resistance is a reason of the low clinical outcomes of the current strategies to mitigate chemotherapy resistance. Recently, some studies showed that genome expression may predict response to drugs and patient outcome [1–3]. Moreover, the study of genes that influence drug activity and toxicity, known as pharmacogenomics, could offer the possibility of tailoring therapy to the specific profile of individual patients and tumors. A pharmacogenetic approach can thus potentially increase response rates and survival outcome while decreasing toxicity and overall treatment costs. Nevertheless, the current pharmacogenomic studies describe genes determining the sensitivity and resistance to chemotherapy drugs, but, at the best of our knowledge, these studies did not address the correlations among these genes. With regard to drug metabolism, experimental studies devoted to the direct measurements of concentration profiles of metabolites and metabolizing enzymes, as well as a plethora of E-mail addresses: [email protected], [email protected].

computational procedure supporting the identification of drug targets, have led to the creation of reliable, even if incomplete, knowledge of the network describing the biotransformations and the mechanism of action of drugs. What is still lacking in the panorama of network pharmacology [4] is an integrative methodology able to (i) infer with a general mathematical proceeding both gene network and metabolic network, (ii) simulate them, and (iii) assemble them into a larger network showing the correlations between genes and metabolism. The use of such a result is of immediate understanding. In fact, integrative network inference methods approaches could support and even guide the experiments devoted to the identification of mechanisms of action of oncological drugs. Integrative computational inference procedures are expressions of the recent paradigm of algorithmic systems biology [5], that propels a computational approach to infer new biological knowledge, and dissect the complexity and understand the time evolution of biological networks. Computational integrative network inference procedures, when applied to unexplored biochemical processes, can also act as generators of hypotheses about the association of genes to metabolizing enzymes and thus can help in saving time and reducing the cost of the wet activities. Many software platforms and applications are facing this challenge. One of the recent large-scale projects in this direction is Ingenuity [6]. Ingenuity systems are web-based application for analyzing and interpreting the biological meaning of huge amount of genomics data. In particular Ingenuity Knowledge Base

1359-6446/06/$ - see front matter ß 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.drudis.2013.10.026

www.drugdiscoverytoday.com

151

REVIEWS

Reviews INFORMATICS

contains biological and chemical interactions and functional annotations created from millions of individually modelled relationships between proteins, genes, complexes, cells, tissue, drugs, and diseases. Furthermore, Ingenuity systems offer also analytical tools for data statistics, data modelling and model inference. The goal of network inference is to deduce topological and – in case of availability of appropriate data – also the outative causal relationships among the components of a biological system from a corpus of data describing such a system [7–10]. Namely, the topological and causal structures of a biological process are often represented by network and visualized as a graph. The nodes of the graph/network represent the components of the biological systems (molecules, metabolites, proteins, etc.). An edge connecting two nodes indicates the existence of an interaction between the biological entities represented by these nodes; the orientation of an edge denotes the cause-effect relationship between two nodes or between two interactions or pathways. The interpretation of the edge orientation depends on the type of network. There are no formal syntax and formal semantics for the specification of a biological systems as a network, however there are some very common symbols familiar to and accepted by the majority of biologists and modelers. For instance, in gene networks, oriented edges are often accompanied by a sign (þ or ) that indicate activation or inhibition respectively, while in metabolic or protein–protein networks, oriented edges represent directions of chemical reactions. Dashed or dotted arrows with differently shaped tips, such as ! or a (or accompanied by a sign þ or ) are frequently used to indicate activation and inhibition of biochemical reactions and pathways. Many cartoons depicting biological networks use this notation, we refer the reader to some of the most known networks and pathways databases, such as Biocarta [11], KEGG Pathways Database [12], REACTOME [13], NCBI BioSystems [14], and Pathway Commons [15]. The data used to infer biological networks concern time series and/or steady state of the abundance of the system’s components such as genes, molecules, proteins, metabolites, enzymes. Namely, metabolic reactome data, measurements of the metabolites concentration and rate of drug influx and efflux, as well as gene expression data are the input of the advanced techniques of network inference applied to pharmacology. Once a network is inferred we dispose of a pharmacokinetics model deeply rooted in real in vivo or in vitro experimental data. The network is a model that can be specified in a mechanistic way and whose dynamic can be simulated through algorithms implementing the kinetics of biochemical interactions in case of metabolic network, or the dynamic of gene interaction in case of gene networks. When applied to pharmacology and in particular to drug discovery and design process, network inference methodologies have to be necessarily integrative [16]. Namely, in order to be predictive and reliable they need to integrate different types of data, mainly gene expression levels and metabolic data, because the response to a pharmacological treatment is directed by the genetic profile of the patient and by its individual metabolic and pharmacokinetic rates. Furthermore, nowadays there is stronger and stronger demand by the biologists and modelers to integrate in the network inference methods also simulation algorithms. Alternating runs of network inference and network simulation allows to iteratively evaluate the effectiveness of the inference methods at recovering 152



models of complex gene networks and complex pharmacokinetics in cases of a limited amount of data [17]. Furthermore, by simulating the inferred network it is possible to get the time behavior also of the unobserved system’s components. In this way the simulated time series of all system components together with the inferred network used as an a priori knowledge, can be used for a further run of network inference, in an iteratively process. Finally, software tools for network inference equipped with procedure for the integration of different types of inferred network are highly demanded [18] by the community of modelers and biologists. The researchers employed in the identification of chemoresistance mechanisms are particularly interested in disposing of procedure and tools for data integration upstream and/or downstream of the network inference procedure. While available sequence of gene expression data are rapidly provided by high-throughput experiments, our current knowledge of the interactions among genes responsible for the tumor chemoresistance and chemosensitivity is still poor. Furthermore, although recent advancements in experimental high-throughput metabolomics and proteomics allowed the collection of a huge amount data both on single proteins or metabolites concentration profiles and on metabolic pathways [19–22] what we currently know about the relationship between gene networks determining drug sensitivity and/or resistance and networks of drug transporters and metabolizing enzymes constitutes only a small fraction of what remains to be discovered. In situations in which there is a large amount of data, but little theory describing the source of those data (i.e. models or hypotheses), Bayesian inference provides a natural and principled way of including into an inferential probabilistic framework prior knowledge embedded in the observed data to deduce models underlying the observed data themselves. In Bayesian inference all forms of uncertainty are expressed in terms of probability. The Bayesian method starts with the formulation of a model (or hypothesis) that we suppose is adequate to describe the background information and the data. A prior probability is assigned to this model. Finally the Bayes’ theorem is used to evaluate a posterior probability (i.e. a degree of belief) for the model in light of the available data [23]. It is worth noting that the Bayesian approach is not directly concerned with the creative process, how to generate new hypothesis or models. It is concerned mainly with assessing the extent to which models reproduce the available knowledge and data. Its use is then particularly suitable to a present context in which the high quality high-throughput data enable the formulation of a priori knowledge able to initialize and guide the inference procedure. In this paper we focus on pancreatic cancer cell lines treated with gemcitabine. We review a general integrative network inference methods we developed in the last two years in close collaboration with experts in pharmacokinetics and pharmacodynamics (Veltkamp et al. [24]) and guided by the methods proposed by mathematicians and computer scientists (Samoilov et al. [25]) that can be used to infer both the network of interactions among genes responsible of the gemcitabine resistance and the network of gemcitabine metabolism in pancreatic cancer cells. The inferential framework combines variational Bayesian methods recently proposed by Lawrence et al. [26] with correlation-based methods presented by Lecca et al. [28,29].

REVIEWS

Time series of gene expression level

dynamics of the cellular biochemical reactions the stochasticity play an significant role [37]. At the end of the run of variational inference step followed by the simulation of the inferred network, the inferred networks (genetic and metabolic) are given as an input to a correlation-based inference proceeding. This proceeding finds out the relationship between the nodes of the two networks and outputs the final integrated network that correlates the genetic network to the metabolic one. A scheme of the whole work-flow of the procedure we are proposing in this study is shown in Fig. 1.

Variational Bayesian inference of stochastic reaction processes In Bayesian inference evidence or observations are used to update or to newly infer the posterior probability that a hypothesis may be true. Variational Bayesian methods are a family of techniques for approximating intractable integrals arising in Bayesian inference. These methods are used in complex statistical models consisting of observed variables as well as unknown parameters and latent variables, with various sorts of relationships among the three types of random variables. Typically, the parameters and latent variables are grouped together as unobserved variables. Variational Bayesian methods are primarily used for two purposes: To provide an analytical approximation to the posterior probability of the unobserved variables, in order to do statistical inference over these variables.

Concentration (nmol/L) 0 100 200 300 400 500 600

This network inference approach includes a simulation step in which the dynamics of the inferred networks are estimated, and then used together with the inferred networks as the input of a correlation-based inference method designed for merging the gene network and the metabolic network into a larger model. The inputs of the inference procedure are gene expression levels of four genes involved in gemcitabine transport and metabolism. These genes are hENT1, dCK, RRM1 and RRM2 and have been identified to be implied in the development of acquired gemcitabine resistance in pancreatic cancer cells lines (see recent studies of Y. Nakano) [30], H. Fujita et al. [31], and S. Ohhashi et al. [32]. The time series of the main gemcitabine metabolites have been measured in recent experiments of S. Veltkamp et al. [24]. We use also synthetic time series of the enzymes dCK and RR obtained by a mechanistic model of gemcitabine metabolism recently developed by the author and co-workers [33,34] and showing a very good agreement with the available experimental data. The core of the inference engine is a variational approximation of Bayesian inference for stochastic processes [26,35,36] described by the diffusion approximation equation of Markov jump processes. The inference engine works in alternation with a simulation algorithm. This algorithm solves the stochastic differential equations of the diffusion processes to obtain the dynamics of each components of the systems, also of those that are not observed. Here, the choice of an inference model for stochastic processes is motivated by the fact that at the single cell level the

Time series of metabolites and enzymes concentrations

0

5

10 15 Time (h)

20

Approximated variational inference Correlation-based inference

Simulation

Gene regulatory and metabolic networks integrate in a larger network

Drug Discovery Today

FIGURE 1

Integrative network inference work-flow: variational approximation of Bayesian inference combined with correlation-based approach is applied to time series input data of gene expression and metabolites concentration profile. Inferred networks are simulated to validated the dynamics of the components of the inferred networks against the experimental data. www.drugdiscoverytoday.com

153

Reviews INFORMATICS


REVIEWS


Input 1 (mandatory) : time series of the expression level of observed genes and metabolizing enzymes associated to the gemcitabine resistance

Reviews INFORMATICS

Rate equations of diffusion approximation of Markov Jump Processes

Plugging priors in variational inference models

f(G1, G2, …, GN)

Input 2 (if available): a priori knowledge about the correlation between chemoresistance (IC50) and a combination of gene expressions f(G1, G2, …, GN)

Building prior models

IC50 of oncological drug

Approximated variationl Bayesian inference

Minimizing

Output: rate equations of

Divergence between trial posterior and true posterior

the posterior processes Drug Discovery Today

FIGURE 2

Scheme depicting inputs, operations and algorithms, and outputs of the approximated variational inference methods used in this study. Adapted from Lecca et al. [28].

To derive a lower bound for the marginal likelihood (some-

times called the ‘evidence’) of the observed data (i.e. the marginal probability of the data given the model, with marginalization performed over unobserved variables). This is typically used for performing model selection, the general idea being that a higher marginal likelihood for a given model indicates a better fit of the data by that model and hence a greater probability that the model in question was the one that generated the data. In the following sections we describe in the detail the mathematical structures and the steps of our inference methods. Figure 2 illustrates its workflow. The inputs obtained in wet lab experiments and discussed in detail later in this article, are of two types: INPUT 1, that is mandatory and consist of the time series of gene expression levels and metabolites and enzymes concentration; and INPUT 2, that is optional and consists of the experimentally observed relationships between the drug 154


concentration that inhibits the cell proliferation by 50% – IC50 – and a mathematical combination of gene expression levels. The inputs are used to build a prior model of the data based on diffusion process equations. Then we fit this model in the variational approximated Bayes’ rule to find the rate equations of the posterior processes through a minimization of the divergence between a trial posterior and the true posterior.

A formulation of variational inference We follow the notation and the proceeding of A. Ruttor et al. in [26] in deriving a formulation of an approximate inference procedure for Markov jump processes modeling biochemical systems. Let discretize the time into small intervals of length Dt and consider this discretized sample path: X0:K fXðt 0 Þ; Xðt 0 þ DtÞ; . . . ; Xðt 0 þ 2DtÞ; . . . ; Xðt 0 þ KDtÞg; where K 2 N .


REVIEWS

so that we obtain log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ

l¼1

where the superscript () denotes the index of the species, for example,

l¼1

ðX0:K ÞlogQðX0:K Þ dX0:K

Yð1Þ ðtÞ fY ð1Þ ðt 1 Þ; Y ð1Þ ðt 2 Þ; . . . ; Y ð1Þ ðt N Þg

(4)

From the Bayes theorem in Eqn 1 we obtain the following equation

is the vector of the measurements of the variable Y for the species 1. According to this notation X(tl), with l = 1, 2, . . ., N is the vectors of the number of individuals of each species at time tl, and thus ð jÞ X(tl) 2 X0:K. Similarly, Yl fYl g with j = 1, 2, . . ., d. From the Bayes theorem, the posterior probability of the trajectory X0:K is Q Pprior ðX0:K Þ PðY j jXðt l ÞÞ Pposterior ðX0:K Þ ¼ Pmarginal ðY1 ; Y2 ; . . . ; YN Þ

(1)

where P(Yj|X(tl)) is the likelihood to observe Yl given the model X(tl), and the marginal likelihood (also termed evidence) is as follows: Z Y N Pmarginal ðY1 ; Y2 ; . . . ; YN Þ ¼ ½ PðYl jXðt l ÞÞPprior ðX0:K Þ dX0:K : (2) l¼1

Except for very simple models, the integration in Eqn 2 (i.e. the marginalization of the integration variable X0:K) is intractable. Dealing with this intractability is the main challenge of Bayesian inference, for estimating both the posterior distribution and the marginal likelihood. When the marginalization is not tractable and we revert to approximation to recover the posterior, we refer to this as approximate inference. There exist many approximation methods and algorithmic procedures (e.g. the Laplace’s approximation, variational inference, Markov Chain Monte Carlo, and some others). We refer the reader to [27] to a comprehensive explanation and comparison of the advantages and limitations of these different techniques and we focus on variational inference. According to variational approaches to approximate inference of the posterior, when the marginal likelihood is analytically intractable, it is common practice to introduce a trial posterior Q(X0:K) that we would like to match with the true posterior. Q(X0:K) is an undefined probability distribution, by which we can multiply and divide the argument of the integral in Eqn 2. Considering the logarithm of the marginal likelihood, we obtain Z

N Y QðX0:K Þlogf½ P

Pprior ðX0:K Þ gdX0:K ðYl jXðt l ÞÞ QðX0:K Þ Z Z N Y ¼ QðX0:K Þlogf½ PðYl jXðt l ÞÞPprior ðX0:K Þg dX0:K Q

YðtÞ fYð1Þ ðtÞ; Yð2Þ ðtÞ; . . . ; YðdÞ ðtÞg

log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ ¼ log

Z

QðX0:K Þ½

N Y P l¼1

Pprior ðX0:K Þ dX0:K : ðYl jXðt l ÞÞ QðX0:K Þ

(3)

We can manipulate Eqn 3 by using the Jensen’s inequality that states the following:

N Y Pprior ðX0:K Þ PðYl jXðt l ÞÞ ¼ Pposterior ðX0:K ÞPmarginal ðY1 ; Y2 ; . . . ; YN Þ: l¼1

(5) Now, substituting 5 in 4 gives log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ Z ¼ QðX0:K Þlog½Pposterior ðX0:K Þ dX0:K Z þ QðX0:K Þlog½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ dX0:K Z QðX0:K ÞlogQðX0:K Þ:

(6)

Note that the second term of the sum on the right-hand side of Eqn 6 is the expectation of [Pmarginal(Y1, Y2, . . ., YN)] under distribution Q(Xo:K), and consequently Z

QðX0:K Þlog½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ dX0:K ¼ log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ:

Therefore log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ log½Pmarginal ðY1 ; Y2 ; . . . ; YN Þ Z þ QðX0:K Þlog½Pposterior ðX0:K Þ dX0:K Z QðX0:K ÞlogQðX0:K Þ dX0:K which means that the difference between the lower bound of the log-marginal likelihood obtained using the Jensen’s inequality and the posterior probability is Z QðX0:K ÞlogQðX0:K Þ dX0:K KL½QðX0:K ÞjjPposterior Z QðX0:K ÞlogPposterior ðX0:K Þ dX0:K :

(7)

This is known as the Kullbach-Leibler (KL) divergence between the two distribution Q(X0:K) and Pposterior(X0:K). Using the Markov property to represent Q(X0:K) and Pposterior(X0:K) we obtain

QðX0:K Þ ¼

K Y

QðXk jXk1 ÞÞ

k¼1

Z log

f ðxÞ dx log f ðxÞ dx;

if

f ðxÞ > 0;

8x2R

Pposterior ðX0:K Þ ¼

K Y

Pposterior ðXk jXk1 Þ

k¼1


155

Reviews INFORMATICS

Let assume to have noisy observations of the systems at N time points for each of the d species composing the system. Let Y(t) be the vector of the experimental observations of the variable Y, that is,

REVIEWS


and therefore the Kullback-Leibler divergence in Eqn 7 can be expressed as follows KL½QðX0:K ÞjjPposterior ðX0:K Þ ¼ ½

Z

K Z X ½ dXk1 QðXk1 Þ k¼1

dXk QðXk jXk1 Þlog

QðXk jXk1 Þ : Pposterior ðXk jXk1 Þ

(8)

Reviews INFORMATICS

The optimal Q is chosen to minimize the KL divergence [35,26]. Once the networks have been inferred by the available data, we simulate their dynamics with Gillespie stochastic simulation algorithm [41] and obtain the time series of all the nodes. The inferred networks along with the time series of all their nodes are the input of a correlation-based network inference method recently developed by the author and consultable at the reference [29]. This method computes a time-lagged correlation between the species of the system and uses this parameter both as a measure of the strength of the relationships between species and/or genes and as a parameter useful to determine the cause–effect relationships between species and interaction. The output of the inference in Lecca et al. [29] is a network showing the correlations and the cause-effect relationships (in the form of oriented edges) between the component of the gene networks and those of the metabolic network.

The case study: the resistance to gemcitabine in pancreas tumors Gemcitabine (20 , 20 -difluorodeoxycytidine, dFdC) is an anti-cancer chemotherapy drug used for the treatment of variuos types of tumors, such as pancreas cancer, non-small cell lung cancer, bladder cancer, soft-tissue sarcoma, metastatic breast cancer [42,24]. Although it has been proposed that the genes for gemcitabine transport and metabolism are involved in the mechanism of cellular resistance to gemcitabine, it is not fully understood if there is an interplay among these genes and how they influence the gemcitabine transport, metabolism, and therefore resistance mechanisms. Few years ago Nakano et al. [30] established various gemcitabine-resistant subclones of human pancreatic cancer cell lines, and investigated changes in gene expression associated with gemcitabine transport and metabolism. Nakano and co-workers identified four genes responsible for the development of gemcitabine resistance: hENT1, dCK, RRM1, and RRM2. Through quantitative real time light cycler PCR analysis, Nakanos et al. [30] showed that the balance of these gene expression levels is associated with inherent and acquired resistance to gemcitabine in pancreatic cancer cells. Moreover, they proved that resistance of cancer cells to gemcitabine is determined by the ratio of these gene expression levels, but not predicted by that of a single gene. Finally, they concluded that the ratio of gene expression might be useful as a predictive marker for the efficacy of gemcitabine therapy in pancreatic cancer patients. Two years later, in 2009, Ashida et al. [43] confirmed the results of Nakano et al [30]. Ashida et al. [43] used a recent technique known as ultrasound-guided fine-needle aspiration biopsy (EUSPNA) to obtain pancreatic tissues. Then they used focused DNA 156


array (FDA) to simultaneously measure the expression of the putative genes associated to the resistance to gemcitabine. In this study we try to address the same questions faced by of Nakano et al. [30] about the relationships between genes and metabolizing enzymes from a network perspective: we investigate the possibility to infer correlations among the genes of the chemoresistance, and then we try to link this network to the metabolism of the drug. In the following subsection we introduce the reader to the metabolism and mechanisms of action of the gemcitabine and to the mechanisms of resistance. In order to not make the notation too dull, we will not use a different script style to indicate the gene and the enzyme that is codified by the gene. The context in which we mention genes and enzymes will make clear the meaning of the notation. Only in formulas in which both gene and enzyme appear, we will add the postfix ‘_mRNA’ and ‘_enzyme’, respectively, to distinguish them.

Gemcitabine metabolism Gemcitabine is an anticancer nucleoside analog in which the hydrogen atoms on the 20 -carbon of deoxycytidine are replaced by fluorine atoms. As with fluorouracil and other analogues of pyrimidines, the triphosphate analogue of gemcitabine replaces one of the building blocks of nucleic acids, in this case cytidine, during DNA replication. The process arrests tumor growth, as only one additional nucleoside can be attached to the ‘faulty’ nucleoside, resulting in termination of DNA replication and ultimately leading the cells to apoptosis. In this study we refer to the recent model proposed by Veltkamp et al. [24]. The model is depicted in Fig. 3. The gemcitabine is transported into cells by equilibrative and concentrative nucleoside transporter hENT1. Then, it is phosphorylated by deoxycytidine kinase (dCK) to become the gemicitabine mono-phosphate (dFdC-MP). dFdC-MP is phosphorylated to its active diphosphorylated (dFdC-DP) and triphosphorylated (dFdC-TP) forms with the intervention of nucleoside monophosphate kinase (NMPK) and nucleoside diphosphate kinase (NDPK), respectively. The triphosphate metabolite (dFdC-TP) competes with the natural nucleoside triphosphate for the incorporation into the DNA and blocks cells in the early DNA synthesis phase. Gemcitabine is also rapidly metabolized by cytidine deaminase to 20 , 20 -difluorodeoxyuridine (dFdU), which is phosphorylated to its monophospate form (dFdU-MP) by thymidine kinase-1 (TK1), to its disphosphate (dFdU-DP) by NMPK and to triphosphate (dFdU-TP) by NDPK. Recently the activity of dFdU-TP has been associated with the cytotoxic effect of the drug [24].

Chemoresistance Many mechanisms contribute to gemcitabine chemoresistance in pancreatic tumor cells [44,43,45,30,46–49]. In recent studies performed on human cancer cell lines, hENT1 was found to be the major gemcitabine transporter [50,51]. If gemcitabine is not transported into the cell via hENT1 it cannot inhibit cell growth [47], but increased hENT1 abundance facilitates efficient cellular entry of gemcitabine and confers increased cytotoxicity [52]. The phosphorylation of gemcitabine by deoxycytidine kinase (dCK) is a rate-limiting step. Deficiency in dCK activity has been considered to be one of the main mechanisms responsible for the develop-


REVIEWS

CDA

Extracellular dFdC

dFdU

hENT1

em

lm Cel

Reviews INFORMATICS

hENT1

ne bra CDA

dFdU

dFdC dCK

TK1

dFdC-MP dCMPD

NMPK

CDP RR dCDP

dFdU-MP

NMPK dFdU-DP

dFdC-DP

NDPK

NDPK

NDPK dFdU-TP

dFdC-TP

Cytoplasm

dCTP

DNA Drug Discovery Today

FIGURE 3

Biotransformation and reaction of the metabolic network of gemcitabine (dFdC) [24].

ment of resistance to gemcitabine. Another factor in gemcitabine resistance is the overexpression of ribonucleotide reductase (RR). Ribonucleotide reductase is mainly responsible for the conversion of ribonucleosides to deoxyribonucleoside triphosphates (dNTPs), which are essential for DNA polymerisation and repair (see Nakano et al. [30] for a quick review about these results). RR consists of the dimerised large and small subunits, M1 and M2, respectively. The M1 subunit possesses a binding site for enzyme regulation (regulatory subunit), and the M2 subunit is involved with RR activity (catalytic subunit). Nakano et al. [30] studied eight human pancreatic adenocarcinoma cell lines, PK1, PCI43,KLM1, PK8, PK9, MIAPaCa2, KP1N, and BxPC3, and generated gemcitabine-resistant cells by exposing the PCI43, PK1, and KLM1 cell lines to incrementally increasing gemcitabine concentrations starting at 3 nM. As the cells adapted to the drug, they doubled the gemcitabine concentration, and they cultured intermediate resistant variants for at least 4weeks. In the Nakano et al. paper [30] the cell lines were named as follows: G for gemcitabine, followed by the nM concentration at which the cell line grew logarithmically. Nakano et al found that most resistant variants were PCI43-G4000, PK1-G4000, and KLM1-G4000 and were resistant to continuous exposure to gemcitabine at 4000 nM. The experiments [30] were performed using cells in the exponential phase of growth. Ashida et al. [43] added three more genes to the list of genes codifying for metabolic and catabolic enzymes: cytidine deaminase (CDA), 50 -nucleotidase (50 -NT), and deoxycytidylate

deaminase (DCD), that, in order to keep as much simple as possible the diagram and our analysis, are not included in Fig. 3 and in this study.

The experimental data Public data of the time courses of mRNA concentrations of the genes implied in resistance to gemcitabine and concentration profiles of the metabolites and the metabolizing enzymes have been processed to obtain the gene network and the metabolic network.

Gene expression level data We used the time series of the expression level of dCK, RRM1 and RRM2 obtained by Whitfield et al. [53] and stored in public database Cyclebase [54]. ‘A priori’ knowledge about the gene network Recent experimental observations [50,30] provide useful information about the correlations between the genes responsible for the resistance to gemcitabine. We can use these information as prior knowledge to incorporate in our inferential framework. Chemoresistance is expressed as the drug concentration that inhibits cell proliferation by 50% (IC50 values) and is usually determined from concentration-effect relationship [30]. Nakano and coworkers [30] found that the expression ratio R¼

hENT1 dCK RRM1 RRM2

(9)


157

REVIEWS


dCK

10 hENT1 8

–

–

6

R

RRM1

Reviews INFORMATICS

4

RRM2 Drug Discovery Today

FIGURE 5

2 0 0

100

200

300

400

IC 50 Drug Discovery Today

Hypothetical correlation network of gene responsible for the resistance to gemcitabine. The gene hENT1 is uncorrelated to all the other gens, whereas both RRM1 and RRM1 are supposed to be negatively correlated (as indicated by sign ‘’) to dCK to reflect the experimental observations of Nakano et al. [30]. According to these observations, in the development of gemcitabine resistance either RRM1 or RRM2 expression increase simultaneously to a decrease of expression level of dCK.

FIGURE 4

Chemorsistance to gemcitabine IC50 is correlated to the expression ratio R in Eqn 9. This plot has been obtained by Nakano et al. [30].

is an informative marker for predicting and monitoring the responses of pancreatic cancer patients to gemcitabine. This expression ratio significantly correlates with the resistance to gemcitabine in pancreatic cancer cells, including also acquired gemcitabine resistant cells (Fig. 4). A decrease of the ratio R indicates inherent and acquired chemoresistance [30]. We can use this result to build an initial hypothesis about the correlations among genes involved in the resistance mechanisms. Eqn 9 tells us that product hENT1 dCK is correlated to the product RRM1 RRM2. Now, to find correlations between each single gene we take into consideration the experimental observations of Nakano et al. [30] about each single gene mRNA profile. We summarize these observations in the following. hENT1 In studies carried on in the first decade of 2000, it is speculated that populations of cells with lower hENT1 abundance may be relatively gemcitabine resistant owing to reduced intracellular accumulation of the drug. In fact, pharmacological inhibition of hENT1 has been reported to render the cells gemcitabine resistant. However, the study of Nakano et al. [30] and the majority of the subsequent researches, show that the expression of the hENT1 gene is not reduced in the development of gemcitabine resistance, and did not correlate with IC50 values of gemcitabine in eight pancreatic cell lines. These data are suggesting that hENT1 expression alone does not reflect inherent and acquired resistance to gemcitabine. In this regards, the increase in hENT1 expression in the PKI-G4000, PCI-G4000 and KLM1-G4000 subclones may be explained as a compensatory adaptation to higher degrees of chemoresistance to gemcitabine, that do not belong to the ranges of the clinical settings. dCK The expression of the dCK gene alone does not correlate with sensitivity to gemcitabine in eight pancreatic cell lines. The dCK deficiency is involved in a higher level of acquired 158


resistance, but not in a lower grade of resistant cells or parental cells.1 RRM1 and RRM2 Either RRM1 or RRM2 expression alone correlates with acquired chemoresistance. In the development of gemcitabine resistance, increases in either RRM1 or RRM2 expression, in addition to a decrease in dCK, may indicate that cells are regulated to avoid damages by gemcitabine incorporation in gemcitabine resistant cells. In the light of these observations we can approximate the formula 9 for the chemoresistance in the following way R

dCK RRM1 RRM2

(10)

and thus we can formulate a prior hypothesis about the correlations between genes responsible for gemcitabine resistance. We hypothesize that the gene hENT1 is uncorrelated to all the other genes, whereas both RRM1 and RRM2 are supposed to be negatively correlated to dCK (Fig. 5). The anticorrelation between RRM1 and dCK and between RRM2 and dCK expresses the fact that either RRM1 or RRM2 expression increase simultaneously to a decrease of expression level of dCK [55]. It is worth to mention here that the edges in the network showed in Fig. 5 do not represent physical interactions or functional interactions, but only the expected anticorrelated time behavior of the expression level observed in the development of resistance.

Metabolic data With regard to the metabolism of gemcitabine, we used the time series data on the concentration of the gemcitabine and its metabolites published by Veltkamp et al. [24]. The concentration of the following metabolites has been measured at four time points (0, 4, 12, and 24 hours): extracellular dFdC (dFdC_out), intracellular dFdC, intracellular dFdC_MP, dFdC_DP, dFdC_TP, dFdU, dFdU_MP, dFdU_DP, and dFdU_TP. The time series of CDA and dCMDP have been simulated by a mechanistic stochastic mass

1 Experiments performed on highly gemcitabine resistant clones are quite different from clinical settings.


REVIEWS

TABLE 1

dCK_mRNA

Gene

Gene

RRM1 dCK dCK RRM2 dCK hENT1 RRM1 RRM2 hENT1 RRM1 RRM1 hENT Marginal log-likelihood

a

k

0.004055018 0.235188159 0 0 0 0 6228.899

0.1757378 0.1074255 0 0 0 0

action model developed in [33,34]. The concentration of TK1 has been kept constant and equal to 1. The unit of measurement of concentration is nM. The time series of the concentration of the enzymes dCK and RR have been derived by the gene expression level of the genes dCK and RR available in Cyclebase [54] using the mass action law as in the following: d½dCK enzyme / dCK mRNA dt d½RR enzyme / RR mRNA dt that imply dCK enzymeðtÞ ¼

Z

dCK mRNAðtÞ dt Z RR enzymeðtÞ ¼ RR mRNAðtÞ dt

where dCK_mRNA(t) and RR_mRNA(t) are polynomial functions of time t obtained by fitting a polynomial to the experimental data of expression level of dCK and RR.

Results Gene network Eqn 10 tells us that the expression level of dCK is correlated to the product of RRM1 and RRM2 expression levels by the value of the chemoresistance, but it does not tell anything about a possible relationship between dCK and RRM1 and RRM2 separately. The value of chemoresistance is mainly specific to a patient, and also the steady-state value of the ratio in Eqn 9 is patient-dependent. What we are looking for is a mathematical quantitative relationship between dCK and RRM1 and RRM2 separately, that is not patient-dependent. In order to make the model as simple as possible we assume that the diffusion process drift is linear in the time derivative of RRM1 and RRM2 expression levels. Namely, d½dCK ¼ c1 G1 þ c2 G2 þ WðtÞ dt

(11)

where G1 / d[RRM1]/dt and G2 / d[RRM2]/dt, and c1 and c2 are constants. In Table 1 we report the estimates of the parameters of the posterior process whose rate has been assumed a priori to be as in Eqn 11 and depicted by the network of Fig. 5. These parameters

hENT1_mRNA 0.004055018

0.235188159

RRM1_mRNA

RRM2_mRNA

RR Drug Discovery Today

FIGURE 6

Network inferred by the variational inference method. The edge do not indicate any sort of physical, biochemical or regulatory interaction. In this context they indicate an algebraic relationship between the expression level of the involved genes. The number close to the edges are the inclusion probability of the edges in the model, that is, the existence and the strength of the relationship in the model specification. We found the inclusion probability of the link between dCK and RRM2 is about 60 times higher than the inclusion probability of the link between dCK and RRM1 in the model.

have been obtained with the variational inference methods described in Variational Bayesian inference of stochastic reaction processes section. a is the variational estimate of the posterior inclusion probabilities, averaged over settings of parameters (as in use by [56]). k is the variational estimate of the posterior mean coefficients, averaged over settings of the parameters (as in use by [56]). The a is higher for the couple dCK–RRM2 than for the couple dCK–RRM1. This indicates that the link between dCK and RRM2 is higher probable than the link between dCK and RRM1. Finally the last row of Table 1 reports the value of marginal log-likelihood (i.e. the logarithm of Eqn 2 for the diffusion model in Eqn 11). The variational inference method returns a gene network as in Fig. 6. The graphical representation of the correlations between genes of chemoresistance is given in Fig. 6. Our integrative network inference method predicts that the modulation of expression of dCK influences the expression both of RRM1 and RRM2. However, the posterior inclusion probability of the variable RRM2 has been found to be about 60 times higher than the posterior inclusion probability of RRM1. Both RRM1 and RRM2 are predictors of dCK expression level, as it is well known and confirmed by experiments [31,32,30] and reported in public pathways databases (see for example STRING database) [38], but our algorithm found that RRM2 should be the most probable predictor for dCK. It is known that the downregulation of dCK expression is strongly correlated with acquired resistance to gemcitabine in pancreatic cancer [32,30,31]. As a consequence, dCK is a promising specific enhancer to increase the sensitivity of pancreatic cancer to gemcitabine. The correlation predicted between dCK, RRM1 and RRM2 is in agreement with the recent experiments suggesting RRM1 and RRM2 are involved in the mechanism of resistance to gemcitabine. The difference in the posterior inclusion probability of the interactions dCK–RRM1 and dCK–RRM2 would reflect the hypothesis that the association between dCk and RRM2 is stronger that the association between dCK and RRM1. Although this result www.drugdiscoverytoday.com

159

Reviews INFORMATICS

Estimates of the gene ‘network’ parameters obtained with variational inference methods. a is the variational estimate of the posterior inclusion probabilities, averaged over settings of the parameters. k is the variational estimate of the posterior mean coefficients, averaged over settings of the parameters. The a is about 60 times higher for the couple dCK-RRM2 than for the couple dCK–RRM1. This indicates that the link between dCK and RRM2 is higher probable than the link between dCK and RRM1.

REVIEWS


TABLE 2

Variational Bayesian estimates of the model parameters [28]. Species

Drift

a

k

dFdC _ out0

k1dFdC_out +k2dFdC_in k2dFdC_out k2dFdC_in +k1dFdC_out k4dFdC_in k6dFdC_in dCK k12dFdC_in k4dFdU_in +k3dFdU_out k14dFdU_in TK1 k34dFdU_in CDA k3dFdU_out +k5dFdC_out k1dFdC_out k7 dFdC_MP k8dFdC_MP +k9dFdC_DP k13dFdC_MP dCMPD k9dFdC_DP k10dFdC_MP +k8dFdC_DP k11dFdC_TP +k10dFdC_DP k15dFdU_MP k16dFdU_MP +k17dFdU_DP k17dFdU_DP k18dFdU_DP +k19 dFDU_TP k19dFdU_TP +k18dFdU_DP k20dCK_mRNA k21dCTP dCK k22RR_mRNA k33RR dFdC_DP

a1 =1 a2 = 1 a3 =1 a4 =1 a5 =1 a6 =1 a7 =1 a8 =1 a9 =1.00000000 a10 =1.00000000 a11 =0.05693185 a12 =0.05693185 a13 =1 a14 =1 a15 =1 a16 =0.0001848516 a17 =0.0001848516 a18 =0.0009929985 a19 =0.0100000000 a20 =0.0008894042 a21 =0.0001446555 a22 =0.0008894042 a23 =0.0099918765 a24 =0.0008731788 a25 =0.0008622886 a26 =0.0008622886 a27 =0.0077662533 a28 =0.007722776 a29 =0.007722776 a30 =0.009999999 a31 =0.009999999 a32 =0.007722547 a33 =0.001048428 a34 =0.009991198 a35 =0.0004318596 a36 =0.0009293393

2.215987058 0.007972246 0.007972246 0.1607113 1.0839144 0.4894519 125.0990550 0.6714444 0.11467949 0.01828059 7.69691549 7.69691549 0.008473699 1.203491124 0.120373355 0.04589507 0.04589507 0.19631441 0.00000000 0.07428256 0.01374268 0.07428255 0.0001288694 0.0002772840 0.06357436 0.06357435 0.36750269 2.666065e02 2.666064e02 1.798279e05 7.275922e09 1.078708e05 9.741248e06 3.530087e05 3.433878e03 1.553951e05

dFdC _ in0 Reviews INFORMATICS

dFdU _ in0

dFdU _ out0

dFdC _ MP0

dFdC _ DP0

dFdC _ TP0 dFdU _ MP0

dFdU _ DP0

dFdU _ TP0 dCK0 RR0

has been inferred from a rich set of data obtained in four experiments of M.L. Whitfield et al. in 2002 [54,28], further (still unavailbale) experiments, collecting time-series of mRNA expression levels of dCK, RRM1, RRM2 certainly would contribute to validate this finding. Nevertheless, the correlation statistical analysis of the data reported in Y. Nakano et al. [30] and showing the relationship between gemcitabine resistance and expression levels of dCK, RRM1, and RRM2 mRNA in the process of acquired gemcitabine resistance in gemcitabine-resistant PCI43, PK1, and KLM1 subclones, revealed a statistical significant increment of RRM2 mRNA level corresponding to a significative decrement of dCK mRAN level in PCI43 and PK1 subclones. The anti-correlation between dCK and RRM1 has been found less statistically significant. Furthermore, the inference algorithm predicts the absence of correlation between dCK and hENT1. On the opposite, in E. Giovannetti et al. [39] and J. Spratlin et al. [40] the level of hENT1 expression was found to be correlated with survival in pancreatic cancer patients treated with gemcitabine. Moreover, hENT1 was reported to be a predictive marker for the efficacy of gemcitabine therapy in pancreatic cancer patients. Ohhashi et al. [32] obseved that these data were obtained on the basis of expression analyses. 160


However experimental studies of these authors found no change in hENT1 expression in gemcitabine-resistant pancreatic cancer cells, and no change in cellular gemcitabine sensitivity after inhibition of hENT1 expression. Ohhashi and co-workers explained the reason for this inconsistency stating that it is possible that the level of hENT1 mRNA expression is not directly correlated with the functional level of hENT1 [32].

Metabolic network In Table 2 the first column contains the names of the drug, metabolites and metabolizing enzymes. a is the variational estimate of the posterior inclusion probabilities. k is the variational estimate of the posterior mean coefficients, averaged over settings of the parameters. The symbol 0 indicates the derivative of the species concentration with respect to the time. The second column reports the terms of the rate equation for each species: this equation is the drift of the prior process, to which we added normal noise. The k0 s are the rate constants, such that k 2 Rþ . The second column is the prior knowledge; we did not include the reactions involving CDP, dCDPT, and RR (Fig. 3) we did not include the reactions involving CDP, dCDPT, and RR since for


+

α1

dCK_mRNA

α4

α33 dCK

–

α34

dCTP

α3

dFdC_out

dFdC_in

dFdU_out α13

CDA α8

α9

dFdU_in TK1

α7

α16

dFdC_MP

α11 dCMPD α19

α25

Reviews INFORMATICS

hENT1_mRNA

REVIEWS

dFdU_MP

RR_mRNA α35 RR

α17

α18

dFdC_DP α21

α23

dFdC_TP

α26

α28

dFdU_DP α29

α30

dFdU_TP Drug Discovery Today

FIGURE 7

Gene and metabolic network. The dashed edge connecting hENT1 and dFdC_in indicates non zero correlation between the nodes. This arc has been obtained with the time-lagged correlation method in Lecca et al. [29]. The dotted doubly oriented edge between dCK_mRNA and RR_mRNA is the arc obtained by the variational inference method applied on the genes of the resistance. The thickness of the arrows represents the value of the inclusion probability. Very thick arrows indicate a 1, medium thick arrow indicate a 0.01, and finally thin arrows represent a 0.0001. Adapted from Lecca et al. [28].

this pathways of the metabolic network the synthetic data are more than the experimental ones, and we think that including a prior knowledge mostly based on synthetic data could jeopardize the reliability of the inference. The third column shows for each term of the rate equation of the prior, the inclusion probability and the estimates of the posterior mean coefficients, averaged over settings of the parameters. For each species we report also the value of the marginal log-likelihood. Applying a time-lagged correlation method [29] taking as an input the gene and metabolic networks, and the time series generated by the stochastic simulation of the dynamics of the networks, we obtained the network shown in Fig. 7. The dashed and the dotted oriented arcs represent correlations. The dashed edge connecting hENT1 and dFdC_in indicates non zero correlation between the nodes. This arc has been obtained with the timelagged correlation method in Lecca et al. [29]. The correlationbased infernece procedure of Lecca et al. [29] finds correlations also between the expression level of hENT1 and the phosphorylated metabolites of dFdC_in; however, these are second order correlations due to the stronger correlation between hENT1_mRNA and dFdC_in. The dotted doubly oriented edge between dCK_mRNA and RR_mRNA is the arc obtained by the variational inference method applied on the genes of the resistance. It represents the correlation pathway in Fig. 5. The result in Fig. 7 is correct since it is confirmed by well established experiments.

Conclusions Bayesian variational inference, traditionally used for variable selection in regression problems has been recently proposed also as method of inference of biological processes at molecular scale

[26,28]. This method holds the promises to handle noisy, incomplete data and to be flexible to include prior knowledge. Nevertheless, nowadays, there is not a single method of network inference, able to address all the main questions challenging the intellectual efforts of biologists and modelers, for example, the accuracy in predicting combinatorial regulations in gene networks and/or cooperativity mechanisms in metabolic pathways and, most importantly the ability to be combine different data from different type of biological process to infer networks integrating gene regulatory process with biochemical metabolic pathways, essentially the ability to provide an integrated systems view on the biological processes under investigation. we prospect that a suitable integration of different methods joined to a technique capturing combinatorial, non-linear relationships among genes could be a way towards a reliable and complete inferential scheme. In this review we presented a summary of our recent experience in the development of integrative network inference method to reverse engineer biochemical reaction and gene regulation pathways. The methods discussed here are integrative because it operates on different kind of data coming from different biological systems (gene and metabolites systems), it integrates model simulation in its framework and it assembles into a larger integrated network all the inferred networks thorough a correlation-based inference method coupled to the variational procedure. The approach is integrative for three reasons: it integrates prior knowledge derived from different types of experimental data and observations and formalizes it in a probabilistic model to guide inference procedure; it integrates model simulation in the inference framework to assess the goodness of the inference www.drugdiscoverytoday.com

161

REVIEWS

results; and finally, it assembles into a larger integrated network all the inferred networks. We showed how this method performs on the case study of the gemcitabine chemoresistance. We inferred the correlations between the genes reported as responsible for the resistance, and the main metabolic biotranformations of drug. Finally we found out the correlations between genes and metabolic enzymes.


Strongly integrative approach to the identification of relationships between genes of resistance and drug metabolizing enzymes is recommended. The result of the inference procedure is a mathematical model deeply rooted in experimental data that is able to work as an in silico lab and as a generator of hypothesis about the association of genes to metabolizing enzymes in the desirable perspective of saving time and reducing the cost of the wet activities.

Reviews INFORMATICS

References 1 Nakai, Y. et al. (2005) Identifying genes with differential expression in gemcitabineresistant pancreatic cancer cells using comprehensive transciptome analysis. Oncol. Rep. 14, 1263–1267 2 Sato, J. et al. (2011) Gene expression analysis for predicting gemcitabine resistance in human cholangiocacinoma. J. Hepatobiliary Pancreat. Sci. 18, 700–711 3 Toffalorio, F. et al. (2010) Expression of gemcitabine- and cisplatin- related genes in non-small-cell lung cancer. Pharmacogenom. J. 10, 180–190 4 Hopkins, A.L. (2008) Network pharmacology: the next paradigm in drug discovery. Nat. Chem. Biol. 4, 682–690 5 Priami, C. (2009) Algorithmic systems biology. Commun. ACM 52, 80–88 6 Ingenuity home page: www.ingenuity.com, 2013. 7 Stelniec-Klotz, I. et al. (2012) Reverse engineering a hierarchical regulatory network downstream of oncogenic KRAS. Mol. Syst. Biol. 8, 601 8 Madhamshettiwar, B.P. et al. (2012) Gene regulatory network inference: evaluation and application to ovarian cancer allows the prioritization of drug targets. Genome Med. 4, 41 9 Tanaka, N. and Papoian, G.A. (2010) Reverse-engineering of biochemical reaction networks from spatio-temporal correlations of fluorescence fluctuations. J. Theor. Biol. 264 (2), 490–500 10 Crombach, A. et al. (2012) Efficient reverse-engineering of a developmental gene regulatory network. PLoS Comput. Biol. 8 (7), e1002589 11 Biocarta: http://www.biocarta.com/, 2012. 12 KEGG Pathways Database: http://www.genome.jp/kegg/pathway.html, 2012. 13 REACTOME: http://www.reactome.org/ReactomeGWT/entrypoint.html, 2012. 14 NCBI Biosystem Database: http://www.ncbi.nlm.nih.gov/Structure/biosystems/ docs/biosystems_about.html, 2012. 15 PC Pathways Commons: http://www.pathwaycommons.org/pc/, 2012. 16 Azmi, A.S. , Wang, Z. , Philip, P.A. , Mohammad, R.M. and Sarkar1, F.H. (2010) Proof of concept: a review on how network and systems biology approaches aid in the discovery of potent anticancer drug combinations. Mol. Cancer Ther. 9 (12), 3137–3144 17 (2002) Evaluating functional network inference using simulation of complex systems. Bioinformatics 18, S216–S224 18 Chao, S.-Y. et al. (2011) An integrative approach to identify cancer chemoresistanceassociated pathways. BMC Med. Genomics 4 19 (2012) Quantitative high-throughput metabolomics: a new era in epidemiology and genetics. Genome Med. 4 20 Han, J. et al. (2008) Towards high-throughput metabolomics using ultrahigh-field Fourier transform ion cyclotron resonance mass spectrometry. Metabolomics 4, 128– 140 21 (2008) Metabolomics: a global biochemical approach to drug response and disease. Annu. Rev. Pharmacol. Toxicol. 48, 653–683 22 Heux, S. et al. (2012) A high-throughput metabolomics method to predict high concentration cytotoxicity of drugs from low concentration profiles. Metabolomics 8, 433–443 23 Baldi, P. and Brunak, S. (2001) Bioinformatics. The Machine Learning Approach. The MIT Press 24 Veltkamp, S.A. et al. (2008) New insights into the pharmacology and cytotoxicity of gemcitabine and 20 , 20 -difluorodeoxyuridine. Mol. Cancer Ther. 7, 2415–2425 25 (2001) On the deduction of chemical reaction pathways from measurements of time series of concentration. Chaos 11, 108–114 26 (2010) Approximate inference for stochastic reaction systems. In Learning and Inference in Computational and Systems Biology (1 edn) (Lawrence, N.D., Girolami, M., Rattray, M., eds), MIT Press 27 Lawrence, N.D., Girolami, M., Rattray, M., eds) (2010) Learning and Inference in Computational and Systems Biology., MIT Press 28 Lecca, P. (2013) An integrative network inference approach to predict mechanisms of cancer chemoresistance. Integr. Biol. 5, 458–473

162


29 Lecca, P. et al. (2012) Inferring biochemical reaction pathways: the case of the gemcitabine pharmacokinetics. BMC Syst. Biol. 6 30 Nakano, Y. et al. (2007) Gemcitabine chemoresistance and molecular markers associated with gemcitabine transport and metabolism in human pancreatic cancer cells. Br. J. Cancer 96, 457–463 31 Fujita, H. et al. (2010) Gene expression levels as predictive markers of outcome in pancreatic cancer after gemcitabine-based adjuvant chemotherapy. Neoplasia 12 (10), 807–817 32 Ohhashi, S. et al. (2008) Down-regulation of deoxycytidine kinase enhances acquired resistance to gemcitabine in pancreatic cancer. Anticancer Res. 28, 2205–2212 33 Lecca, P. et al. (2011) Modelling and estimating dynamics of tumor shrinkage with Blenx and kinfer. In Proceedings of UKSim 2011 – 13th International Conference on Modelling and Simulation, Cambridge, UK, August. pp. 75–80 34 Kahramanogullari, O. et al. (2012) Algorithmic modeling quantifies the complementary contribution of metabolic inhibitions to gemcitabine efficacy. PloS One (in press) 35 Archambeau, C. and Opper, M. (2011) Approximate Inference for Continuous-Time Markov Processes. Cambridge University Press 36 (2008) The variational approximation for Bayesian inference. Life after the EM algorithm. Signal Process. Mag., IEEE 25, 131–146 37 Wilkinson, D.J. (2006) Stochastic Modelling for Systems Biology. CRC Press/Taylor & Francis Group 38 STRING database home page: http://string-db.org/. 39 Giovannetti, E. et al. (2006) Transcription analysis of human equilibrative nucleoside transporter-1 predicts survival in pancreas cancer patients treated with gemcitabine. Cancer Res. 66, 3928–3935 40 Spratlin, J. et al. (2004) The absence of human equilibrative nucleoside transporter 1 is associated with reduced survival in patients with gemcitabine-treated pancreas adenocarcinoma. Clin. Cancer Res. 10, 6956–6961 41 Gillespie, D.T. (1976) A general method for numerically simulating the stochastic time evolution of coupled chemical species. J. Comp. Phys. 22, 403–434 42 Chu, E. and DeVita, V.T. (2007) Physicians’ Cancer Chemotherapy Drug Manual. Jones & Bartlett 43 Ashida, R. et al. (2009) Gemcitabine sensitivity-related mRNA expression in endoscopic ultrasound-guided fine-needle aspiration biopsy of unresectable pancreatic cancer. J. Exp. Clin. Cancer Res. 28 44 Andersson, R. et al. (2009) Gemcitabine chemoresistance in pancreatic cancer: molecular mechanisms and potential solutions. Scand. J. Gastroenterol. 44, 782–786 45 Duxbury, M. et al. (2003) Identification of novel channel and transporter genes associated with acquired gemcitabine resistance in pancreatic adenocarcinoma cells. Poster session presented at Pharmaconference, Switzerland 46 Ina, S. et al. (2010) Identifying molecular markers for chemosensitivity to gemcitabine in pancreatic cancer: increased expression of interferon-stimulated gene 15 kd is associated with intrinsic chemoresistance. Pancreas 39, 473–485 47 Rauchwerger, D.R. et al. (2000) Equilibrativesensitive nucleoside transporter and its role in gemcitabine sensitivity. Cancer Res. 60, 6075–6079 48 Kim, M.P. and Gallick, G.E. (2008) Gemcitabine resistance in pancreatic cancer: picking the key players. Clin. Cancer Res. 14, 1284 49 Huanwen, W. et al. (2009) Intrinsic chemoresistance to gemcitabine is associated with constitutive and laminin-induced phosphorylation of fak in pancreatic cancer cell lines. Mol. Cancer 8, 125 50 Garcia-Manteiga, J. et al. (2003) Nucleoside transporter profiles in human pancreatic cancer cells: role of hcnt1 in 20,20-difluorodeoxycytidine-induced cytotoxicity. Clin. Cancer Res. 9, 5000–5008 51 Spratlin, J.L. and Mackey, J.R. (2010) Human equilibrative nucleoside transporter 1 (hENT1) in pancreatic adenocarcinoma: towards individualized treatment decisions. Cancers 2, 2044–2054

52 Ritzel, M.W. et al. (2001) Recent molecular advances in studies of the concentrative Na+dependent nucleoside transporter (cnt) family: identification and characterization of novel human and mouse proteins (hcnt3 and mcnt3) broadly selective for purine and pyrimidine nucleosides (system cib). Mol. Membr. Biol. 18, 65–72 53 Whitfield, M.L. et al. (2002) Identification of genes periodically expressed in the human cell cycle and their expression in tumors. Mol. Biol. Cell 13, 1977–2000

REVIEWS

54 Cyclebase: http://www.cyclebase.org/, 2012. 55 Re´jiba, S. et al. (2009) Gemcitabine-based chemogene therapy for pancreatic cancer using ad-dck::umk gdept and tsrr sirna strategies. Neoplasia 7, 637–650 56 Carbonetto, P. and Stephens, M. (2012) Scalable variational inference for Bayesian variable selection in regression, and its accuracy in genetic association studies. Bayesian Anal. 7, 73–108

Reviews INFORMATICS



163

Network-based inference methods for drug repositioning.

ARACNe-AP: gene network reverse engineering through adaptive partitioning inference of mutual information.

Accelerated parallel algorithm for gene network reverse engineering.

Minimum network constraint on reverse engineering to develop biological regulatory networks.

Molecular mechanisms of chemoresistance in gastric cancer.

Untangling statistical and biological models to understand network inference: the need for a genomics network ontology.

Biological network inference using low order partial correlation.

Enabling network inference methods to handle missing data and outliers.

On protocols and measures for the validation of supervised methods for the inference of biological networks.

Reverse engineering cellular networks with information theoretic methods.

Systems biology: the role of engineering in the reverse engineering of biological signaling.

Matrix Factorization for Transcriptional Regulatory Network Inference.

Instrumental variable methods for causal inference.

Systematic reverse engineering of network topologies: a case study of resettable bistable cellular responses.

Reverse engineering highlights potential principles of large gene regulatory network design and learning.

The feasibility of genome-scale biological network inference using Graphics Processing Units.

The molecular mechanisms of chemoresistance in cancers.

Dynamical and topological robustness of the mammalian cell cycle network: a reverse engineering approach.

Molecular mechanisms of chemoresistance in osteosarcoma (Review).

Assessment of network inference methods: how to cope with an underdetermined problem.

Reverse engineering transcriptional gene networks.

Reverse engineering of multilayer coatings for ultrafast laser applications.

Elastic-Net Copula Granger Causality for Inference of Biological Networks.

Network biology elucidates metastatic colon cancer mechanisms.