The complexity of cancer signaling networks limits the efficacy of most single-agent treatments and brings about challenges in identifying effective combinatorial therapies. In this study, we used chronic active B-cell receptor (BCR) signaling in diffuse large B-cell lymphoma as a model system to establish a computational framework to optimize combinatorial therapy in silico. We constructed a detailed kinetic model of the BCR signaling network, which captured the known complex cross-talk between the NFκB, ERK, and AKT pathways and multiple feedback loops. Combining this signaling model with a data-derived tumor growth model, we predicted viability responses of many single drug and drug combinations in agreement with experimental data. Under this framework, we exhaustively predicted and ranked the efficacy and synergism of all possible combinatorial inhibitions of eleven currently targetable kinases in the BCR signaling network. Ultimately, our work establishes a detailed kinetic model of the core BCR signaling network and provides the means to explore the large space of possible drug combinations. Cancer Res; 77(8); 1818–30. ©2017 AACR.

Major Findings

Using chronic active B-cell receptor signaling in diffuse large B-cell lymphoma as a model system, we developed a kinetic-modeling based computational framework for predicting effective combination therapy in silico. By integrative modeling of signal transduction, drug kinetics, and tumor growth, we were able to directly predict drug-induced cell viability responses at various dosages, which were in agreement with published cell line experimental data. We implemented computational screening methods that identified potent and synergistic combinations in silico and validated our independent predictions in vitro.

Quick Guide to Equations and Assumptions

A kinetic model was constructed for the BCR signaling network according to the following rules. For protein–protein binding interactions, we assumed that this type of reaction were under equilibrium and solved the steady-state level of protein complex analytically in the model,

where [TA] and [TB] stand for the total concentration of protein species A and B, respectively, [AB] represents the concentration of the protein complex, |K$| is the inverse of the dissociation constant Kd. For kinase catalyzed reactions, we adopted the classic Michaelis–Menten kinetics,

where |\frac{{dP}}{{dt}}$| is the rate of catalytic product formation, kcat is the enzyme turnover rate, Km is the Michaelis–Menten constant, [E] and [S] are concentration of enzyme and substrate, respectively. Each interaction in the BCR signaling network was depicted by corresponding terms according to above rules. The full kinetic model consists of 28 state variables, each representing concentration of a specific form of a protein species, depicted by 10 steady-state equations and 18 ordinary differential equations (see Supplementary Text S1).

A data-driven tumor growth model was formulated to link signaling response to cell viability output. We assume that tumor cell population is at exponential growth phase,

where N0 and N are cell number at time 0 and time t, |{r_d}$| is basal death rate, while growth rate |r$| is dependent on three downstream survival and proliferation signals NFκB, pAKT, and pERK (normalized by untreated control) through a Hill function,

Here, r* is basal growth rate, n is Hill coefficient, w1, w2, w3 are the estimated weights of three signaling outputs NFκB, pAKT, and pERK. Therefore, viability response defined as the ratio of cell number monitored under treated condition |N$| (for a time span of |T$|⁠) and untreated control |{N_c}$| can be predicted as following,

The activation of intracellular signaling pathways in response to environmental stimulus leads to important cell decisions such as proliferation. The amplitude and duration of pathway activation are precisely and robustly controlled by complex regulatory loops to maintain cellular homeostasis. In cancer, activating mutations or deletion of signaling repressors frequently result in sustained and exaggerated pathway activation that drives uncontrolled tumor survival and proliferation. Targeted therapies that use small-molecule inhibitors to repress specific signaling pathway members, for example, kinases, can directly block oncogenic pathway activation and lead to tumor cell death. These targeted therapies are expected to provide improved efficacy and reduced toxicity compared with chemotherapy.

However, clinical application of targeted therapy is facing several challenges such as low response rate and frequently acquired drug resistance. The limited efficacy of single-agent–targeted therapy is at least partially due to pathway cross-talks and compensatory circuits within signaling networks targeted by these agents (1). Cross-talks and compensatory circuits allow signals to bypass drug inhibition and reactivate downstream effectors. By simultaneously repressing multiple nodes in a signaling network, combination therapy has the potential to completely extinguish oncogenic signaling and induce more potent and durable treatment response. Thus, novel drug combinations where two or more drugs work cooperatively to suppress corrupted signaling networks need to be identified to achieve maximum therapeutic efficacy. The complexity of signaling networks makes it difficult to simply guess which combinations will be effective and synergistic and which ones will not. Moreover, given the large number of possible drug combinations against complex signaling networks, comprehensive experimental screening, including exploration of multiple dosages, is not practically feasible. Besides, results acquired from such screening may be specific to the cell line tested, thus lacking general applicability to highly variable primary tumors found in patients.

Computational models of signaling networks that can accurately reconstruct signaling dynamics in silico may represent a useful alternative to experimental screening and trial-and-error experimental investigation. Once proven reliable, these models can be used to exhaustively test the efficacy of a large number of single drug and drug combinations by quantifying signaling output under corresponding network perturbations. Even though computational modeling has been widely used to study the dynamics of signaling network in the past decades, the development of cancer signaling models and its application to predicting effective combinatorial therapies is still lacking. Here we demonstrate the feasibility of this approach using chronic activation of B-cell receptor (BCR) signaling in diffuse large B cell lymphoma (DLBCL) as a model system. We adopted a systems biology approach and established a computational framework to optimize anti-DLBCL combinatorial therapy in silico. The proposed approach is broadly applicable and can be used for other malignancies driven by aberrantly active signaling pathways.

The deregulation of BCR signaling is central to the pathogenesis of many B-cell malignancies. It is especially central in the activated B-cell–like subtype of DLBCL (ABC DLBCL). ABC DLBCLs exhibit chronic active BCR signaling, and are addicted to constitutive activation of downstream survival and proliferation signals such as NFκB (2). It has recently been found that a subset of the germinal center B-cell like subtype of DLBCLs are also dependent on BCR signaling through activation of the PI3K/AKT pathway (3). Multiple small-molecule inhibitors against BCR signaling were developed and proved effective in killing BCR-dependent DLBCLs in vitro and in vivo (4–6). However when tested in clinical trials, single-agent treatments again demonstrated limited responsiveness and efficacy (7), suggesting an urgent need for the design of effective combination therapies.

In this work, we present a kinetic modeling–based computational framework for predicting and optimizing combinatorial therapy against chronic active BCR signaling (Fig. 1). We constructed a detailed kinetic model of the BCR signaling network parameterized by published signaling responses and protein concentrations. Mathematical models of proximal BCR signaling and downstream transcriptional network have been reported (8–10). But, to our knowledge, this is the first kinetic model to reconstruct the entire core BCR signaling network in silico. Using published drug response data in a BCR signaling–dependent cell line, we trained a tumor growth model, which in combination with the kinetic model, allowed us to simulate viability response upon various targeted treatments. Under this framework, we exhaustively tested the efficacy and synergism of all possible combinations of inhibition of eleven currently targetable kinases in the BCR signaling network. We discuss how these results pave the way for the discovery of effective drug combinations.

