Sources of inter-individual variability leading to significant changes in anti-PD-1 and anti-PD-L1 efficacy identified in mouse tumor models using a QSP framework

While anti-PD-1 and anti-PD-L1 [anti-PD-(L)1] monotherapies are effective treatments for many types of cancer, high variability in patient responses is observed in clinical trials. Understanding the sources of response variability can help prospectively identify potential responsive patient populations. Preclinical data may offer insights to this point and, in combination with modeling, may be predictive of sources of variability and their impact on efficacy. Herein, a quantitative systems pharmacology (QSP) model of anti-PD-(L)1 was developed to account for the known pharmacokinetic properties of anti-PD-(L)1 antibodies, their impact on CD8+ T cell activation and influx into the tumor microenvironment, and subsequent anti-tumor effects in CT26 tumor syngeneic mouse model. The QSP model was sufficient to describe the variability inherent in the anti-tumor responses post anti-PD-(L)1 treatments. Local sensitivity analysis identified tumor cell proliferation rate, PD-1 expression on CD8+ T cells, PD-L1 expression on tumor cells, and the binding affinity of PD-1:PD-L1 as strong influencers of tumor growth. It also suggested that treatment-mediated tumor growth inhibition is sensitive to T cell properties including the CD8+ T cell proliferation half-life, CD8+ T cell half-life, cytotoxic T-lymphocyte (CTL)-mediated tumor cell killing rate, and maximum rate of CD8+ T cell influx into the tumor microenvironment. Each of these parameters alone could not predict anti-PD-(L)1 treatment response but they could shift an individual mouse’s treatment response when perturbed. The presented preclinical QSP modeling framework provides a path to incorporate potential sources of response variability in human translation modeling of anti-PD-(L)1.


Plasma (Central) Compartment
The antibodies are cleared from the bloodstream at the rates k el,A P and k el,A L for anti-PD-1 and anti-PD-L1 respectively. Distribution to the peripheral compartment occurs at rates k cp,A P and k cp,A L while flow back from the peripheral into the central compartment occurs at rates k pc,A P and k pc,A L . Flux between the central and tumor compartments (k ct,A P and k ct,A L ) are described in Section 1.1 and is based on the interstitial tumor volume (V t ).

Transfer to tumor compartment
Antibody distribution to the tumor is based on the work of Li et al. [6] which models movement of the antibodies as a function of the permeability of the capillary and diffusion into the tumor tissue. The capillary radius (R cap ), the average radius of tissue surrounding each blood vessel (R Krogh ), and the permeability constant (P A P ) determine movement through the capillary membrane while the diffusion constant (D A P ) and the tumor radius (R tumor ) determine the rate of diffusion. Both permeability and diffusion constants are estimated based on the molecular weights of a bivalent antibody of 150 kDa [6]. The fraction of tumor volume that is outside of cells (the intersitial volume) is denoted by .

Peripheral Compartment
Antibody amounts in the peripheral compartment (denoted with a p subscript) are purely functions of influx from and outflow to the central compartment. No synthesis or degradation is modeled.

Tumor Compartment
We assume the tumor is well mixed and that T cell concentrations are constant as tumor volume changes. The fraction of tumor volume not comprised of cells (the interstitial volume) is .

