## Abstract

Combining radiotherapy with immune checkpoint blockade may offer considerable therapeutic impact if the immunosuppressive nature of the tumor microenvironment (TME) can be relieved. In this study, we used mathematical models, which can illustrate the potential synergism between immune checkpoint inhibitors and radiotherapy. A discrete-time pharmacodynamic model of the combination of radiotherapy with inhibitors of the PD1–PDL1 axis and/or the CTLA4 pathway is described. This mathematical framework describes how a growing tumor first elicits and then inhibits an antitumor immune response. This antitumor immune response is described by a primary and a secondary (or memory) response. The primary immune response appears first and is inhibited by the PD1–PDL1 axis, whereas the secondary immune response happens next and is inhibited by the CTLA4 pathway. The effects of irradiation are described by a modified version of the linear-quadratic model. This modeling offers an explanation for the reported biphasic relationship between the size of a tumor and its immunogenicity, as measured by the abscopal effect (an off-target immune response). Furthermore, it explains why discontinuing immunotherapy may result in either tumor recurrence or a durably sustained response. Finally, it describes how synchronizing immunotherapy and radiotherapy can produce synergies. The ability of the model to forecast pharmacodynamic endpoints was validated retrospectively by checking that it could describe data from experimental studies, which investigated the combination of radiotherapy with immune checkpoint inhibitors. In summary, a model such as this could be further used as a simulation tool to facilitate decision making about optimal scheduling of immunotherapy with radiotherapy and perhaps other types of anticancer therapies. *Cancer Res; 76(17); 4931–40. ©2016 AACR*.

To help the design and efficacy analysis of combined anticancer therapies, this article proposes a set of mathematical equations that describe the pharmacodynamics of radiotherapy in combination with two paradigmatic immunotherapies, namely the blockers of the PD1–PDL1 axis and of the CTLA4 pathway. These equations explain several experimental results reported in preclinical and clinical settings. Together with published pharmacokinetic models, they pave the way for the efficient *in silico* design and optimization of combined anticancer therapies mixing immunotherapy and radiotherapy; this approach is currently under active investigation because of its strong potential synergism.

## Introduction

Immune checkpoint inhibitors are considered as a major advance in a variety of cancer diseases with dismal prognosis (1). Despite outstanding progress in melanoma and lung cancers, a majority of patients will nevertheless fail in achieving long-lasting responses. Associating radiotherapy with immunotherapy has shown synergistic effect and is now considered as a promising strategy to stretch both response and survival at bedside (2). However, combining radiotherapy and immunotherapy can be done following a variety of different modalities. Determining the exact scheduling of the radioimmunotherapy association is therefore critical, and standard empirical approaches (trial-and-error) will hardly meet the current requirements of fast-track approvals of novel therapies. In this respect, developing modeling and simulation tools able to describe biological and pharmacodynamics processes should help to better determine an optimal solution among numerous possibilities.

Radiotherapy is widely used in the treatment of primary and metastatic tumors. The biological effects of radiation on tumors include DNA damage, modulation of signal transduction and alteration of the inflammatory tumor microenvironment (TME). The induction of irreparable DNA damage and cell-cycle arrest can induce lethal mechanisms, including mitotic catastrophe, apoptosis, necrosis, autophagy, and senescence (3). Traditionally, the DNA damages and the subsequent tumor cell kill have been explained with the classical 4 R's of radiobiology: repopulation, reoxygenation, repair of sublethal damage, and redistribution within the cell cycle (4). These mechanisms have long been described by widely used radiotherapy models such as the linear-quadratic (5–7).