Figure 1.

Outline of the approach taken in the current study. The central BCR signaling network was constructed on the basis of validated protein–protein interactions from experimental literature. Parameters of molecular reaction kinetics were estimated from phosphorylation time course data and protein concentrations were retrieved from MOPED protein expression database. A phenotypic tumor growth model was trained on cell viability assays of inhibitor treatments to link signaling response to viability outcome. In the end, simulation of the signaling model in combination with the tumor growth model was conducted to optimize treatment strategy. The model's prediction was compared with published drug response data and new prediction-driven hypotheses were tested independently in vitro.

Figure 1.

Outline of the approach taken in the current study. The central BCR signaling network was constructed on the basis of validated protein–protein interactions from experimental literature. Parameters of molecular reaction kinetics were estimated from phosphorylation time course data and protein concentrations were retrieved from MOPED protein expression database. A phenotypic tumor growth model was trained on cell viability assays of inhibitor treatments to link signaling response to viability outcome. In the end, simulation of the signaling model in combination with the tumor growth model was conducted to optimize treatment strategy. The model's prediction was compared with published drug response data and new prediction-driven hypotheses were tested independently in vitro.

Close modal

Protein–protein interactions in the BCR signaling network

Antigen-induced BCR crosslinking allows SRC family kinases, mainly LYN, to phosphorylate the immuno-receptor tyrosine–based activation motifs (ITAM) of the intracellular BCR subunits Igα (CD79A) and Igβ (CD79B; ref. 11). Dually phosphorylated ITAM motifs then recruit SYK and activate it via tyrosine phosphorylation (12). Activated SYK phosphorylates adapter BLNK, which recruits BTK to the plasma membrane to facilitate its phosphorylation and subsequent activation by SYK and LYN (13). Activated BTK further phosphorylates PLCγ2, which catalyzes the hydrolysis of phosphatidylinositol-4,5-bisphosphate (PI(4,5)P2) into diacylglycerol (DAG) and inositol trisphosphate (IP3; ref. 14). DAG together with elevated intracellular calcium induced by IP3 triggered endoplasmic reticulum (ER) calcium-release activates PKCβ (15), which then stimulates two divergent pathways that activate NFκB and ERK, respectively. Phosphorylation of CARD11 by PKCβ leads to the assembly of the CBM complex composed of CARD11, BCL10, and MALT1 (16). CBM acts as a scaffolding complex that facilitates IKK phosphorylation by TAK1, which in turn phosphorylates IKB and induces its degradation, releasing NFκB into the nucleus to elicit transcriptional activity (17). In addition, protease activity of MALT1 positively regulates NFκB signaling by cleaving and inactivating inhibitors against NFκB activation such as A20 and RELB (18, 19). In the meantime, PKCβ and DAG activate RASGRP, which triggers the canonical MAPK signaling cascade, leading to eventual phosphorylation and activation of ERK (20). On the other hand, SYK and LYN phosphorylate BCAP and CD19, respectively, which activate PI3K by membrane recruitment (21, 22). PI(3,4,5)P3 synthesized by PI3K further facilitates PDK1-catalyzed AKT phosphorylation by binding to both proteins via their plextrin homology (PH) domains (23). Importantly, LYN negatively regulates PI3K signaling by activating SHIP1, which hydrolyzes PI(3,4,5)P3 into PI(4,5)P2 (24).

Besides major signal transduction pathways as described above, our model includes key regulatory interactions in the BCR signaling network such as pathway cross-talks and feedback loops. The PI3K pathway positively regulates NFκB and ERK signaling by enhancing BTK membrane recruitment via PI(3,4,5)P3 binding. At the same time, it conversely attenuates ERK signaling via AKT catalyzed RAF phosphorylation (25). It has recently been found that MEK negatively regulates PI3K/AKT signaling by recruiting PTEN to the plasma membrane (26), which dephosphorylates PI(3,4,5)P3 into PI(4,5)P2. BTK amplifies BCR signaling by two coupled positive feedback loops. It recruits PIP5K to the plasma membrane, which produces PI(4,5)P2 to sustain both PI(3,4,5)P3 synthesis and PI(4,5)P2 hydrolysis (27). In addition, BTK phosphorylates BCAP, further facilitating the membrane recruitment of PI3K (22). The activity of BTK is attenuated by active PKCβ via disruption of its membrane localization, constituting a negative feedback loop (28). Besides, another indirect feedback from PKC to SYK was added into the model as knockdown of PKCδ was shown to mediate hyperphosphorylation of SYK (29). Furthermore, multiple negative feedback loops exist within the MAPK signaling cascade to fine-tune its activation amplitude and duration (30).

Cell viability assay

In the published drug response experiment (31), cells were treated at time 0 and incubated for 48 hours. Viability response was normalized to the plate positive control (bortezomib) and negative control(DMSO) as described previously (31). These normalized data were used for comparison with simulation results.

For our independent validation of synergistic and antagonistic drug combinations, the DLBCL cell lines TMD8, HBL1, and OCI-LY10 were generously donated from the Melnick laboratory at Weill Cornell Medicine (New York, NY) in 2014. These cell lines were verified using the cell line authentication service at Biosynthesis (BSI) using their STR Profiling and Comparison Analysis Service (32) in 2015. TMD8, HBL1 cell lines were cultured using RPMI1640 media supplanted with 10% FBS, 5% HEPES, 5% l-glutamine, and 5% penicillin/streptomycin. OCI-LY10 cell line was cultured using Iscove's Modified Dulbecco's Media supplanted with 10% FBS, 5% l-glutamine, and 5% PS. All cell culture media and reagents were from Life Technologies. R406, dasatinib, and MI-2 were purchased from Selleck Chemicals. Cells were grown at concentrations sufficient to keep untreated cells in exponential growth during the time of drug exposure. Cells were treated with 6 doses of each drug or combination in triplicate. Drug combinations were administered in constant ratio. Cell viability was determined by an ATP luminescent method (CellTiter-Glo, Promega). Luminescence was measured with the Syngery4 microplate reader (BioTek). Cell viability in drug-treated cells was normalized to vehicle-treated controls.

Kinetic model of BCR signaling network

Upon antigen-induced BCR crosslinking, most kinases involved in the BCR signaling pathway are recruited to the plasma membrane through activation of adaptor proteins such as BLNK, BCAP, and CD19 to form an integrated signalosome (33). Therefore, we expect most interactions depicted in our kinetic model to happen in a localized fashion. Thus, we do not consider any compartmentalization and model the cell as a well-mixed system.