Cells
In the tumor compartment the concentrations of tumor cells (T i ), PD-1+CD8+ T cells (E i ), and other PD-1+ T cells (H) were modeled. Tumor growth is modeled using a logistic growth curve due to its ability to describe tumor growth past the initial stage of exponential growth [5] and for consistency with previous mouse tumor modeling publications [2,7] with growth rate k pfr,T and carrying capacity k limit,T . Tumor cells (T V ) grow according to this curve until they are damaged either through ADCC (f ADCC () ) due to binding with anti-PD-L1 or through interaction with a CTL with rate k kill . Damaged tumor cells (T M ) undergo apoptosis with rate k apop . The total tumor volume (V t ) is the sum of the undamaged and damaged tumor cells divided by the constant tumor cell concentration (Eq. 16). Additional volume from immune cells is assumed negligible due to their smaller size and lower concentration (5-7 um, 2.15 × 10 2 cells per mm 3 , vs. 12-25 um [4], 10 8 cells per mm 3 [1] for tumor cells). For numerical stability, a minimum tumor size of 1 mm 3 was implemented.
Anti-PD-L1 molecules capable of ADCC damage tumor cells according to the number bound per tumor cell ( A:L V T V ) with a maximum rate of E ADCC max and achieving half that rate at EC50 ADCC .
Other PD-1 positive CD8+ T cells (H) cells serve only as a source of PD-1 to competitively bind with PD-L1 on tumor cells, and are synthesized at the rate k s,H and degrade at the rate k d,H . We assume that the immune cell concentration is constant as the tumor volume changes, as such the synthesis rate is multiplied by the tumor volume (V t ).
Inside the tumor compartment, PD-1+CD8+ T cells are synthesized as inactive cells (E I ). The cells go through eight proliferation stages (E A , E P i for i = 1, ..., 7) before becoming CTLs (E D i for i = 1, ..., 10) [7]. Eight proliferation stages was chosen to match the magnitude of T cell expansion seen in [8]. Each CTL can damage ten tumor cells before becoming exhausted (E D 0 ) [3]. Inactive, active, and proliferating CD8+ T cells degrade at a rate of k d,CD8 while CTLs degrade at a rate of k d,CT L .
The binding of PD-L1 expressed on tumor cells to PD-1 expressed on inactive PD-1+CD8+ T cells (P I : L V ) controls T cell behavior. This binding inhibits the influx (f influx ) of inactive PD-1+CD8+ T cells into the tumor microenvironment, and the activation (f activate ) of PD-1+CD8+ T cells in the tumor compartment. In the absence of antibody treatments PD-1:PD-L1 per inactive CD8+ T cell is at its steady state concentrations (approximately 7.1) causing inactive PD-1+CD8+ T cells to enter the tumor compartment at the rate of k s,I V t (Eq. 56). This influx rate ensures that the CD8+T cells reach a steady state concentration during control simulations.
Both the activation and influx functions reach their maximum (E i max ) when there is no bound PD-1:PD-L1. As the bound per inactive T cell increases, the activation and influx decreases. Since both the activation and influx functions are dependent on PD-1:PD-L1 binding on inactive CD8 cells, we use the same EC50 and Hill coefficient (c=1) for both.
In the presented model set up, during control simulations, nearly all tumor cells are undamaged and the inactive CD8+ T cell concentration converges to steady state value from data (section 4) leading to a consistent value for P I :L V E I while fitting the remaining parameters. While for each individual the exact receptor occupancy per inactive CD8+ T cell may not be 6.6, the value works well for maintaining influx equal to k s,I V t at baseline and higher when receptor occupancy drops with treatment.
Activated CD8+PD-1+ T cells (E A ) follow the same activation rates as described above. Cells deactivate at the rate k A2Ibasal and double with rate k pfr ; for each cell leaving E A , two cells enter E P 1 .
After proliferating, the cells are now cytotoxic lymphocytes and have the ability to damage tumor cells. Every time a cell damages a tumor cell it moves from stage In stage E D 0 , the cell is considered exhausted, cannot damage any tumor cells, and is left to degrade. CTLs damage tumor cells at the rate k kill .

Free Antibodies
Antibody amounts in the tumor (denoted with a t subscript) are dependent on flux from the central compartment (see Section 1.1) and binding to and unbinding from receptors on the various cells within the tumor compartment at rates k on,AP , and k off,AP for anti-PD-1 and k on,AL , and k off,AL for anti-PD-L1. PD-1 on H or E cells are denoted P i for i ∈ {H, I, A, P, D} while PD-L1 on tumor cells are denoted L V and L M for tumor cells and damaged tumor cells respectively. No degradation of free antibodies happens within the tumor compartment, but antibodies bound to receptors on cells can be internalized (described in Section 3.5). Antibodies bound to cells are denoted as A:P i for i ∈ {H, I, A, P, D} for anti-PD-1 bound to H or E cells while A:L V and A:L M denote anti-PD-L1 to PD-L1 complexes on tumor cells and damaged tumor cells, respectively.