Cancer immunotherapy aims at stimulating the ability of the host immune system to eliminate tumor cells through the recruitment and activation of cytotoxic effector cells. A close relationship exists between the tumor cells, the innate and the adaptive immune cells (8); this interplay may be explained by the theory of immunosurveillance or the 3 E's: elimination, equilibrium, and tumor escape. The cancer cells are first eliminated by an innate immune response involving natural killer cells, then by a specific immune response provoked by antigen-specific cytotoxic T cells. But the strength of the immune response gradually decreases during the equilibrium phase, as a result of several distinct and nonexclusive mechanisms, which eventually help the tumor to escape the immune system. Even though these complex biological mechanisms could not be known at the time, the gradual loss of tumor immunogenicity had been established in the 1970s by pioneering studies on the link between the size of a tumor and its ability to induce an immune response (9–11).

This escape mechanism could be, at least partly, reversible: several preclinical and clinical studies have reported the ability of radiation to activate the immune system, to enhance the immune response and to contribute to the control of cancer (12–15). The CD8^{+} T lymphocytes play a major role in antitumor immunity, hence radiotherapy requires their presence for postradiation tumor control (16) and their activity is stimulated by irradiation (17). Wattenberg and colleagues (18) reported on several cases where radiation sensitizes the tumor cells to immune-mediated killing, through a process called immunomodulation (18, 19). Importantly, they showed that the increase in the strength of the immune response was dose related.

Furthermore, the effect of radiation on the immune response may not be only one way: radiotherapy may also induce immunosuppressive effects by increasing the tumor cell expression of PD-L1, the secretion of TGFβ, the expansion of CD4^{+} T regulatory (Treg) cells, of myeloid-derived suppressor cells (MDSC), the maturation of type M2 macrophages (20, 21), and the death of the immune effectors located in the irradiation field (22, 23). Hence, the effects of irradiation could involve a balance between immunostimulation and immunosuppression.

Enhancing the antitumor immune response has long been a field of intensive research. As of today, the most promising drugs having reached clinical approval belong to two types of immune checkpoint inhibitors. They consist in two families of mAbs: the first one (such as ipilimumab) targets the cytotoxic T-lymphocyte–associated protein 4 (CTLA-4); the second one (such as nivolumab, pembrolizumab, and lambrolizumab) targets the programmed cell death 1 (anti-PD1) or its ligand (anti-PDL1). Preclinical studies have provided evidence of the synergistic effect between radiation, anti–CTLA-4 mAbs and/or anti PD-1 or PD-L1 mAbs in different *in vivo* models (12, 22, 24).

The abscopal effect of radiotherapy is defined as a regression of tumor or metastases in a site outside the field of irradiation, occurring at such a distance that the effect cannot be explained by direct or scattered radiation; it is immune-mediated and probably caused by the release of danger signals that activate CD8+ cytotoxic T cells, which participate in the control of the nonirradiated residual tumor and metastases (15). The abscopal effect of radiation in combination with ipilumumab was reported in different case reports (14, 25–27). These data from preclinical and clinical studies have confirmed the concept that radiotherapy could lead to the conversion of the tumor into an *in situ* vaccine.

This diversity of medical options paves the way to multimodal anticancer therapies that could unleash synergies between the different modalities. However, the determination of the most beneficial combination of immune therapies with radiation, in terms of doses, scheduling and synchronism for each modality, is crucial to the induction and persistence of a local and systemic antitumor immunity. To find such an optimal combination, a trial-and-error process would require too much effort and would most probably yield to the suboptimal selection of the best performing strategy among a very small set of candidates.

Herein, a mathematical model is proposed to help in the search of such a combined treatment. It contains a simplified model of tumor growth coupled with antitumor immune response. These equations correctly describe several experimental results about the relationship between the immunogenicity of a tumor and its mass dynamics. In particular, this framework explains how the tumor immunogenicity initially increases with its mass, but later vanishes as the TME becomes progressively more immunosuppressive. Then, by deriving a modified version of the linear-quadratic model, the effect of irradiation on the antitumor immune response are simply modeled by the release of tumor antigens out of irradiated tumor cells. Finally, pharmacodynamic equations are derived to describe the effects of the two families of immune checkpoints inhibitors (targeting either the PD1–PD-L1 axis or the CTLA4 pathway). These equations explain the experimental synergies between irradiation and immune checkpoint inhibitors that have been reported by several authors. Because this model is capable of describing a variety of experimental results, it could probably be used as a tool to simulate multimodal anticancer therapies and also help in the selection of the most efficient ones.