As protein–protein binding is in general a very fast process at the second time scale, we assumed that this type of reaction were under equilibrium and solved the steady-state level analytically in the model. For a reversible protein–protein binding reaction,

at steady-state,

where [A] and [B] stand for the concentration of freed form of A and B; [TA] and [TB] stand for the total concentration of A and B; [AB] represents the concentration of the complex. By solving the above three equations, we have

where |K = \displaystyle{\frac{{{k_ + }}}{{{k_ - }}}$|⁠, is the inverse of the dissociation constant Kd.

Under a few circumstances where a protein may bind to more than one partner, the interactions were considered independently for simplification.

For kinase catalyzed reaction, we adopted the classic Michaelis–Menten kinetics,

where |\frac{{dP}}{{dt}}$| is the rate of catalytic product formation, kcat is the enzyme turnover rate, Km is the Michaelis–Menten constant, [E] and [S] are concentration of enzyme and substrate, respectively.

Reactions in the BCR signaling network were written into corresponding equations according to rules discussed above. The full model consists of 28 state variables each representing concentration of a specific form of a protein species, depicted by 10 steady-state equations and 18 ordinary differential equations (ODE; see Supplementary Text S1). Kinetic parameters in the model are summarized in Supplementary Table S1. Bounds for kinetic parameter estimation are listed in Supplementary Table S2. Parameters of total protein concentrations are summarized in Supplementary Table S3. We performed parameter sensitivity analysis where each parameter was perturbed independently across four orders of magnitudes and viability response was recorded. We found overall robustness and identified the most sensitive parameters as parameters regulating main axis of the NFκB pathway (see Supplementary Fig. S1).

Input signal

We imposed temporal pLYN and pSYK stimulus in normal BCR signaling modeled by two double exponential functions,

where |{T_{LYN}}$| and |{T_{SYK}}$| are total concentration of LYN and SYK, respectively. |{f_{LYN}}$| and |{f_{SYK}}$| are the phosphorylated fraction of LYN and SYK, respectively. |{f_{LYN}}$| and |{f_{SYK}}$| were estimated by fitting to the phosphorylation time course data using a genetic algorithm. A,B,C,D |{\tau _1},{\tau _2},{\tau _3},{\tau _4}$|were estimated by fitting to the normalized pLYN and pSYK time course.

In contrast, we imposed constitutive pLYN and pSYK stimulus in simulations of diseased ABC DLBCLs, with negative feedback from SYK and PKC, respectively.

The 0.2 coefficient is to account for LYN attenuation effect due to CD79B mutation in TMD8 (2). The parameter for negative feedback is estimated by fitting to single-drug viability response.

Tumor growth model

Assume a tumor cell population is at exponential growth phase (34, 35),

where N0 and N are cell number at time 0 and time t, |{r_d}$| is basal death rate, while growth rate |r$| is dependent on three downstream survival and proliferation signals NFκB, pAKT, and pERK (normalized by untreated control) through a Hill function,

Here r* is basal growth rate, n is Hill coefficient, w1, w2, w3 are the estimated weights of three signaling outputs NFκB, pAKT, and pERK. Therefore, viability response defined as the ratio of cell number monitored under treated condition |N$| (for a time span of |T$|⁠) and untreated control |{N_c}$| can be predicted as following,

Parameters required in this function were trained with viability data of three single-drug viability responses, namely NFκB, AKT, and MEK inhibitor, respectively. First the level of the three downstream survival and proliferation signals were predicted via simulation of the signaling model, and then input into the tumor growth functions to compute the viability output. Parameters were chosen by minimizing the sum of residuals between the viability predictions and experimental data.

Drug kinetics

To simulate an inhibitor's effect at a given dosage, percent activity of the targeted kinase was calculated via the medium effect equation,

The drug's IC50 was taken from literature (Supplementary Table S4), while m was assumed to be 1 under a first-order approximation. Then, the activity of the targeted kinase (i.e., parameters representing catalytic or activation rate of targeted kinase) was reduced to the corresponding percentage in the kinetic model. We list perturbed parameters in each simulated inhibitor treatment in Supplementary Table S4.

Synergy quantification

Under the Bliss Independence model, the additive effect of two inhibitors is computed as the multiplication of the effect of individual inhibitors,

where |{F_{UA}}$| indicates the fraction unaffected. To evaluate mode of interaction between two inhibitors, we computed viability response at 10×10 virtual dosages by varying the percent inhibition of each targeted kinase independently from 0% to 90% at 10% interval. These viability values were used to estimate parameter that minimizes the following metric,

where x,y are virtual dosages for inhibitor 1 and 2, respectively. |\beta \langle {1,\beta = 1,\beta } \rangle 1$| indicates synergism, additive, and antagonism, respectively.

Computational software

The numerical simulations of ODEs were performed in Matlab _R2012a while data analysis was performed in R. Specifically, we used ode15s function in Matlab for solving ODEs and ga function in Matlab for parameter estimation using genetic algorithms. We have uploaded scripts for ODE simulation and parameter estimation in github (https://github.com/weiduweillcornell/BCR-signaling-model).

Kinetic modeling of BCR signaling network reproduces normal BCR signaling in silico

We first curated the central BCR signaling network by gathering experimentally validated protein–protein interactions from literature. The reconstructed network is shown in Fig. 2, and includes three major signaling pathways downstream of BCR, namely NFκB, PI3K/AKT, and RAF/RAS/ERK. We chose to include these three pathways because they have been reported to closely regulate cell survival and proliferation in B cells and B-cell malignancies (36). Detailed information on the interactions included in the constructed BCR signaling network is provided in Materials and Methods.

Figure 2.

The central BCR signaling network constructed from literature. Antigen binding induces BCR aggregation and subsequent phosphorylation, which further triggers a complex signaling cascade initiated by phosphorylated LYN and SYK. The BTK–PLCγ2–PKCβ pathway activates downstream NFκB and ERK through divergent paths, while membrane recruitment of PI3K leads to AKT activation. Pathway cross-talks and feedback regulations are highly abundant in the network.

Figure 2.

The central BCR signaling network constructed from literature. Antigen binding induces BCR aggregation and subsequent phosphorylation, which further triggers a complex signaling cascade initiated by phosphorylated LYN and SYK. The BTK–PLCγ2–PKCβ pathway activates downstream NFκB and ERK through divergent paths, while membrane recruitment of PI3K leads to AKT activation. Pathway cross-talks and feedback regulations are highly abundant in the network.

Close modal

Instead of directly applying mass action kinetics to characterize elementary reactions in the network, we chose to adopt more streamlined mathematical representations derived from mass action law under reasonable assumptions (see Materials and Methods). This strategy greatly reduced the number of variables, equations, and most importantly parameters required in the mathematical model. As elementary protein–protein binding reactions generally reach equilibrium within seconds, we modeled them by deriving the steady-state relationships from mass action law (see Materials and Methods). For enzymatic reactions such as phosphorylation or dephosphorylation, we adopted a classic Michaelis–Menten kinetic framework. To parameterize the model, we first retrieved protein concentrations in B lymphocytes quantified by mass spectrometry from the MOPED protein expression database (37). We modeled LYN and SYK as two independent input signals that triggered a downstream response. Their activation kinetics were approximated by double exponential functions (see Materials and Methods) where parameters in the functions were estimated by fitting to the phosphorylation time course data (38). We then used genetic algorithms to optimize the remaining 72 kinetic parameters within bounded biologically reasonable ranges (Supplementary Table S2; refs. 39–41) by minimizing residual sum of squares between simulated phosphorylation time courses and published Western blot data (38). Experimental data and simulated results were each normalized to their respective maximum value for comparison. Ten sets of kinetic parameters were identified from 5,000 independent runs that fit almost equally well (Supplementary Table S5). Simulated trajectory under these 10 parameter sets together with phosphorylation time course data are shown in Fig. 3A. To examine whether our fitting method is prone to overfitting problem, we manually injected different levels of white noise into the training time course data and refitted the model using our genetic algorithm. The simulation was done using optimized parameter values as initial population of parameter choice and 10 independent optimization processes were performed for each noise level. We showed that the residual sum of square greatly increases as the noise level increases (Supplementary Fig. S2). This indicates that the model is not prone to fitting to noise. Parameter sensitivity analysis was performed as described in Materials and Methods (Supplementary Fig. S1). We note that predicted activation time course of BLNK, PLC, and PKC deviate from experimental data. We believe that the most likely reason behind this relatively poor fit may be currently unknown regulatory interactions involving additional proteins in the BCR signaling pathway. However the relatively poor fit in these components in the upstream of NFκB pathway would not significantly influence the model's predictions as our parameter sensitivity analysis indicates viability response is most sensitive to parameters in the downstream part of the NFκB pathway (Supplementary Fig. S1). Of note, we were independently able to find 39 kinetic parameter values from the literature, and we compared these values with the range of estimated 10 parameter sets (Fig. 3B). By shuffling the literature-retrieved parameters 10,000 times, we found that the literature-retrieved parameter values fall within the estimated parameter ranges significantly more often than random (P = 0.05; Fig. 3C). We note, however, that many discrepancies were found between estimated and published parameters (Fig. 3B). We speculate that many of these discrepancies are likely due to in vitro nature of the experiments used to quantify kinetic parameters.

Figure 3.

Simulation of normal BCR signaling and estimation of kinetic parameters. A, Simulated trajectory of ten parameter sets in comparison with published phosphorylation time course data. B, Comparison between literature-retrieved parameter values with simulation-estimated parameter ranges. Box plot indicates the simulation-estimated parameter ranges, while black dots represent literature-retrieved parameter values. C, Empirical cumulative distribution of the number of parameters that would fall within simulation-estimated parameter ranges by chance derived by shuffling the literature-retrieved parameters 10,000 times.

Figure 3.

Simulation of normal BCR signaling and estimation of kinetic parameters. A, Simulated trajectory of ten parameter sets in comparison with published phosphorylation time course data. B, Comparison between literature-retrieved parameter values with simulation-estimated parameter ranges. Box plot indicates the simulation-estimated parameter ranges, while black dots represent literature-retrieved parameter values. C, Empirical cumulative distribution of the number of parameters that would fall within simulation-estimated parameter ranges by chance derived by shuffling the literature-retrieved parameters 10,000 times.

Close modal

Combining BCR signaling model with a tumor growth model predicted cell viability response upon single and combinatorial drug treatments in a BCR signaling–dependent ABC DLBCL cell line

We next sought to simulate the effect of various small-molecule inhibitors on ABC DLBCL cell viability and to compare simulation results with published drug response data in a BCR signaling–dependent ABC DLBCL cell line TMD8 (31). We selected TMD8 because of the extensive drug combinatorial data available on this cell line (31). We first made several modifications to the model to accommodate the differences between normal BCR signaling and aberrant BCR signaling in ABC DLBCL. Instead of applying a temporal upstream stimulus, we assumed constitutive LYN and SYK phosphorylation as observed both in ABC DLBCL cell lines and in primary DLBCL patient samples (see Materials and Methods; refs. 2, 42). In addition, we accounted for genetic alterations in members of the BCR signaling network in TMD8 compared with normal B cells. Specifically, TMD8 was shown to carry CD79B mutation that attenuates LYN activity by approximately 80% (2). Correspondingly, we decreased the enzymatic activity of LYN in the model to the same extent (see Materials and Methods).

To predict cell viability response from signaling output, we formulated a tumor growth model in which the growth rate of tumor cells is dependent on the weighted sum of the three downstream survival and proliferation signals NFκB, ERK, and AKT through a Hill function (see Materials and Methods). A similar formalism has been used to characterize tumor growth of ERBB-amplified breast cancer driven by ERK and AKT activation (43). We used published viability response data in TMD8 to parameterize the tumor growth model, where cells were treated with IKK, AKT, and MEK inhibitors at multiple dosages (31). Specifically, using the median effect equation (44), we calibrated the percent activity left on the targeted kinase for each inhibitor at a given dosage based on the inhibitor's IC50 value (see Materials and Methods; Supplementary Table S4). We then reduced the activity of the targeted kinase to the same level in the model and simulated steady-state signaling output. Parameters in the tumor growth model were estimated by minimizing residual sum of squares between predicted viability response and experimental data.

We first simulated single-drug viability response of inhibitors covering the NFκB, PI3K/AKT, and MAPK pathway and compared with experimental data. We observed that in silico simulation with the BCR signaling model and the tumor growth model recapitulated the viability response of the three training single-drug response, namely IKK, AKT, MEK inhibitors (Fig. 4A; ref. 31). This is not surprising as the growth model was fitted on the basis of training data. As independent predictions, we also simulated drug response of inhibitors targeting other kinases in the network, for example, CAL-101 against PI3K, ibrutinib against BTK, and found that predicted viability response matched favorably with TMD8 drug response data (Fig. 4B; ref. 31). At the same time, we found simulated viability response of SYK inhibition to deviate from experimental data (gray line), yet this discrepancy can be partially rescued by adding a negative feedback from SYK to LYN (blue line). It has been reported that SYK functions as a negative regulator of BCR signaling by phosphorylating Ig-α (29, 45, 46). As Ig-α primarily interacts with LYN, we assumed in the model that SYK indirectly negatively regulates LYN.

Figure 4.

Training and prediction of single drug viability response in ABC DLBCL cell line TMD8. A, Tumor growth model parameterization using single drug viability response of inhibitors targeting NFκB, AKT, and MEK. Gray dashed lines correspond to simulation results of model without SYK to LYN-negative feedback, while brown dashed lines correspond to simulation results of model without PI3K-NFκB cross-talk. B, Single drug viability response of inhibitor targeting various kinases against BCR signaling network.

Figure 4.

Training and prediction of single drug viability response in ABC DLBCL cell line TMD8. A, Tumor growth model parameterization using single drug viability response of inhibitors targeting NFκB, AKT, and MEK. Gray dashed lines correspond to simulation results of model without SYK to LYN-negative feedback, while brown dashed lines correspond to simulation results of model without PI3K-NFκB cross-talk. B, Single drug viability response of inhibitor targeting various kinases against BCR signaling network.

Close modal

Beyond single-drug viability response, we also simulated combinatorial drug response of ibrutinib in combination with various other kinases targeting the BCR pathway, and observed the predicted response contour to match favorably with experimental results (Fig. 5). These results demonstrate that our model can correctly capture the interaction between inhibitors as well.

Figure 5.

Combinatorial drug viability responses of BTK inhibitor ibrutinib in combination with additional inhibitors targeting BCR network intermediates were predicted and compared with experimental data.

Figure 5.

Combinatorial drug viability responses of BTK inhibitor ibrutinib in combination with additional inhibitors targeting BCR network intermediates were predicted and compared with experimental data.

Close modal

Overall, these results suggest that viability response of small-molecule inhibitors targeting the BCR signaling network can be predicted via in silico simulation of the BCR signaling model in combination with the tumor growth model.

Cross-talk between PI3K and NFκB pathway mediates efficacy of PI3K inhibition in TMD8

In both the drug response data and model's simulation, we observed that PI3K inhibition is significantly more effective at inhibiting tumor growth than blockage of its downstream effector AKT. A similar phenomenon was reported in other studies, where PI3K inhibition was shown to attenuate NFκB transcriptional activity (3, 47). We hypothesized that the efficacy of PI3K inhibition is primarily attributed to suppression of NFκB signaling, which is mediated by upstream cross-talk between the PI3K and NFκB pathways. To test this hypothesis, we abolished the cross-talk between PI3K and NFκB by knocking out in silico PI(3,4,5)P3-mediated membrane recruitment of BTK in the signaling model. Under this condition, we resimulated the viability response of PI3K inhibition, which showed significant reduction compared with both experimental data and simulation with the full signaling model (Fig. 4, brown line). This result supports the notion that the upstream crosstalk between PI3K and NFκB pathway is critical in mediating tumor growth inhibition by PI3K inhibitor. It also provides further rational support for the clinical use of PI3K inhibitors in DLBCL that are dependent in NFκB signaling (3, 47).

Computational optimization of targeted therapy against chronic active BCR signaling

Using the above modeling framework, we sought to identify targeted therapies against the BCR signaling network that most effectively inhibit tumor growth. We exhaustively tested all drug pairs based on 11 small-molecule inhibitors currently available that target various kinases in the network, yielding 55 treatment strategies in total. In each scenario, viability response was simulated at 10×10 virtual dosages where each targeted kinase was inhibited at 0% to 99% evenly spaced in log10 space. We calculated area under the combinatorial viability response surface as an overall indicator of drug combination potency. The smaller the value is, the more potent the drug target combination is (Fig. 6A). We found that under the same inhibition potency, efficacy of different treatment strategies was highly variable, ranging from almost no growth inhibition to up to 80% reduction (Fig. 6B). Specifically, inhibiting downstream of the NFκB signaling pathway, especially through MALT1 and IKK inhibitor, exhibited the most prominent efficacy, and combined MALT1 and IKK blockage yielded highest tumor growth inhibition. In comparison, tumor cell growth was relatively insensitive to blockage of MAPK pathway in our simulations. In summary, this computational screening result suggests that various treatment strategies against a signaling network can yield highly variable therapeutic responses and that in silico simulation can help identify targets that confer intrinsic vulnerability.

Figure 6.

A, Computational optimization of treatment strategy against chronic active BCR signaling. Viability response surface of three drug target combinations. Two horizontal axes correspond to virtual dosage of two different inhibitors, while the vertical axis indicates cell viability normalized to untreated condition. For each drug target combination, 10 virtual dosages (between 0% to 99% inhibition evenly spaced in log10 space) of each single drug are tested. B, Bar plot of simulated viability response of all possible dual inhibition on 11 kinases in the BCR signaling network that are currently targetable. Here viability responses are calculated area under the combinatorial viability response surface (as shown in Fig. 6A) as an overall indicator of drug combination potency. Binary codes on the bottom indicate the treatments applied (black, targeted inhibition).

Figure 6.

A, Computational optimization of treatment strategy against chronic active BCR signaling. Viability response surface of three drug target combinations. Two horizontal axes correspond to virtual dosage of two different inhibitors, while the vertical axis indicates cell viability normalized to untreated condition. For each drug target combination, 10 virtual dosages (between 0% to 99% inhibition evenly spaced in log10 space) of each single drug are tested. B, Bar plot of simulated viability response of all possible dual inhibition on 11 kinases in the BCR signaling network that are currently targetable. Here viability responses are calculated area under the combinatorial viability response surface (as shown in Fig. 6A) as an overall indicator of drug combination potency. Binary codes on the bottom indicate the treatments applied (black, targeted inhibition).

Close modal

We then sought to identify drug combinations that are synergistic via computational simulations. For a given two-drug combination, the combinatorial drug response at 10× 10 virtual dosage as discussed above were used to estimate mode of drug interaction under the Bliss independence model (see Materials and Methods; Fig. 7A). Computational screening predicted dual blockage of LYN and SYK as the most synergistic combination. To test this prediction, we treated ABC DLBCL cell lines TMD8, HBL1, and OCI-LY10 with LYN inhibitor dasatinib and SYK inhibitor R406, at multiple doses. Comparing combinatorial drug response data to theoretical additive response predicted by the Bliss independence model (see Materials and Methods), we confirmed synergism between dasatinib and R406 (Fig. 7B). We also tested and confirmed the predicted antagonism between dual SYK and MALT1 inhibition using SYK inhibitor R406 and MALT1 inhibitor MI-2 across three ABC DLBCL cell lines (Fig. 7B).

Figure 7.

A, Modes of interaction of all pairwise inhibitions under Bliss Independence model. β < 1, β = 1, β > 1 correspond to synergism, additive, and antagonism, respectively. B,In vitro validation of predicted synergistic and antagonistic drug combination in ABC DLBCL cell line TMD8, HBL1, and OCI-LY10. R406, dasatinib, and MI-2 are inhibitors against SYK, LYN, and MALT1, respectively.

Figure 7.

A, Modes of interaction of all pairwise inhibitions under Bliss Independence model. β < 1, β = 1, β > 1 correspond to synergism, additive, and antagonism, respectively. B,In vitro validation of predicted synergistic and antagonistic drug combination in ABC DLBCL cell line TMD8, HBL1, and OCI-LY10. R406, dasatinib, and MI-2 are inhibitors against SYK, LYN, and MALT1, respectively.

Close modal

It is increasingly acknowledged that aberrant BCR signaling plays a central role in the development and maintenance of many B-cell malignancies (48). Although a large panel of small-molecule inhibitors against BCR signaling have been developed, rational methodologies that can predict effective combinatorial therapy and guide the design of specific treatment strategy in individual patients have been lacking. To bridge this gap, we aimed to construct the first kinetic model of the core BCR signaling network and use this model to investigate targeted therapy against BCR signaling. We showed that simulations with the signaling model reconstructed dynamics of normal B-cell signaling in silico. Combining the signaling model with a data-trained tumor growth model successfully predicted viability response of multiple drug combinations, and identified novel synergistic drug combinations such as LYN and SYK inhibitor.

In this work, we chose to develop a highly detailed kinetic model of the BCR signaling network, modeling direct protein–protein interactions using highly quantitative ODEs wherever possible. We note that there are other simpler modeling techniques for analysis of signaling pathways, for example, simplifying signaling cascades as coupled Hill functions (43) or using Boolean network models. To investigate the predictive power of a simpler model, we constructed a Boolean model of a simplified BCR signaling network consisting of only the 11 targetable nodes together with signaling pathway end points and a node representing cell viability (Supplementary Fig. S3). In this model, each node has only two states (0, inactive; 1, active), and the state value in a particular step is determined by the values of all its regulators in the previous step. A node will be active if the majority of its regulators are activating, except that the Viability node is active when any of the signaling outputs is active. Assuming the two input nodes LYN and SYK are constitutively active, all possible initial states (2^12) are exhaustively simulated until reaching attractors. These initial states ended up in two attractors, both attractors consisting of states in which Viability is on (Supplementary Fig. S4). Then all two drug combinations were tested. To simulate drug-mediated inhibition, each targeted node was made constitutively inactive. There are four drug combinations that result in a single global attractor in which Viability is off. These four drug combinations are BTK-MEK, BTK-RAF, SYK-MEK, and SYK-RAF. These four drug combinations ranked low in the drug efficacy predictions made by the full ODE model (Fig. 6B). Moreover, experimental data in TMD8 cell line (31) indicates that the BTK-MEK and BTK-RAF combinations are not very effective in decreasing cell viability—the MEK and RAF inhibitors are not responsive (along the y-axis) and synergism is lacking (Supplementary Fig. S5). Overall, these results suggest that a simple Boolean network model is not able to capture the same results as the full ODE model and that its prediction accuracy is lacking when compared with experimental data.

The clinical application of targeted therapy is frequently challenged by highly variable drug response among cancer patients. Heterogeneous response to BCR signaling–targeted therapy was observed in ABC DLBCL cell lines (47) and in clinical trials of ABC DLBCL patients (7). This is likely due to differential expression of proteins in the BCR signaling pathway that impact pathway activities in individual patient. In this study, signal transduction is explicitly modeled on the molecular level, which provides a straightforward framework to incorporate protein level variation to develop patient-specific predictive models. Specifically, expression levels of proteins within the BCR signaling network can be measured experimentally using protein expression profiling techniques such as reverse-phase protein array (RPPA) or other approaches. Then protein expression changes in patients relative to cell line can be incorporated into the model to predict optimal treatment strategy for individual patients. We believe that the use of patient-specific predictive models can greatly improve the performance of targeted therapy in cancer.

In this study, we focused on the exploration of combinatorial efficacy of kinase inhibitors against BCR signaling pathway components. In the future, it may be compelling to combine kinase inhibitors with traditional chemotherapeutic drugs in the clinical setting. Regrettably, capturing the combinatorial effect of kinase inhibitors and chemotherapeutic drugs likely requires a complex theoretical framework that may involve a more detailed mechanistic characterization of cell cycle and DNA damage response than the simple growth model we use here.

Finally, we note that besides DLBCL, aberrant BCR signaling was shown to play a role in other B-cell malignancies such as chronic lymphocytic leukemia (CLL; ref. 49) and mantle cell lymphoma (MCL; ref. 50). In phase II studies of BTK inhibitor ibrutinib, 71% and 68% overall response rate (ORR) was reported in CLL and MCL patients, respectively (51), suggesting targeting BCR signaling as promising treatment strategy. Correspondingly, the overall framework and predictions reported in this work may also be useful to identify drug combinations for CLL and MCL-targeted treatment.

L. Cerchietti reports receiving a commercial research grant from Celgene. A.M. Melnick is a consultant at Epizyme and reports receiving a commercial research grant from Janssen, GSK, and Eli Lilly. No potential conflicts of interest were disclosed by the other authors.

Conception and design: W. Du, A.M. Melnick, O. Elemento

Development of methodology: W. Du, O. Aly

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): R. Goldstein, Y. Jiang, O. Aly

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): W. Du, R. Goldstein, L. Cerchietti, A.M. Melnick

