Understanding the Mechanism of Cell Death in Gemcitabine Resistant Pancreatic Ductal Adenocarcinoma: A Systems Biology Approach.

BACKGROUND: Gemcitabine is the standard chemotherapeutic drug administered in advanced Pancreatic Ductal Adenocarcinoma (PDAC). However, due to drug resistance in PDAC patients, this treatment has become less effective. Over the years, clinical trials for the quest of finding novel compounds that can be used in combination with gemcitabine have met very little success.
OBJECTIVE: To predict the driving factors behind pancreatic ductal adenocarcinoma, and to understand the effect of these components in the progression of the disease and their contribution to cell growth and proliferation.
METHODS: With the help of systems biology approaches and using gene expression data, which is generally found in abundance, dysregulated elements in key signalling pathways were predicted. Prominent dysregulated elements were integrated into a model to simulate and study the effect of gemcitabine-induced hypoxia.
RESULTS: In this study, several transcription factors in the form of key drivers of cancer-related genes were predicted with the help of CARNIVAL, and the effect of gemcitabine-induced hypoxia on the apoptosis pathway was shown to have an effect on the downstream elements of two primary pathway models; EGF/VEGF and TNF signalling pathway.
CONCLUSION: It was observed that EGF/VEGF signalling pathway played a major role in inducing drug resistance through cell growth, proliferation, and avoiding cell death. Targeting the major upstream components of this pathway could potentially lead to successful treatment.

INTRODUCTION

Pancreatic ductal adenocarcinoma (PDAC) is the most common form of pancreatic cancers, with a high mortality rate of over 95% worldwide [1] with approximately 459,000 new registered cases, of which about 432,000 deaths occured. According to GLOBOCAN [1], the incidence rate for this cancer is highest in Europe and North America, with equal incidences in both male and female population. This makes pancreatic cancer the fourth leading cause of cancer-related death in the world [2, 3]. Besides the lack of early detection and diagnosis, a major obstacle that lies in the pathway of PDAC treatment is the drug-resistant phenotype, which develops overtime [4]. Current forms of treatment including traditional chemotherapy and radiation therapy do little to improve the survival rate, which is the lowest for all cancers at only 9% for the 5-year relative survival rate [5]. A comprehensive understanding of how the tumor microenvironment operates may help enhance therapeutic development and improve the quality of survival.

Gemcitabine is the first line of chemotherapy administered in patients with advanced PDAC [6, 7]. It is a deoxycytidine analogue that inhibits ribonucleotide reductase M1 and M2, preventing the replication of DNA, and thus inhibiting cancer growth. The problem with this form of treatment arises when most of the patients develop resistance towards gemcitabine overtime [4, 8]. Several clinical trials have since been conducted to come up with compounds that can be used in combination with gemcitabine. This effort has, however, been futile as significant improvement has not been observed on overall survival rates. Several clinical trials have since been conducted to come up with compounds that can be used in combination with gemcitabine. This effort has, however, been futile as significant improvement has not been observed on overall survival rates. Some of the drugs currently used in combination with gemcitabine include 5-Fluorouracil (5-FU), irinotecan and folfirinox, which provide around 4.3 months of overall increased survival for the patient [9]. This, however, came with adverse effects. Another combination of gemcitabine with nab-paclitaxel was reported to increase the overall survival by up to 2 months when compared to treatment with gemcitabine alone [10], however this form of treatment does not change the overall survival by a significant amount. Even with the knowledge of the mechanism of action of gemcitabine [11], the overall effect that it has on other important signalling pathways such as the TNF, EGF, and VEGF signalling pathways, and how the mechanism of apoptosis is regulated under the influence of gemcitabine is not understood yet.

Systems biology is an emerging branch of bioinformatics that combines experimental data with mathematical modelling to analyze complex and dynamic biological systems [12]. The main purpose of systems biology is to study the interaction between genes and proteins, and understand how these interactions give rise to specific behavior. Several biological problems can be solved by systems biology through the application of graph theoretical concepts. Thus, by applying such approaches, we aim to provide a better understanding of the mechanism underlying gemcitabine drug resistance in PDAC, which could lead to the development of potential therapeutics for successful treatment in the future.

MATERIALS AND METHODS

Gene Expression Dataset and Analysis