As a general modeling principle, no attempt was made to include all the currently available knowledge on immunology and radiotherapy. Trying to build a biologically exhaustive model would have resulted in an unnecessary complex mathematical framework of reduced practical interest. Instead, a limited number of biological hypotheses have been considered to derive a reduced set of equations. As a consequence, the model uses only a small number of parameters although sufficient to describe the experimental results of interest.

The following simplifying assumptions have been made:

About the immune escape mechanism, we assumed that the growth of the tumor mass induces a progressive state of immune tolerance; this process has been reported in the 1970s by several authors and later attributed to a number of biological processes such as immune cell exhaustion, the activation, and accumulation of inhibitory immune cells such as Tregs, MDSCs, and M2-type tumor-infiltrating macrophages (8–11); other interesting contributing factors to immune escape, such as the selection of poorly immunogenic variants, have not been modeled because their additional explanatory power has been estimated relatively modest in comparison with their mathematical complexity;

As a consequence of the immune escape mechanism, the activity of the immune response is a non-monotonic (or biphasic) function of the tumor mass; it increases with tumor mass as long as the volume remains below a certain threshold, then it decreases; this is coherent with results in refs. 9–11, 28–32. The coupling between the tumor dynamics and the immune response is obtained mathematically by the product of two mechanisms. In the first mechanism, the tumor produces antigens that stimulate the immune response in a tumor size-dependent manner; this explains the increasing part of the immune response curve as a function of tumor size. The second mechanism is a tumor size-dependent neutralization of the immune response, which produces the decreasing part of the immune response curve after the tumor-size has reached a critical threshold, where the effect of immune neutralization exceeds the effect of immune stimulation;

The nonimmunologic effects of radiotherapy have been described by a linear-quadratic model (5–7), which we have modified to describe more finely the tumor mass dynamics after irradiation;

On the stimulation of the immune response by irradiation, the different explanations given by several authors (15, 33–37), even though they differ by their biological nature, have been assumed to produce additive effects that can be approximately described by a unique equation, which involves the release of tumor antigens by radiotherapy;

The downregulation of the immune response by irradiation is modeled by two mechanisms. First, radiotherapy may kill the immune effectors in the irradiation field. This is a mathematical interpretation of preclinical data published in ref. 37. Second, the immune effectors downregulate the immune response; because effectors are in relation with tumor antigens, and therefore the radiation dose, this creates a dose-dependent immune downregulation, in coherence with the preclinical data published in ref. 22, involving the PD1–PD-L1 axis, or with other biological models involving the TGFβ signaling pathway (20, 21);

Because the blockade of the PD1–PDL1 axis promotes the antitumor activity of the T-cell response by reducing the expansion of tumor-infiltrating MDSC and T regulatory lymphocytes (1), we made the anti-PD1 mAb serum concentration an inhibitor of the immune response downregulation;

Because the immune response in anti-CTLA4 mAbs treatment is thought to be mostly caused by an increased activity of CD4

^{+}T cells and the subsequent stimulation of the secondary (or memory) immune response (1), we made the anti-CTLA4 mAb concentration a stimulating factor of the memory response and we neglected its direct contribution to the primary immune response (even though CTLA4 blockade does also have a modest effect on the primary immune response, which is caused by the increased inhibition that CTLA4 blockade exerts on the growth of the tumor mass; this smaller tumor mass reduces the inhibition of the primary immune response by the tumor mass according to Eq. D). Hence, our mathematical formulation of the immune response has been separated into two additive components: the primary immune response and the secondary (or memory) immune response, each part corresponding to its eponymous biological counterpart; the secondary response is assumed to appear progressively after the primary response, as a result of the cumulative effect of the primary response. The secondary response appears relatively slowly but has a long-lasting effect and it is not sensitive to tumor downregulation, contrary to the primary response. This is our mathematical interpretation of immunologic processes for which one can find a comprehensive review in ref. 1;

