Tyrosine kinase inhibitors (TKIs) constitute a class of cancer therapeutics, many of which are known to cause cardiotoxicity as a major adverse event. Reported cardiotoxicities include heart failure, cardiomyopathy, conduction abnormalities, QT prolongation, and myocardial injury. The most common toxicity is systolic dysfunction or cardiomyopathy, potentially leading to heart failure, which is most likely mediated through direct toxicity of cardiomyocytes (Albini et al., 2010; Eschenhagen et al., 2011; Force and Kolaja, 2011; Raschi and De Ponti, 2012; Ewer and Ewer, 2015).

Trastuzumab, an inhibitor of the HER2 receptor tyrosine kinase (Slamon et al., 2001; Piccart-Gebhart et al., 2005; Romond et al., 2005; Force et al., 2007) was both the first monoclonal antibody TKI given FDA approval (in 1998) and the first TKI reported to cause cardiotoxicity (Wu et al., 2016). Since the reports of trastuzumab-induced toxicity, several additional targeted cancer therapeutics have been classified as cardiotoxic, observations that have contributed to the emergence of a new research field, cardio-oncology (Albini et al., 2010; Bellinger et al., 2015).

Previous studies have shown that TKI-related cardiotoxicity, as seen with trastuzumab, is mostly due to the targeting of pathways that are shared between malignancies and cardiovascular cells (De Keulenaer et al., 2010; Bellinger et al., 2015). Investigations of these adverse events revealed that many of the tyrosine kinases targeted by TKIs serve critical roles in survival and maintenance of cardiomyocytes, leading to unintended on-target toxicity. At the same time, many TKIs inhibit multiple kinases simultaneously, which can cause off-target toxicity (Chen et al., 2008; Force and Kerkelä, 2008; Force and Kolaja, 2011).

Despite the risk of cardiotoxicity, TKIs are still one of the highly effective and favored cancer therapeutics on the market (Eschenhagen et al., 2011; Force and Kolaja, 2011). The success of drugs such as trastuzumab and imatinib, a small molecule inhibitor used to treat chronic myeloid leukemia (CML), has inspired the development of additional TKIs. As of April, 2015, 25 small molecule TKIs have entered the market (Shah and Morganroth, 2015), with many more under development (Bellinger et al., 2015). Given the booming research in the development of TKIs, it would be beneficial to develop a systematic strategy to: (1) evaluate and predict how new TKIs will affect signaling networks in cardiomyocytes; and (2) identify interventions that can reverse and/or mitigate any associated cardiotoxicity. These questions are well-suited to be addressed using a quantitative systems pharmacology (QSP) approach that combines large-scale measurements with mechanism-based mathematical modeling. The diversity of TKI targets and the complexity of cellular mechanisms responsible for cardiotoxicity mean that two drugs with similar targets may operate through different mechanisms, and the effects of two TKIs with different targets may converge on a common pathway. Untangling this type of complexity generally requires computational approaches that are based on biological mechanisms. Therefore, our aims in this Perspective are to review the progress that systems approaches have made in predicting TKI-induced cardiotoxicity and to offer suggestions for how mathematical modeling can be applied to elucidate mechanisms and predict potential adverse events caused by new drugs.

The canonical roles of tyrosine kinases are found in mitogenesis and related processes such as differentiation, metabolism, and migration. Constitutive activation of tyrosine kinase (TK) signaling, via either gain-of-function (GOF) mutations or overexpression due to gene amplification, is found in about 70% of malignancies (Blume-Jensen and Hunter, 2001; Chen et al., 2008). Well-understood examples include overexpression of ERBB2 in HER2⁺ breast cancer (Force et al., 2007) and the constitutively active oncogenic fusion protein BCR-ABL, which can cause CML (Force et al., 2007; Chen et al., 2008; Force and Kolaja, 2011). This dependency of tumor formation and proliferation on TK signaling led to the rise of TKIs as promising anti-cancer therapeutics.

Currently, there are two chemical classes of TKIs: (1) humanized monoclonal antibodies (mAbs) and (2) small molecule inhibitors (Force et al., 2007; Chen et al., 2008; Force and Kolaja, 2011). Small molecule TKIs can be further subcategorized based on whether they compete with ATP for the binding pocket or interact with other regions of the protein (Force and Kolaja, 2011). Additionally, TKIs are often identified by the intended target(s) or the target specificity (Force et al., 2007; Bellinger et al., 2015; Gharwan and Groninger, 2015). The most common target groups that are used to classify TKIs include EGFR/ERBB2 inhibitors, VEGFR inhibitors, ABL inhibitors, and multi-targeted drugs that are designed to inhibit at least two different target groups such as VEGFR and ABL (see Table 1 for descriptions of the important cellular signaling proteins mentioned in the manuscript). Figure 1A shows the currently-approved TKIs, grouped by the published targets, and indicates how these classifications frequently overlap.