Writing, review, and/or revision of the manuscript: W. Du, O. Aly, A.M. Melnick, O. Elemento

Study supervision: L. Cerchietti, O. Elemento

We thank all Elemento and Melnick laboratory members for productive discussions.

This work was supported by the CAREER grant from National Science Foundation (DB1054964), NIH grant R01CA194547, the Starr Cancer Consortium and Hirschl Trust. A.M. Melnick is supported by NCI R01 CA143032, R01 CA104348, the Burroughs Welcome Foundation, and the Chemotherapy Foundation.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

1.
Du
W
,
Elemento
O
. 
Cancer systems biology: embracing complexity to develop better anticancer therapeutic strategies
.
Oncogene
2015
;
34
:
3215
25
.
2.
Davis
RE
,
Ngo
VN
,
Lenz
G
,
Tolar
P
,
Young
RM
,
Romesser
PB
, et al
Chronic active B-cell-receptor signalling in diffuse large B-cell lymphoma
.
Nature
2010
;
463
:
88
92
.
3.
Chen
L
,
Monti
S
,
Juszczynski
P
,
Ouyang
J
,
Chapuy
B
,
Neuberg
D
, et al
SYK inhibition modulates distinct PI3K/AKT- dependent survival pathways and cholesterol biosynthesis in diffuse large B cell lymphomas
.
Cancer Cell
2013
;
23
:
826
38
.
4.
Advani
RH
,
Buggy
JJ
,
Sharman
JP
,
Smith
SM
,
Boyd
TE
,
Grant
B
, et al
Bruton tyrosine kinase inhibitor ibrutinib (PCI-32765) has significant activity in patients with relapsed/refractory B-cell malignancies
.
J Clin Oncol
2013
;
31
:
88
94
.
5.
Fontan
L
,
Yang
C
,
Kabaleeswaran
V
,
Volpon
L
,
Osborne
MJ
,
Beltran
E
, et al
MALT1 small molecule inhibitors specifically suppress ABC-DLBCL in vitro and in vivo
.
Cancer Cell
2012
;
22
:
812
24
.
6.
Naylor
TL
,
Tang
H
,
Ratsch
BA
,
Enns
A
,
Loo
A
,
Chen
L
, et al
Protein kinase C inhibitor sotrastaurin selectively inhibits the growth of CD79 mutant diffuse large B-cell lymphomas
.
Cancer Res
2011
;
71
:
2643
53
.
7.
Roschewski
M
,
Staudt
LM
,
Wilson
WH
. 
Diffuse large B-cell lymphoma-treatment approaches in the molecular era
.
Nat Rev Clin Oncol
2014
;
11
:
12
23
.
8.
Barua
D
,
Hlavacek
WS
,
Lipniacki
T
. 
A computational model for early events in B cell antigen receptor signaling: analysis of the roles of Lyn and Fyn
.
J Immunol
2012
;
189
:
646
58
.
9.
Ravichandran
S
,
Rao
KV
,
Jain
S
. 
Bistability in a model of early B cell receptor activation and its role in tonic signaling and system tunability
.
Mol Biosyst
2013
;
9
:
2498
511
.
10.
Martinez
MR
,
Corradin
A
,
Klein
U
,
Alvarez
MJ
,
Toffolo
GM
,
di Camillo
B
, et al
Quantitative modeling of the terminal differentiation of B cells and mechanisms of lymphomagenesis
.
Proc Natl Acad Sci U S A
2012
;
109
:
2672
7
.
11.
Gauld
SB
,
Cambier
JC
. 
Src-family kinases in B-cell development and signaling
.
Oncogene
2004
;
23
:
8001
6
.
12.
Rolli
V
,
Gallwitz
M
,
Wossning
T
,
Flemming
A
,
Schamel
WW
,
Zurn
C
, et al
Amplification of B cell antigen receptor signaling by a Syk/ITAM positive feedback loop
.
Mol Cell
2002
;
10
:
1057
69
.
13.
Baba
Y
,
Hashimoto
S
,
Matsushita
M
,
Watanabe
D
,
Kishimoto
T
,
Kurosaki
T
, et al
BLNK mediates Syk-dependent Btk activation
.
Proc Natl Acad Sci U S A
2001
;
98
:
2582
6
.
14.
Kim
YJ
,
Sekiya
F
,
Poulin
B
,
Bae
YS
,
Rhee
SG
. 
Mechanism of B-cell receptor-induced phosphorylation and activation of phospholipase C-gamma2
.
Mol Cell Biol
2004
;
24
:
9986
99
.
15.
Spitaler
M
,
Cantrell
DA
. 
Protein kinase C and beyond
.
Nat Immunol
2004
;
5
:
785
90
.
16.
Shinohara
H
,
Yasuda
T
,
Aiba
Y
,
Sanjo
H
,
Hamadate
M
,
Watarai
H
, et al
PKC beta regulates BCR-mediated IKK activation by facilitating the interaction between TAK1 and CARMA1
.
J Exp Med
2005
;
202
:
1423
31
.
17.
Wegener
E
,
Krappmann
D
. 
CARD-Bcl10-Malt1 signalosomes: missing link to NF-kappaB
.
Sci STKE
2007
;
2007
:
pe21
.
18.
Coornaert
B
,
Baens
M
,
Heyninck
K
,
Bekaert
T
,
Haegman
M
,
Staal
J
, et al
T cell antigen receptor stimulation induces MALT1 paracaspase-mediated cleavage of the NF-kappaB inhibitor A20
.
Nat Immunol
2008
;
9
:
263
71
.
19.
Hailfinger
S
,
Nogai
H
,
Pelzer
C
,
Jaworski
M
,
Cabalzar
K
,
Charton
JE
, et al
Malt1-dependent RelB cleavage promotes canonical NF-kappaB activation in lymphocytes and lymphoma cell lines
.
Proc Natl Acad Sci U S A
2011
;
108
:
14596
601
.
20.
Coughlin
JJ
,
Stang
SL
,
Dower
NA
,
Stone
JC
. 
RasGRP1 and RasGRP3 regulate B cell proliferation by facilitating B cell receptor-Ras signaling
.
J Immunol
2005
;
175
:
7179
84
.
21.
Okada
T
,
Maeda
A
,
Iwamatsu
A
,
Gotoh
K
,
Kurosaki
T
. 
BCAP: the tyrosine kinase substrate that connects B cell receptor to phosphoinositide 3-kinase activation
.
Immunity
2000
;
13
:
817
27
.
22.
Fujimoto
M
,
Fujimoto
Y
,
Poe
JC
,
Jansen
PJ
,
Lowell
CA
,
DeFranco
AL
, et al
CD19 regulates Src family protein tyrosine kinase activation in B lymphocytes through processive amplification
.
Immunity
2000
;
13
:
47
57
.
23.
Ding
Z
,
Liang
J
,
Li
J
,
Lu
Y
,
Ariyaratna
V
,
Lu
Z
, et al
Physical association of PDK1 with AKT1 is sufficient for pathway activation independent of membrane localization and phosphatidylinositol 3 kinase
.
PLoS ONE
2010
;
5
:
e9910
.
24.
O'Neill
SK
,
Getahun
A
,
Gauld
SB
,
Merrell
KT
,
Tamir
I
,
Smith
MJ
, et al
Monophosphorylation of CD79a and CD79b ITAM motifs initiates a SHIP-1 phosphatase-mediated inhibitory signaling cascade required for B cell anergy
.
Immunity
2011
;
35
:
746
56
.
25.
Zimmermann
S
,
Moelling
K
. 
Phosphorylation and regulation of Raf by Akt (protein kinase B)
.
Science
1999
;
286
:
1741
4
.
26.
Zmajkovicova
K
,
Jesenberger
V
,
Catalanotti
F
,
Baumgartner
C
,
Reyes
G
,
Baccarini
M
. 
MEK1 is required for PTEN membrane recruitment, AKT regulation, and the maintenance of peripheral tolerance
.
Mol Cell
2013
;
50
:
43
55
.
27.
Saito
K
,
Tolias
KF
,
Saci
A
,
Koon
HB
,
Humphries
LA
,
Scharenberg
A
, et al
BTK regulates PtdIns-4,5-P2 synthesis: importance for calcium signaling and PI3K activity
.
Immunity
2003
;
19
:
669
78
.
28.
Kang
SW
,
Wahl
MI
,
Chu
J
,
Kitaura
J
,
Kawakami
Y
,
Kato
RM
, et al
PKCbeta modulates antigen receptor signaling via regulation of Btk membrane localization
.
EMBO J
2001
;
20
:
5692
702
.
29.
Limnander
A
,
Zikherman
J
,
Lau
T
,
Leitges
M
,
Weiss
A
,
Roose
JP
. 
Protein kinase Cdelta promotes transitional B cell-negative selection and limits proximal B cell receptor signaling to enforce tolerance
.
Mol Cell Biol
2014
;
34
:
1474
85
.
30.
Reth
M
,
Brummer
T
. 
Feedback regulation of lymphocyte signalling
.
Nat Rev Immunol
2004
;
4
:
269
77
.
31.
Mathews Griner
LA
,
Guha
R
,
Shinn
P
,
Young
RM
,
Keller
JM
,
Liu
D
, et al
High-throughput combinatorial screening identifies drugs that cooperate with ibrutinib to kill activated B-cell-like diffuse large B-cell lymphoma cells
.
Proc Natl Acad Sci U S A
2014
;
111
:
2349
54
.
32.
Cabrera
CM
,
Cobo
F
,
Nieto
A
,
Cortes
JL
,
Montes
RM
,
Catalina
P
, et al
Identity tests: determination of cell line cross-contamination
.
Cytotechnology
2006
;
51
:
45
50
.
33.
Dal Porto
JM
,
Gauld
SB
,
Merrell
KT
,
Mills
D
,
Pugh-Bernard
AE
,
Cambier
J
. 
B cell antigen receptor signaling 101
.
Mol Immunol
2004
;
41
:
599
613
.
34.
Assanga
I
,
Lujan
L
. 
Cell growth curves for different cell lines and their relationship with biological activities
.
Int J Biotechnol Mol Biol Res
2013
;
4
:
60
70
.
35.
Rodriguez-Brenes
IA
,
Komarova
NL
,
Wodarz
D
. 
Tumor growth dynamics: insights into evolutionary processes
.
Trends Ecol Evol
2013
;
28
:
597
604
.
36.
Niiro
H
,
Clark
EA
. 
Regulation of B-cell fate by antigen-receptor signals
.
Nat Rev Immunol
2002
;
2
:
945
56
.
37.
Kolker
E
,
Higdon
R
,
Haynes
W
,
Welch
D
,
Broomall
W
,
Lancet
D
, et al
MOPED: model organism protein expression database
.
Nucleic Acids Res
2012
;
40
:
D1093
9
.
38.
Kumar
D
,
Srikanth
R
,
Ahlfors
H
,
Lahesmaa
R
,
Rao
KV
. 
Capturing cell-fate decisions from the molecular signatures of a receptor-dependent signaling response
.
Mol Syst Biol
2007
;
3
:
150
.
39.
Chen
WW
,
Niepel
M
,
Sorger
PK
. 
Classic and contemporary approaches to modeling biochemical reactions
.
Genes Dev
2010
;
24
:
1861
75
.
40.
Stroppolo
ME
,
Falconi
M
,
Caccuri
AM
,
Desideri
A
. 
Superefficient enzymes
.
Cell Mol Life Sci
2001
;
58
:
1451
60
.
41.
Peters
GH
,
Branner
S
,
Moller
KB
,
Andersen
JN
,
Moller
NP
. 
Enzyme kinetic characterization of protein tyrosine phosphatases
.
Biochimie
2003
;
85
:
527
34
.
42.
Bogusz
AM
,
Baxter
RH
,
Currie
T
,
Sinha
P
,
Sohani
AR
,
Kutok
JL
, et al
Quantitative immunofluorescence reveals the signature of active B-cell receptor signaling in diffuse large B-cell lymphoma
.
Clin Cancer Res
2012
;
18
:
6122
35
.
43.
Kirouac
DC
,
Du
JY
,
Lahdenranta
J
,
Overland
R
,
Yarar
D
,
Paragas
V
, et al
Computational modeling of ERBB2-amplified breast cancer identifies combined ErbB2/3 blockade as superior to the combination of MEK and AKT inhibitors
.
Sci Signal
2013
;
6
:
ra68
.
44.
Chou
TC
,
Talalay
P
. 
Quantitative analysis of dose-effect relationships: the combined effects of multiple drugs or enzyme inhibitors
.
Adv Enzyme Regul
1984
;
22
:
27
55
.
45.
Sidorenko
SP
,
Law
CL
,
Klaus
SJ
,
Chandran
KA
,
Takata
M
,
Kurosaki
T
, et al
Protein kinase C mu (PKC mu) associates with the B cell antigen receptor complex and regulates lymphocyte signaling
.
Immunity
1996
;
5
:
353
63
.
46.
Heizmann
B
,
Reth
M
,
Infantino
S
. 
Syk is a dual-specificity kinase that self-regulates the signal output from the B-cell antigen receptor
.
Proc Natl Acad Sci U S A
2010
;
107
:
18563
8
.
47.
Kloo
B
,
Nagel
D
,
Pfeifer
M
,
Grau
M
,
Duwel
M
,
Vincendeau
M
, et al
Critical role of PI3K signaling for NF-kappaB-dependent survival in a subset of activated B-cell-like diffuse large B-cell lymphoma cells
.
Proc Natl Acad Sci U S A
2011
;
108
:
272
7
.
48.
Young
RM
,
Staudt
LM
. 
Targeting pathological B cell receptor signalling in lymphoid malignancies
.
Nat Rev Drug Discov
2013
;
12
:
229
43
.
49.
Burger
JA
,
Chiorazzi
N
. 
B cell receptor signaling in chronic lymphocytic leukemia
.
Trends Immunol
2013
;
34
:
592
601
.
50.
Boukhiar
MA
,
Roger
C
,
Tran
J
,
Gressin
R
,
Martin
A
,
Ajchenbaum-Cymbalista
F
, et al
Targeting early B-cell receptor signaling induces apoptosis in leukemic mantle cell lymphoma
.
Exp Hematol Oncol
2013
;
2
:
4
.
51.
Hendriks
RW
,
Yuvaraj
S
,
Kil
LP
. 
Targeting Bruton's tyrosine kinase in B cell malignancies
.
Nat Rev Cancer
2014
;
14
:
219
32
.