From these biological hypotheses, we derived the set of equations presented in the next section.

The model is described by 5 discrete-time equations:

Here, these equations describe the changes over a one-day period, between day *n* and day *n* + 1, though any time step could be chosen. A visual representation of the links between the variables is provided in Fig. 1, while the detailed explanations are provided in the next section.

## Materials and Methods

### Tumor dynamics

The compartment *T _{n}* is the tumor mass on day

*n*, expressed in grams. The parameter μ > 0 is constant and expresses a constant exponential growth. Hence, over a one-day period, in the absence of treatment and assuming no activity of the immune system, the tumor mass is multiplied by exp(μ). This simplification is adequate to describe the experimental results presented here, involving small-size murine tumors, even though other tumor mass dynamics could be used, such as the Gompertz growth model for example, because the cell loss factor is probably volume dependent in human tumors as explained in refs. 38 and 39.

The immune system is the only antitumor agent explicitly considered: Depending on its intensity, it can merely slow down tumor growth or can provoke a shrinkage of the tumor mass and eventually eradicate it. The tipping point between these two behaviors is Z_{prm,n}. + Z_{sec,n}. = μ. If Z_{prm,n}. + Z_{sec,n}. < μ the immune system slows down tumor growth, if Z_{prm,n}. + Z_{sec,n}. > μ the immune system makes the tumor mass shrink.

The constant μ implicitly includes other cell loss factors such as, for example, the necrosis of tumor cells resulting from inadequate vascularization and the resulting lack of nutrients or oxygen. Hence, μ could be made a dynamic parameter in a more complex model involving, for example, tumor angiogenesis.

The quantity *S*_{n}(** d**) denotes the survival probability at time

*n*for a cumulative radiation dosing regimen

**(**

*d*=*d*

_{1},…,

*d*

_{n}), where

*d*

_{i}is the radiation dose received on day i (

*d*

_{i}being possibly zero). Its expression is derived from our modified formulation of the linear-quadratic model (see the Supplementary Material section for explicit mathematical derivation and expression).

### Antigen dynamics

The compartment A_{n} models the amount of tumor antigens. Antigens are released either naturally by tumor cells at a constant rate ρ or after irradiation at a rate ψ.[1−*S*_{n}(** d**)], where

*S*

_{n}(

**) is the radiation survival probability as defined previously. The constant rate ρ controls the intrinsic (nonradiation-induced) immunogenicity of the tumor, while the constant ψ controls the radiation-induced immunogenicity. The constant parameter λ∈[0,1] controls the rate at which these antigens disappear from the system and produce immune effector cells L**

*d*_{n}.

Note that despite the use of the word antigen, the proposed model is actually agnostic about the precise biological mechanism by which radiation stimulates an immune response: As long as this stimulation is approximately proportional to the number of cells killed, the model holds and the compartment A_{n} may describe a time delay between irradiation and immune stimulation.

### Lymphocytes (or immune effectors) dynamics

The compartment L_{n} models the immune effectors at the tumor site: tumor-infiltrating lymphocytes (TIL), tumor-associated macrophages (TAM), etc. This compartment is fed at a constant rate λ > 0 from the compartment of the tumor antigens A_{n}, that are either produced naturally or released by irradiation. Immune effectors exit the system at a constant rate ϕ > 0.

The immune cells that are located in the tumor when radiation is administered are killed and this is modeled by the survival probability δ_{n}, for which we might have used a linear-quadratic model. However, for simplicity we have assumed that no immune effector in the tumor survives the doses used in clinical practice, hence: *δ*_{n} = 0 if radiation is applied at day *n* on the tumor, *δ*_{n} = 1 otherwise.