The initial discovery of TKI-induced cardiotoxicity was made during the groundbreaking clinical trials of trastuzumab, the first such drug to be marketed (Seidman et al., 2002; De Keulenaer et al., 2010). However, estimated cardiotoxicities of 3–7% with trastuzumab alone and 25% when the drug was administered with an anthracycline (Slamon et al., 2001) were only determined during a post hoc analysis. Similar retrospective analyses have been performed to estimate that sunitinib causes left ventricular dysfunction with an incidence of 4–11% (Yeh and Bickford, 2009; Lenneman and Sawyer, 2016) and the VEGFR inhibitor bevacizumab induces either cardiomyopathy or heart failure in 1.5–3% of patients (Yeh and Bickford, 2009). These examples demonstrate the difficulties associated with identifying cardiotoxicity during drug development. Because clinical trials are primarily focused on evaluating efficacy, they often lack appropriate safety screening measures to identify side effects (Force et al., 2007).

Overall, of the 30 TKIs currently marketed for use in the United States, 26 list serious cardiac side effects as a “black box warning” in their prescription information (FDA and CDER, 2012; Boehringer Ingelheim International GmbH, 2014; Gharwan and Groninger, 2015). The cardiac related black box warnings of TKIs can be categorized into: cardiomyopathy, arrhythmia, myocardial infarction, hypertension, and pericardial effusion, based on the specific potential adverse events listed in the package insert. In Figure 1B, which indicates both the adverse events and the target class for each TKI, we observe no obvious association between the intended primary target, and the reported cardiac risks.

The initial discovery of TKI-induced cardiotoxicity was met with surprise due to the fact that cardiomyocytes are non-dividing and terminally differentiated (Force et al., 2007). Since TKs were mostly known for their role in proliferation and their association with cancer, these kinases were not expected to have any essential role in cardiomyocytes, and toxicity in heart was not anticipated (Chen et al., 2008; Bellinger et al., 2015). The discovery of TKI-induced cardiotoxicity, therefore, became a driving force for uncovering the roles of tyrosine kinases in heart. The research spurred by these adverse events has allowed us to appreciate that many of the pathways responsible for proliferation in malignant cells also play important roles in cardiomyocytes in: (1) survival signaling; (2) mitochondrial and sarcoplasmic reticulum (SR) homeostasis; and (3) electrical and contractile function.

From survival and homeostasis to contractile function, tyrosine kinases perform a wide variety of important roles in the health and function of cardiomyocytes. Although considerable progress has been made to decipher the roles of individual TKs, the breadth of the different mechanisms involved makes it difficult to draw general conclusions about TKI-induced cardiotoxicity. Moreover, even when some mechanistic details have been uncovered, our understanding is primarily qualitative, and biological mechanisms cannot usually be connected to factors such as dosing and the physiological characteristics of individual patients. Given the past success of TKIs as cancer therapeutics and the drive to develop additional TKIs, it would be beneficial to develop a systematic strategy to: (1) evaluate the potential cardiotoxicity of new TKIs; (2) predict the mostly likely mechanisms involved; and (3) suggest strategies for mitigating and/or reversing toxicity. Systems approaches that perform large-scale measurements and quantitatively compare responses to multiple drugs are likely to be extremely useful for understanding the common and distinct features of cardiotoxicity caused by diverse TKIs.

A common approach in systems-level pharmacology studies is to utilize a cell based, high-throughput drug screening assay. Although cardiovascular pharmacology has traditionally not been well-suited for high throughput studies, the development of hiPSC-CMs has expanded the possibilities. For instance, a recent study described the development of a comprehensive assay that evaluated cellular effects of TKIs in hiPSC-CMs (Sharma et al., 2017). Using hiPSC lines from 13 individuals, the investigators examined how 21 FDA-approved small molecule TKIs affected cell viability, contractility, and gene expression. By integrating the results with literature-reported TKI serum levels in patients, the authors developed a novel cardiac safety index for TKIs (Sharma et al., 2017).