RNA-seq data for the effect of gemcitabine on the PDAC cell line was obtained from publically available dataset [13] at NCBI GEO (GEO ID: GSE105083). PANC-1, a well-known cell line obtained from the ductal region of the pancreas, treated with 100nm and 1µm of gemcitabine was compared against normal, untreated PANC-1 cell line. Gene count data from HTSeq [14] was provided as the input for differential gene expression using DEGUST (http://degust.erc.monash.edu/), an online tool for RNA-Seq exploration analysis, and visualization. Differential gene expression analysis was performed using edgeR method [15], keeping the false discovery rate (FDR), or adjusted p-value cutoff and log FC (fold change) at <0.1 and 1, respectively. A lower cutoff was set for both log FC and FDR to retain more genes since a higher cutoff gave an insufficient number of genes to work with.

Prediction of Transcriptional Factors and Pathway Activities Using CARNIVAL

CAusal Reasoning pipeline for Network identification using Integer VALue programming (CARNIVAL) [16] is a tool that combines gene expression data with prior knowledge network (PKN) to provide a hypothesis on how upstream signalling regulators control downstream targets in a signalling network. Prediction is made based on a casual reasoning approach using Integer Linear Programming (ILP) [17].

To begin with, a PKN was constructed using the resource available from OmniPath [18], a comprehensive database for literature curated signalling pathways. Signalling pathway information for humans was selected, keeping only interactions that were directed and signed. Based on the type of interaction, i.e., activation or inhibition function, a reaction sign of 1 or -1 was assigned, respectively. Using the gene expression data obtained from DEGUST, the activity of transcriptional factors on direct gene targets was estimated using DoRothEA v2 [19]. The top 50 transcriptional factor scores were then stored for further use. For the prediction of pathway activities, PROGENy [20], a tool for inferring activities of signalling pathways from transcriptome data, was used. This method makes use of a linear regression model to predict perturbed pathways. The predicted pathway and transcriptional factor scores were used as input for running CARNIVAL. All parameters were set to default for this step, besides using inverse CARNIVAL. The inverse feature allows the prediction of sub networks when the targets of perturbation are unknown.

Construction of Signalling Pathway Model

For the construction of a logic model to predict the influence of gemcitabine on apoptosis, pathway models were obtained from KEGG [21] and merged together in Cytoscape [22]. This model was further modified by pruning non-essential nodes, and the apoptosis model described by Mai and Liu [23] was integrated into the study. The network consists of extrinsic and intrinsic pathways for apoptosis to evaluate the process of signal transduction. Tumor necrosis factor (TNF) generally activates the extrinsic pathway by binding to the tumor necrosis factor receptors (TNFRs). This leads to a cascade of reaction, wherein TNF receptor-associated death domain protein (TRADD) and Fas-associated death domain protein (FADD) are activated. This results in the activation of caspase8 (CAS8), and subsequently caspase 3 (CAS3). The pro-survival signalling components of the extrinsic pathway inhibits the activity of CAS8 through TNF receptor-associated factors (TRAF) and the cellular form of FLICE-inhibitory protein (cFLIP). However, for the cause of simplicity, the cascade reaction from TNFR to CAS8 was represented by two nodes; one for activation through FADD, and another for inhibition through cFLIP. TNFR mediated nuclear factor kappa-light-chain-enhancer of activated B cells (NF-κB), and TNFR mediated interaction with TP53 was included in the model. Transduction of extrinsic pathway by epidermal growth factor (EGF) through the epidermal growth factor receptor (EGFR), and vascular endothelial growth factor A (VEGFA) through the vascular endothelial growth factor receptor (VEGFR) occurs through the phosphatidylinositol3-kinase (PI3K) pathway, which in turn activates the AKT protein. Activated AKT inhibits caspase 9 (CAS9), an essential component of the apoptosis machinery. Due to the indication that hypoxia plays an unknown, yet major role in the progression against gemcitabine [24-26], Hypoxia-inducible factor 1 (HIF1), a transcription factor that acts as the master regulator of hypoxia was included in the model. HIF1 overexpression plays a significant role in promoting tumor growth, along with angiogenesis [27, 28]. Thus, HIF1 modulated activation of VEGFA was modelled.

The intrinsic pathway consists of a pro-apoptotic Bcl-2-associated death promoter (BAD), Bcl-2-associated X protein (BAX), and anti-apoptotic BCL-2 protein family protein (BCLX). BAD protein receives inhibitory signals through the AKT pathway, further inhibiting its downstream effector, BCLX. BCLX prevents the mitochondria from processing pro-apoptotic signals by inhibiting its function. On the other hand, the BAX protein receives activating signals from tumor protein p53 (TP53). This signal is passed into the mitochondria, which then releases cytochrome C molecules in response to BAX protein, triggering the activation of APAF1, and subsequently, CAS9 and CAS3, leading to cell death. The cell can, however, halt the process of cell death by activating an apoptosis-inhibiting protein called X-linked inhibitor of apoptosis protein (XIAP). This protein is regulated by the nuclear factor kappa-light-chain-enhancer of activated B cells (NF-κB), a protein complex that influences cell survival. However, pro-apoptotic molecules released by the mitochondria can inhibit the activity of XIAIP. The overall setup for this study was simplified, which may vary with some experimental conditions. The model presented in this work simply illustrates the mechanistic signalling process rather than providing a comprehensive aspect into the detailed working process, which is beyond the scope of this study.

Tools and Parameters for Pathway Simulation

Analysis of Networks with Interactive Modeling (ANIMO) [29] is a tool for modelling and simulating signalling networks, where elements of the network are represented in the form of nodes and edges. The activity of each node is specified by the user based on experimental data. The user can also specify whether the nodes activate or inhibit the downstream node, and the rate at which each reaction occurs (very slow, slow, medium, fast, very fast). The activity level can be expressed between 0 and 100, where 0 signifies an inactive node, while 100 represents full activity. For this study, a simple model for the apoptosis pathway described above was used for simulation.

For this study, the protein activity levels used in the model were curated from the literature. For consistency, all data associated with protein activities were mined from experiments where PANC-1 cell lines were treated with 100nm of gemcitabine. Due to large numbers of unknown parameters, components with insufficient experimental data were set to 0. Transcription factors, on the other hand, were assigned parameters 0 or 1 to signify “off” or “on” state, respectively.

RESULTS AND DISCUSSION

Inferring Differential Gene Expression and CARNIVAL Results

A total of 1147 differentially expressed genes were obtained from DEGUST, where 905 genes were upregulated and 242 genes were downregulated. Using the list of differentially expressed genes, a casual network was generated using inverse CARNIVAL. The computed network is represented in Fig. ().

The network consisted of 39 transcription factors, of which several were actively involved in the progression of pancreatic cancer in one form or the other. Some notable transcriptional factors which were predicted to be overexpressed were nuclear factor NF-kappa-B p65 subunit (RELA), Tumor protein p53 (TP53), nuclear factor kappa-light-chain-enhancer of activated B cells (NF-κB), and Signal transducer and activator of transcription (STAT). TP53 is an important regulator of apoptosis, which also served as machinery for DNA repair. This transcription factor became active when the cell was subjected to stress. In this study, the expression of TP53 was positively correlated with the increasing dose of gemcitabine, indicating that gemcitabine triggered the activity of TP53, thus paving a path for cell death. NF-κB, on the other hand, acted as a master regulator for cell proliferation and survival, preventing cells from dying, and thus led to uncontrollable cell growth. RELA is a protein which is essential for forming the NF-κB complex, and is positively associated with several types of cancer. STATs are transcription factors found in the intercellular region and their dysregulation is often related to cancer survival and angiogenesis. Another important transcription factor that plays a critical role in gemcitabine resistance is the hypoxia-inducible factor 1-alpha (HIF1A). This transcription factor is known to induce the transcription of several essential genes involved in cell proliferation and survival. HIF1A is very often upregulated in many types of cancer, and this overexpression may be attributed to hypoxic conditions in the tumor microenvironment. Studies have suggested that hypoxia results in gemcitabine resistance in human pancreatic cell lines. Moreover, comparison of expression data between normal PANC-1 cell line, and cell lines treated with 100nm and 1µm gemcitabine indicated an increase in genes regulated by HIF1A. The pathway scores from PROGENy indicated an increase in the level of the hypoxic pathway with increasing gemcitabine dose. Thus, the results obtained from this study support the fact that hypoxia may play a major role in the resistance towards gemcitabine. Although some transcription factors promoted overexpression of genes, others were observed to negatively impact gene regulation. MYC belongs to a family of transcriptional factors that regulate cell proliferation, cell growth, and even apoptosis. While the presence of this protein regulates apoptosis positively, its negative regulation may contribute to anti-apoptotic effects.

Effect of Gemcitabine on the Apoptosis Signaling Pathway

To gain insight into the possible effects of gemcitabine on the apoptosis pathway, a simplified model consisting of components relevant to apoptosis signalling pathway was constructed based on literature review. According to the findings of this study and several other sources [24-26], hypoxia seems to play a central role in helping the tumor cell confer resistance towards gemcitabine. It can be speculated that the activity of components essential for blocking the process of cell death is severely upregulated under this condition. To decipher the underlying process, the constructed model was simulated for a period of 240 minutes. The resulting simulation provided a glimpse into the possible mechanism of cell death regulation in the presence and absence of hypoxia (Fig. ).

It was observed that under active hypoxia, the activity of VEGFA increased rapidly, causing downstream elements AKT and MTOR to remain active indefinitely. HIF1 regulates the production of VEGFA, allowing it to be secreted continuously under the influence of hypoxia [30]. High expression of AKT and MTOR proteins was observed in the PANC-1 cell line treated with 100nm gemcitabine [31]. The activation of AKT led to the induction of IKK, thus keeping the NF-ΚB complex in a perpetual “on” state. This led to the understanding that the activity of IKK prevented the oscillatory signalling of NF-κB. However, the exact mechanism behind this process was not deduced. Since NF-κB regulates the activity of XIAP, the latter became saturated, only slightly inhibited by a weakly expressed SMAC. These findings were found to be consistent with reports found in the literature [32]. On contrast, the activity of CAS8 was noticeably absent, indicating that the extrinsic TNF pathway inhibited the expression of CAS8. The activity of cleaved CAS3, which is the key component for cell death, was also found to be severely down-regulated due to the inhibitory action of XIAP [31, 33, 34]. The increased anti-apoptotic and decreased pro-apoptotic activity of proteins can be correlated with increased hypoxia in gemcitabine treated cells, essentially making them drug-resistant not only due to the absence of cell death, but also due to the induction of epithelial-to-mesenchymal transition (EMT) that helps in cell proliferation. Based on the literature findings, it was evident that hypoxia increased the expression of both HIF1 and NF-κB in PANC-1 cells, making them resistant to gemcitabine treatment [35].

Given the importance of NF-κB in promoting cell survival, and the fact that the expression of NF-κB is correlated with HIF1 due to hypoxia, a hypothesis was put forth to validate the fact that the inhibition of HIF1 could lead to the shutting down of the NF-κB transcription factor, thus effectively leading to cell death in the absence of XIAP. However, the signal transduction for NF-κB to be regulated could occur through the EGF/VEGF pathway, or the TNF pathway. To understand the mechanism behind this operation, the model was divided into two distinct sections in order to test the hypothesis described above. In the first section (Fig. ), only elements that were directly connected to the TNF pathway were included, including the intrinsic pathway. After simulating the model for 200 minutes, it was observed that the level of XIAP was drastically reduced, while the increase in CAS3 level was inversely correlated with that of XIAP. This was due to normal functioning of IKK in this signal transduction model that led to the self-inhibition of NF-κB as described in the canonical model [36].

In the alternate signaling model, EGF/VEGF model was considered as the primary source of signal transduction, where HIF1 influences the production of VEGFA. After simulating the model for 200 minutes, the level of AKT, IKK, XIAP, and MTOR was highly saturated, however CAS3 activity level was still low due to the strong inhibition by XIAP (Fig. and 4B). When hypoxia node was disabled/silenced and the simulation was visualized, normal signaling of apoptosis resumed along with normal regulation of MTOR (Fig. and 4D). This suggests that hypoxia not only regulated the level of apoptosis, but it leads to an increase in cell growth and proliferation by influencing the VEGF pathway.

CONCLUSION

In this study, using gene expression data, the important regulators of signalling pathways were highlighted. Although gemcitabine mechanism of action has been reported in several studies, the explicit effect that it has on other signalling pathways is not fully explored yet. With the application of systems biology, important perturbed elements were identified in this study which could help understand crosstalk between cellular signalling. Studies have reported that NF-κB may be activated under hypoxia, however, the underlying mechanism remains under speculation. In this work, the relation between hypoxia and apoptosis in gemcitabine treated cell lines has been implicated with the help of pathway simulation. Following a cascade of reaction, IKK activates NF-κB, which in turn induces the activation of XIAP, thus preventing the cell from attaining apoptosis. However, in the absence of hypoxia, NF-κB is switched off, leading to the inhibition of XIAP followed by normal cell death (Fig. ).

In conclusion, this study provides an insight into the effect of gemcitabine on gene regulation by identifying several potential components that partake in cellular signalling. These components were essential for the regulation of genes that drive PDAC. It was further concluded that hypoxia acts as a key driver in preventing tumor cells from achieving cell death through NF-κB regulation, and promotes cell growth and proliferation through MTOR. All these processes are thought to be regulated through the EGF/VEGF pathway. Thus, future experimental studies can be carried out in the wet lab to establish the findings presented in this study.