### Primary immune response

The variable Z_{prm,n} models the primary activity of the immune system. As its name indicates, it is our mathematical formulation of the biological concept of primary immune response to the tumor, that is, it is biologically composed of an innate and of an adaptive immune response (but this breakdown is not reflected in this model).

According to our model, the primary immune response is proportional (with the constant ω) to the number of immune effectors L_{n} multiplied by the ratio of immune effectors that are active against the tumor, defined as 1/[1+(κ.T_{n}^{2/3}.L_{n})/(1+p_{1})]:

The constant κ models the intrinsic ability of the tumor to downregulate the immune system. For example, it could contain the propensity of the tumor to express immunomodulatory proteins along with the propensity of immune cells to express receptors for these factors, such as PD-L1 and PD1, respectively. *Ceteris paribus*, the higher the tumor volume, the higher the amount of immunomodulatory factors produced by the tumor, hence the term T_{n}^{2/3}. The power 2/3 is used because it provided a good fit of experimental data on tumor immunogenicity; it links the amount of immunomodulatory factors to the surface of the tumor rather than its volume, which makes sense if one considers that the center of a tumor may contain necrosis and should be less infiltrated by the immune effectors. *Ceteris paribus* again, the higher the number of immune effectors infiltrated in the tumor site, the higher the amount of immunomodulatory factors produced by the tumor cells as a result of increased interactions with the immune system, hence the term L_{n}. *In fine*, the term κ.T_{n}^{2/3}.L_{n} can be thought as an amount of immunomodulatory factors produced by the tumor (or by the immune effectors themselves), in response to the attack by the immune effectors, multiplied by a factor that models the affinity for their receptors.

The variable *p*_{1} is the concentration of an anti-PD1 or an anti–PD-L1 immunotherapy drug, it depends from the dosing regimen and pharmacokinetics, which precise description is out of the scope of this article.

A more detailed explanation of the mathematical formulation is provided in the Supplementary Material, along with a discussion about the concurrent blockade of PD1 and PD-L1, and of the consequences of the existence of several concurrent ligands such as PD-L1 and PD-L2.

### Secondary (or memory) immune response

The variable *Z*_{sec,n} models the secondary activity of the immune system, that is, it corresponds to the biological concept of the memory response after a first encounter with a pathogen. The mechanism by which an anticancer secondary (or memory) immune response establishes itself following a primary immune response is complex and not fully understood. However, we used a small number of probable facts to base our mathematical formulation. For a comprehensive review of these mechanisms, the reader is referred to ref. 1.

The secondary response is probably mediated by auxiliary or helper T cells, or CD4^{+} T cells. The secondary response is conditional to the existence of a primary response and is a progressive process that takes time to complete. Once it is established, the secondary response is long-lasting and relatively insensitive to downregulation by the PD1–PD-L1 axis. Because the inhibitory receptor CTLA4 is expressed mostly on CD4^{+} T cells, we decided that the effect of CTLA-4 blockade would be to inhibit solely the secondary immune response (even though this is probably a simplification). Then, we assumed that for each time step, a constant proportion γ of the activated (or non-exhausted) effectors of the primary immune response will incrementally contribute to the secondary immune response. To secure the long-lasting memory effect, we further assumed that the intensity of the secondary response could not decrease, hence we did not use a clearance rate.

Hence, the secondary immune response follows the primary response, and its intensity is related to the cumulative effect of the primary response. Mathematically we defined the secondary response as the integral of the primary response over time; this integral is simplified into a discrete sum because of the discrete-time formalism.

The variable *γ* is the propensity of the primary response to induce a secondary response per time-step. Biologically, this propensity could be related, among other things, to the number of CD4 lymphocytes (excluding the regulatory ones).