Although this study represents a significant milestone in that it integrates cutting edge technologies such as deep sequencing and high-throughput imaging, room for improvement nonetheless remains. Specifically, the experiments performed in this study represent snapshots of cell state after TKI treatment, and the data can therefore provide only limited insight into the dynamics of toxicity development. Moreover, experiments performed using individual drugs cannot predict how either a second drug or a circulating hormone (e.g., adrenaline, angiotensin) might either exacerbate or protect against cardiotoxicity. Although it is of course possible to expand the assay to treat cells with drug pairs and/or add relevant physiological stimuli, it is not clear a priori which additional perturbations might be informative or relevant. Finally, even when unambiguous results are seen in cellular high-throughput assays, the mechanistic details often remain hidden. It would therefore be helpful to couple such high-throughput measurements with integrative computational analyses that can potentially overcome these limitations.

One way to address the aforementioned limitations is to use mathematical models that mechanistically describe biological dynamics. When the processes simulated by these models overlap with toxicity mechanisms, the simulations can be used to generate testable predictions, to guide experimental studies, and ultimately to make decisions about new drugs based on a quantitative understanding of benefits and risks. In this context, models that describe biological mechanisms through differential equations are frequently referred to as QSP models (Leil and Bertz, 2014; Gadkar et al., 2016). Although precise definitions remain a matter of debate, QSP models are generally distinguished from both purely empirical, statistical approaches such as computing a risk score for a drug based on a series of measurements (Kramer et al., 2013; Mistry et al., 2015), and pharmacokinetic models that can predict the effects of dosing on cardiotoxicity (van Hasselt et al., 2012) but generally offer only limited mechanistic insight. Although QSP models have been exploited to understand cardiotoxicity caused by anthracyclines (de Oliveira et al., 2016), the application of QSP to TKI-induced cardiotoxicity is still in its early stages. Given the recent development of QSP modeling, it is instructive to consider examples in which mechanistic models have been successfully applied to the study of adverse events, in particular drug-induced liver injury (DILI) (Huang et al., 2013; Shoda et al., 2014; Yang et al., 2015), and drug-induced arrhythmias (Moreno et al., 2011; Sarkar and Sobie, 2011; Britton et al., 2013; Cummins et al., 2014; Grandi and Maleckar, 2016; Yang et al., 2016). Specifically, we emphasize how mechanistic models can be integrated with large in vitro data sets, as these studies may provide an important blueprint for future research on cardiotoxicity caused by kinase inhibitors.

An example of the success of QSP models for toxicity applications can be found in the development of DILIsym®, a mathematical model and software package used for predicting DILI. DILIsym® comprises multiple sub-models describing relevant biological processes involved in hepatotoxicity such as drug distribution in the liver, bile acid homeostasis, reactive metabolite generation and disposition, oxidative stress, immune responses, and the hepatocyte life cycle (Woodhead et al., 2017). The value of this approach was recently demonstrated in studies that examined differences in hepatotoxicity between acetaminophen and its less toxic isomer 3′-hydroxyacetnilide. Although the former drug can cause toxicity across many species, the latter has been shown to cause DILI in humans and rats, but not in mice. Using the mechanistic DILIsym® model, a testable hypothesis was generated in which the amount of reactive metabolite produced from each isomer was identified as an important contributor to the observed species differences, and this prediction was confirmed experimentally in a later study (Kyriakides et al., 2016). In addition to the consortium that has developed the DILIsym® package, other groups have gained important insight into DILI through mathematical modeling (Smith et al., 2016; Blais et al., 2017; Thiel et al., 2017).

To understand and predict drug-induced arrhythmias, e.g., Torsades de Pointes (TdP), the Comprehensive in vitro Proarrhythmia Assay (CiPA) initiative is highly relevant. This effort, part of the FDA's Critical Path Initiative, aims to improve the accuracy and cost effectiveness of screening for TdP risk. Whereas current in vitro methods for predicting TdP risk focus almost exclusively on block of K_v11.1 (i.e., the hERG channel), an approach that is often inadequate, CiPA intends to both assess how drugs block multiple ion channels and to combine these measurements with recordings in hiPSC-CMs and mechanistic simulations (Sager et al., 2014; Fermini et al., 2016). A couple of recently-published studies highlight the value that is gained by utilizing QSP models. Lancaster and Sobie, for example, used models of human ventricular myocytes to simulate physiological changes caused by 67 unique drugs, some that are known to cause TdP, and others that are apparently safe. In addition to providing a classification model that was superior to K_v11.1 block alone, the simulations provided testable predictions about the most informative assays to perform in cellular experiments and the specific ion transport pathways that, when affected by a drug, may contribute to TdP risk (Lancaster and Sobie, 2016). More recently, Li et al. showed that incorporating the kinetics of K_v11.1 block into simulations provides superior identification of TdP risk than simply considering steady-state block measurements (i.e., an IC₅₀-value), and the study suggested an experimental protocol for measuring drug block kinetics (Li et al., 2017).