PD-1 Receptors
PD-1 receptors on cells are tracked according to the following assumptions: • New cells are created with the average number of surface receptors per cell (S H , S E for H and E cells respectively).
• PD-1 receptors degrade at the rate k d,P .
• Cells synthesize PD-1 at a rate k s,P,H = S H × k d,P and k s,P,CD8 = S E × k d,P for H and E cells, respectively.
• When a cell degrades, all receptors on the cell are eliminated.
• When a tumor cell that is bound to an PD-1+ cells dies (through apoptosis or the death term from the logistic growth equation), the PD-1 receptors return to the unbound pool.
• To minimize the number of equations, PD-1 receptors on proliferating (E P i ) and CTL (E D i ) cells are gathered into their respective pools (P P and P D ). The PD-1 receptors are assumed to be evenly distributed between the cells in the group.
• When a cell moves one stage downstream (for example, E P 7 to E D 10 ), the receptors on the cell also move one stage downstream.
Comments are provided for the first two equations to explain each term.
+ k off,PL (P P :L V + P P :L M ) − k on,PL (P P L V + P P L M )/(V t ) + k off,AP A:P P − k on,AP A P P P /(V t ) + k pfrCD8 (2P A + A:P A + P A :L V + P A :L M ) + k pfrCD8 (P P + A:P P + P P :L V + P P : + k pfrCD8 (2P P + A:P P + P P :L V + P P :L M ) • When a cell degrades, all the receptors on the cell are eliminated • When a PD-1+ cell that is bound to a tumor cell dies, the PD-L1 receptors return to the unbound pool

Bound antibodies on cells
The following equations follow the same assumptions described in sections 3.3 and 3.4

PD-L1 bound to tumor cells (A:L V ) and damaged tumor cells (A:L M ).
dA: − k A2Ibasal A:P A + f activate P I :L V E I A:P I − k pf rCD8 A:P A dA:P P dt = − k off,AP A:P P + k on,AP A P P P /(V t ) − k d,AP A:P P − k d,CD8 A:P P (42)

Bound PD-1:PD-L1 on cells
The following equations follow the same assumptions described in sections 3.3 and 3.4 Parameter fitting using baseline intra-tumoral T cell concentration data T cell kinetic parameters were constrained using average concentrations of CD8+ T cells from untreated CT26 tumors in Balb/c mice ( Figure 1C, Table S1) that were derived from flow cytometry data. The five CD8+ T cell populations were defined by cell surface markers and their concentrations were derived as follows: marker concentration = viable cell count (10 6 per ml) percentviable × total events × CD45 events × CD8 as %CD45 × marker as %CD45 The model set up for the CD8+ T cell stages has the potential to accurately capture the magnitude and speed of T cell expansion after anti-PD-(L)1 treatment, but has a number parameters not readily available in the literature. These include k pfrCD8 , k d,CD8 , k d,CTL , E activate max , and k s,I . Our goal is to constrain these parameters so that the model predicts accurate CD8+ T cell concentrations for each stage during control simulations.
The known concentrations are See T cell concentration data reference for the specific numbers. Note that these are concentrations, while our original model equations are in amount. We assume that as the tumor grows, the concentrations of CD8+ T cells remains constant. This is true in model simulations due to the k s,I V t term in equation 56. We write them now as concentrations for use in this exercise.
Since we are interested in the total number of CD8 cells, and that does not change as they travel through the D i compartments, we can use [D] = 10 i=0 [E D i ] to simplify equations 60-61 as follows: Step Step 2: Constrain k d,CD8 in terms of k pfr,CD8 From equation 65 we get: 2k pfr,CD8 k pfr,CD8 + k d,CD8 [E A ] is known, so use a numerical solver to get k d,CD8 = k pfr,CD8 /2.3581.
Step 3: Constrain k d,CTL in terms of k pfr,CD8 From equation 63: 2k pfr,CD8 k pfr,CD8 + k d,CD8 We solved for k d,CD8 k pfr,CD8 numerically in step 2. Combine these two equations for the result k d,CTL = 11.7403 × k pfr,CD8 . Step Parameters k s,I and k A2Ibasal cannot be uniquely identified with high confidence given the limited data available. Thus we fixed k A2Ibasal to the value used in [7] and calculated k s,I to match. Step We know [E D ] and [E D 0 ] and we know k d,CTL from step 3. Note that [ The exact ratio of undamaged to damaged tumor cells varies based on many parameters, but in the majority of simulations, almost all tumor cells are undamaged, so we can approximate [T V ] ≈ 0.9999 × C T C . Thus, we can find k kill with a numerical solver.