The variable *c*_{4} is the concentration of an anti-CTLA4 immunotherapy drug, it depends from the dosing regimen and pharmacokinetics, that are out of the scope of this article.

The constant *r* is chosen greater than 1. Hence, in the absence of the anti-CTLA4 drug, the propensity to induce a secondary immune response per time step is γ/*r* < γ. If an anti-CTLA4 drug is given, this propensity tends to γ by inferior values.

The Supplementary Material contains more details on the mathematical formulation of the effect of the anti-CTLA4 agent on the secondary immune response.

## Results

### Abscopal effects and probability of metastatic rejection

Abscopal effects have been investigated for decades and there are numerous experimental data to choose from. We studied the results in ref. 9 because they provide very detailed data on the link between a tumor size and its propensity to induce an immune response.

After tumor implantation in mice, the author has quantified experimentally the development and decline of immune resistance (or abscopal effect) under the influence of the progressively growing tumor. Within the context of our model, we could replicate his findings and we correctly explain that, as the tumor grows, the percentage of tumor rejection initially increases, then reaches a plateau, ultimately followed by a rapid decline. As shown in Fig. 2, despite the simplifying assumptions of the proposed model, the theoretical tumor rejection probabilities are relatively close to the experimental data for a wide range of tumor volumes; hence, the proposed model was able to describe correctly the immunogenicity of a tumor as a function of its size. The reader is referred to the Supplementary Material section for more details on the mathematical approach that we used to derive a probability of tumor rejection from our dynamic model.

### Maintenance of a sustained response after treatment discontinuation

#### Special case: no secondary immune response.

This case applies to a patient whose immune system has been so deeply depleted by disease or by chemotherapy that the remaining count of CD4^{+} T cell is too low to generate any significant secondary immune response. For results on the depth and duration of CD4 lymphopenia induced by chemotherapy, see for example ref. 40. It is also well known by clinicians that the disease itself, even before chemotherapy, may cause a marked lymphopenia. In that case, Eq. E can be entirely removed since *Z*_{sec,n} ≈ 0.

In this special case, tumor control cannot be maintained after treatment discontinuation and recurrence is highly expected (see Fig. 3B). According to our model, disease control for this type of patient requires theoretically a perpetual blockade of the PD1–PDL1 axis. We will see that this is not the case if a secondary immune response can be stimulated by the primary response. In clinical practice, it is reasonable to assume that at least some residual level of secondary immune response persists, which motivates the more general model in the next section.

#### General case: primary and secondary immune responses.

In the general case, disease recurrence may or may not occur after treatment discontinuation depending on its duration and the propensity to induce a secondary immune response (see Fig. 4).

In the proposed model, the concurrent blockade of CTLA4 accelerates the appearance of a long-lasting secondary response. The long-term persistence of the immune response after treatment discontinuation has been noted in several clinical trials of the blockade of the PD1–PDL1 axis, see for example ref. 41.

### Synergy of radiotherapy and immunotherapy

#### Radiotherapy and PD-L1 blockade.

Our model explains the synergistic effects of immunotherapy in combination with radiotherapy. In Fig. 4, we provide a numerical example inspired from various recent publications, such as ref. 12 or 22; this result shows that the model explains why combined treatment can result in sustained remission, whereas radiotherapy or immunotherapy alone only provide a temporary relief followed by recurrence.

#### Synchronization of radiotherapy and PD1-PDL1 blockade.

As shown experimentally in ref. 22, the model predicts that radiotherapy and PD1 or PD-L1 blockade should have sufficient synchronization for a significant synergistic effect to exist. In the proposed model, the rationale for synchronization relies on the finite quantity of tumor antigens released by irradiation, and the finite lifetime of immune effectors. As time after radiation passes, the quantity of antigens and the immune effectors both decrease and after a while, this reduces the efficacy of immunotherapy that does not produce a substantial synergistic effect (Fig. 4). However, in clinical practice, the concomitant administration of several treatments may result in increased toxicity and adverse events; therefore, more analysis is necessary to derive optimized dosing schedules that are both efficient and have an acceptable tolerability profile.

#### Radiotherapy and concurrent PD1 and CTLA4 blockades.

Results published recently in ref. 14 provide many experimental data on radiotherapy with the concurrent blockades of the PD1–PDL1 axis and/or the CTLA4 pathway. These works explore the abscopal immune effect of radiotherapy if only a fraction of the total tumor burden is irradiated, and how this immune effect may be magnified by blockades of checkpoint inhibitors. A numerical example directly inspired from their results is provided (Fig. 5, left), along with a theoretical simulation with a less aggressive tumor (with a longer doubling time, all other parameters being identical: Fig. 5, right).

The model was able to produce survival curves similar to those produced experimentally by these authors, hence showing that it is capable of modeling the synergy between radiotherapy and PD1/PD-L1 or CTLA4 blockades. The simulated example with the (theoretically) less aggressive tumor shows how the model can be used to compare different modalities while assuming various hypotheses about the underlying tumor.

The Kaplan–Meier survival curves have been produced numerically by assuming (i) that the death automatically results from the tumor mass exceeding a fixed threshold (set to 10 times the mass on day 1), and (ii) that the tumor growth rate is a random variable with a log-normal probability density function[whose SD has been set such that 95% of the tumors have a growth rate μ within (80%–120%) of the peak value]. Despite these simplifications, the model adequately describes the experimental survival curves. The possibility to simulate mathematically the Kaplan–Meier curves is both practical and intuitive for most clinicians.

## Discussion

For decades, improving cancer therapy has been achieved by trial-and-error practices. Most phase III registration studies or multicentric clinical trials calling for a change in practice have been designed on empirical basis. In some cases, such empirical designs have led to approving new drugs with suboptimal dosing (for instance, sunitinib as a single agent has been a major breakthrough for treating renal carcinoma, but it took years of bedside practice to learn that the official dosing was probably not the optimal one; ref. 42). Defining optimal modalities of treatment turns to be even more complicated when several strategies or multitherapies are to be combined. Today, the ever-growing complexity of cancer biology and the plethora of new agents being regularly approved makes determination of optimal dosing and scheduling a challenging task that cannot be addressed anymore by empirical approaches. Consequently, developing computational-aided strategies to better design combined therapies is attractive, because *in silico* simulation can help to test numerous modalities in a time and cost-effective fashion.

The recent achievements in the field of cancer immunotherapy have been recognized as a major breakthrough. However, despite very encouraging results, improvements of the response rates and survival are still needed because many patients will fail in maintaining sustained clinical benefit. In this respect, association of immune checkpoint inhibitors with other strategies (e.g., chemotherapy, radiotherapy) is a promising perspective, as the most striking past improvements in cancer therapy have been achieved thanks to the combination of a variety of strategies. Indeed, it is widely acknowledged now that beating cancer will not come from a single miraculous agent, but rather from the rational combination of a number of complementary strategies that could be chosen from various types of immunotherapy, radiotherapy protocols, antiangiogenic agents, targeted therapies, and different chemotherapy regimen. Ideally the dosing, duration, and relative scheduling of the different components should be optimized by taking into account the potential synergies as well as the possible antagonisms and toxicities. Given the countless possibilities of combinational therapy and the critical challenge of picturing all the interactions between complex and intricated pharmacologic processes, we believe that an integrated mathematical model could be of great help to optimize these combinations *in silico* prior to test them in patients. In addition, such a mathematical model could help to choose strategies that are both medically and economically justified because they are the most likely to eventually yield clinical benefit.

For the first time (to our knowledge), our study introduces such a model that integrates radiotherapy combined with two important types of immunotherapy: the blockades of the PD1–PDL1 axis and of the CTLA4 pathway. This model has been designed to rely on a limited number of parameters that have a biological relevance (if not an exact biological equivalent) and it can be calibrated on experimental data. The mathematical complexity has been purposely kept to a minimum by the use of a small number of discrete-time equations. Indeed, we do believe that phenomenologic modeling will be more easily transposed to routine clinical practice than bottom-up, systems biology approaches (43).

The scientific relevance of the proposed model could be checked by its ability to explain and to reproduce *in silico* a number of experimental results, from decades-old data describing the link between tumor immunogenicity and its volume, to recent and cutting-edge results on the combination between different types of immunotherapy and radiotherapy protocols. Also, the model has a certain degree of robustness with respect to biological hypotheses because the complexity of the immune system has been summarized as an interaction between antigens and a compartment of immune effectors. Though this simplification does not capture all the complexity of the interplay between the many actors of the innate and the acquired immune response, its low level of granularity is sufficient to model the effect of PD-L1 or CTLA4 blockades. The practical interest of this model relies on its capacity to compare the outcome of different strategies and to allow optimizing dosing regimen *in silico*. For example, if the disease is deemed so aggressive that it mandates the use of very risky strategies that might result in poor tolerance or even death by treatment, such an integrated model can be of tremendous help to (i) quantify and justify the benefit-risk ratio of such an option with respect to other less risky strategies, and (ii) optimize the aggressive strategy to minimize its associated risks while preserving the chances of success.

### Conclusion and perspectives

The mathematical model presented can describe several experimental results on the combination between radiotherapy and immune checkpoint inhibitors targeting the PD1–PD-L1 axis and the CTLA4 pathway. Beyond the intellectual interest of achieving a mathematical description of experimental results, this model could have a practical interest by predicting the clinical outcome of new combination of anticancer therapies and to test different scenario. In that respect, the ability to simulate Kaplan–Meier curves is a practical and intuitive way to compare different designs *in silico*, whereas the search of optimized protocols could be partially automated by using mathematical tools.

Future works will focus on the integration of such mathematical tools into clinical practice, and will therefore explore how to manage drug and radiation-induced toxicities in the optimization process. Dose–effect relationships will be investigated and a pharmacokinetic model will be integrated to better describe the relation between the experimental results and the actual dosing protocols. The combination with other therapeutic modalities such as metronomic chemotherapy or targeted agents (such as antiangiogenics) will also be the subject of future research.

## Disclosure of Potential Conflicts of Interest

J. Ciccolini is a consultant/advisory board for Celgene, Roche, Pierre Fabre. F. Barlesi is a consultant/advisory board member for BMS, Roche/Genentech, Astra-Zeneca. No potential conflicts of interest were disclosed by the other authors.

## Authors' Contributions

**Conception and design:** R. Serre, N. André, X. Muracciole, D. Barbolosi

**Development of methodology:** R. Serre, C. Meille, N. André, F. Barlesi, D. Barbolosi

**Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.):** R. Serre, L. Padovani, X. Muracciole, D. Barbolosi

**Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis):** R. Serre, L. Padovani, N. André, F. Barlesi, D. Barbolosi

**Writing, review, and/or revision of the manuscript:** R. Serre, S. Benzekry, C. Meille, N. André, J. Ciccolini, F. Barlesi, X. Muracciole, D. Barbolosi

**Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases):** R. Serre, J. Ciccolini, D. Barbolosi

**Study supervision:** R. Serre, X. Muracciole, D. Barbolosi

## Acknowledgments

The authors thank the Ligue Contre le Cancer de Corse du Sud for its support.

## Grant Support

This work was partly supported by the A*MIDEX Mars Project Funding. This study was supported by a public grant from the French State, managed by the French National Research Agency (ANR) in the frame of the "Investments for the future" Program IdEx Bordeaux - CPU (ANR-10-IDEX-03-02).

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.

## References

^{+}t cells: changing strategies for cancer treatment

^{+}t-lymphocyte regeneration after intensive chemotherapy