The examples of both DILI and drug-induced TdP demonstrate the value that can be added when mechanistic modeling is used to address toxicity. The simulations can uncover the reasons for counterintuitive results, such as drugs that block K_v11.1 but are nonetheless safe (Kramer et al., 2013; Lancaster and Sobie, 2016; Li et al., 2017), or drugs that only cause hepatotoxicity in some species (Kyriakides et al., 2016; Smith et al., 2016; Blais et al., 2017; Thiel et al., 2017; Woodhead et al., 2017). The simulations can also suggest the prioritization of experiments that are most likely to provide additional insight.

TKI-induced cardiotoxicity is a problem that seems well-suited to a QSP approach because tyrosine kinase signaling encompasses large, complex networks with numerous feedback loops, and understanding how a drug alters TK cascades is therefore extremely complicated. Although mechanistic modeling to predict TKI-induced cardiotoxicity is much less well-developed than for DILI or TdP, efforts that may yield breakthroughs in the next few years are currently underway.

As noted above, an important recent study by Sharma et al. (2017) derived a “toxicity index” by examining effects of TKIs in hiPSC-CMs through several assays. Although the large-scale nature of this study justifies the “systems” label, and the toxicity index is a quantitative risk score, it's important to emphasize that QSP modeling, which can provide mechanistic insight and actionable predictions, is complementary to a high-throughput, data-driven strategy such as that used in this study (Sharma et al., 2017).

Another notable recent effort is a study by Shin et al. (2014). These investigators combined experimental measurements with simulations to uncover mechanisms by which high and sustained doses of the β-adrenergic receptor agonist isoproterenol could increase myocyte susceptibility to apoptosis (Shin et al., 2014). This study is an important example of how simulations often generate novel experimentally-testable predictions, and the work can potentially be extended to examine TKI-induced toxicity.

To fill current gaps in knowledge and obtain new mechanistic insight, the Drug Toxicity Signature Generation (DToxS) Center at Mount Sinai has initiated a large-scale project to advance cellular assays and computational approaches that can improve our understanding of TKI-induced cardiotoxicity. One aspect of this project involves employing QSP models that describe biological processes potentially involved in this cardiotoxicity. We outline here an approach by which whole transcriptome expression assays can be integrated with mechanistic models to classify drugs and generate novel, experimentally-testable predictions.

For this part of the analysis, the project has designed a standard experimental protocol (DtoxS—Drug Toxicity Signature Generation Center—SOP, https://martip03.u.hpc.mssm.edu/sop.php) that captures early effects of drugs, as reflected in gene expression changes. The experiments treat cultured cardiomyocyte-derived cells with potentially cardiotoxic TKIs as well as drugs from different classes that are presumably safe. After 48 h, mRNA is harvested, and sequencing is performed to quantify drug-induced changes in gene expression (compared with vehicle-treated controls). Data are released at the DToxS website and can be freely-downloaded (DtoxS—Drug Toxicity Signature Generation Center—Data & Resources, https://martip03.u.hpc.mssm.edu/data.php). The pipeline for integrating these released data with mechanistic mathematical models is shown in Figure 2, top. Changes in mRNA levels in drug treated cells (Figure 2A) can be translated into parameter alterations in models that describe processes potentially relevant to the toxicity (Figure 2B), and simulations are then performed with these models. This workflow assumes that before overt toxicity is induced, drugs can alter the cellular state in relatively subtle ways. Mechanistic simulations may then allow one to predict how this drug-induced altered cellular state influences the response to various physiological stimuli.

For example, results shown here are obtained with a QSP model that describes signaling events relevant to cardiac hypertrophy through a system of 106 ordinary differential equations, each one describing activity of a signaling component (Ryall et al., 2012). This model was chosen for initial simulations because the progression to heart failure due to pathological remodeling often includes a transient induction of hypertrophy. In addition, many known TKI targets (e.g., ERBB2, Raf-1) and critical nodes in cardiac survival signaling (e.g., PI3K, Akt, and ERK) are included. Using this model, simulations were performed to predict how drug-induced network alterations affected 7 hypertrophy biomarkers (the model's “outputs”), under conditions meant to simulate a variety of physiological or pharmacological stimuli (e.g., stretch, angiotensin, EGF, phenylephrine). For instance, time courses in Figure 2C show simulated normalized levels of Brain Natriuretic Peptide (BNP) under six conditions: before and after isoproterenol stimulation, and in three groups of cells: untreated (control), nilotinib-treated, and dasatinib-treated. BNP is an appropriate output to consider because it is both measured in patients with hypertrophy and heart failure and has been shown to be relevant for drug-induced cardiotoxicity (Nousiainen et al., 2002; Sandri et al., 2005; Skovgaard et al., 2014).

These example simulation results, summarized in Figure 2D, suggest that nilotinib leads to increases in BNP, both before and after isoproterenol, whereas dasatinib may reduce the upregulation of BNP that isoproterenol normally causes. The hypertrophy index, a summary statistic, condenses results by summing drug-induced changes across seven biomarkers in response to different stimuli applied in the model. The hypertrophy index confirms the impression from the simulated time courses, namely a pro-hypertrophic response to nilotinib contrasted with a slight anti-hypertrophic response to dasatinib (Figure 2E).

By simulating responses to stimuli, using gene expression changes induced by all drugs, patterns begin to emerge. Specifically, in myocytes that are not exposed to a physiological stimulus, TKIs and non-TKIs cause similar changes in hypertrophy biomarkers (Figure 2F). However, simulations predict that when TKI-treated myocytes are also exposed to isoproterenol or endothelin-1, agonists that are used experimentally to induce hypertrophy (Ichikawa et al., 1996; Yamazaki et al., 1996; Shohet et al., 2004; Ryall et al., 2012), the pro-hypertrophic response is exaggerated compared with non-TKI-treated cells.

These preliminary simulation results indicate the potential strengths of combining large scale measurements with mechanism-based mathematical models. First the simulations do not simply describe existing data—they can predict how drug-treated cells will respond to an additional stimulus that has not yet been applied experimentally. These predictions can be subsequently tested. Second, the simulation approach does not merely generate qualitative predictions; because the quantitative models predict that some drugs and/or stimuli may cause large effects whereas others cause only minor effects, this provides a means to prioritize experimental tests and use resources efficiently. Third, when clear differences are observed, for instance between individual drugs or between drug classes, the simulations predict the mechanisms responsible for the differences.

Although the preliminary simulation results shown here are encouraging, they also hint at the future research that must be performed to fully realize the potential of this approach. First, although the simulations predict how drug-induced changes in gene expression may influence both baseline signaling and cellular responses to stimuli, they do not describe direct inhibition of kinase activity by drugs, which is of course the more traditional and straightforward method for simulating drug effects. We excluded these effects from initial simulations because many TKI targets are not included in the model, but future work will expand the model by systematically adding drug targets based on published protein-protein interaction databases (Warde-Farley et al., 2010) and large-scale kinase inhibition assays (Anastassiadis et al., 2011; Davis et al., 2011). For such work, a promising way to expand the model will be to use large-scale, data-driven network identification algorithms (Thiagarajan et al., 2017) that can provide an unbiased approach for identifying potential off-targets.

A second important extension of the work will be to simulate additional biological processes potentially involved in cardiotoxicity. For instance, cell death via apoptosis, which may be important in the development of toxicity caused by some TKIs, has been described mathematically by many previous models (Schleich and Lavrik, 2013; Shin et al., 2014). Once these are tuned to reflect apoptotic signaling in cardiac myocytes, the models can be integrated with the experimental gene expression data to generate novel predictions. Similarly, models describing mitochondrial function, including the production of reactive oxygen species (Aon and Cortassa, 2012; Bazil et al., 2016; Wacquier et al., 2016), and the coupling of electrical excitation, and contractile function (Rice et al., 2008; Tewari et al., 2016), are also likely to be relevant. Finally, once a number of QSP models, describing additional processes, have been added, further secondary analyses can be performed. These include sensitivity analysis to identify the most important nodes in each model (Sobie, 2009), systematic simulations to potential targets for toxicity mitigation or reversal, and effects of combination therapy.

Here, we have discussed contemporary challenges in understanding TKI-induced cardiotoxicity and have illustrated how a QSP approach can be used to address unresolved questions and improve understanding. Previous successes of QSP in illuminating and predicting other forms of drug toxicity, including hepatotoxicity and drug-induced arrhythmia, demonstrate its potential utility for other drug toxicities. The initial results presented here show how mechanistic models can be integrated with “omics” measurements such as mRNA-seq, generating simulations that can suggest underlying mechanisms and help in prioritizing costly experiments. In the coming years, future work along these lines can be used to develop strategies to mitigate or reverse TKI-induced cardiotoxicity, thereby contributing to the development of therapeutic regimens that are both effective and safe.

JVS and BC performed an extensive literature review and analysis of published clinical and experimental data. JVS performed simulations and analysis. JVS and ES wrote the initial draft of the manuscript. BC, JvH, MB, and JJS provided critical feedback and editorial suggestions.