Mechanistic modeling predicts efficacy of CISH knockout in tumor-infiltrating lymphocytes with synergistic gene editing

Tumor-infiltrating lymphocyte (TIL) therapy is a type of adoptive cell therapy, where the lymphocytes of a cancer patient’s tumor are harvested, expanded in vitro using IL-2 stimulation, and then infused back into the patient[1], [2]. However, even with the use of TIL therapy, cancer cells can survive for various reasons, such as poor lymphocyte infiltration into tumors, chronic activation of the T cell receptor and the immunosuppressive tumor microenvironment [3]. Cytokine-inducible SH2-containing (CISH) protein is a negative regulator of T cell activation, and in a recent clinical trial was knocked out in TILs to improve TIL therapy efficacy[4]. In this study, we developed a mechanistic signaling pathway model to theoretically evaluate the efficacy of CISH knockout (CISH KO) in T cell activation and examine potential alternative target genes that can theoretically be targeted using multiplex gene -editing or drugs to further improve T cell activation and function[5]. Our modeling results demonstrate that CISH knockout increases the transcription of activation biomarkers IL -2 and TNF-α, but also inhibitory biomarkers such PD1 and FasL. Using global sensitivity analysis, we also found that GSK3B, which is responsible for the deactivation of the NFAT, is also predicted to further increase T cell activation when knocked out. In addition, we predict that PDCD1, FAS and CTLA4 can be knocked-out in combination with CISH to further enhance T cell activation and prevent exhaustion and apoptosis.

Adoptive cell therapy (ACT) is a form of personalized cancer treatment that consists of harvesting immune cells, expanding them in vitro and transferring them back to the patient to eliminate tumor cells[1]. There are two main approaches; one where immune cells are harvested (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint directly from cancer patients to treat their own disease, i.e. autologous cell therapy, or another where immune cells are obtained from healthy donors, i.e. allogeneic cell therapy[6]. A major advancement in the efficacy of ACT was in 1976 when it was discovered that the administration of interleukin-2 (IL-2) can expand T cells ex vivo[2]. Moreover, Rosenberg et al. showed that high doses of IL-2 result in the generation of lymphokine-activated killer (LAK) cells that mediate the regression of tumor s in mice [7]. Donohue et al. examined the effect of administration of IL-2 after adoptive cell transfer [8]. Despite their initial promising results in a clinical trial, where 11 of 25 patients with different types of metastatic cancer , including melanoma, colorectal, and renal-cell, showed tumor regression after the administration of LAK, a subsequent clinical trial did not demonstrate a statistically significant benefit compared to the administration of IL2 alone in patients with renal cell cancer[9], [10]. Tumor-infiltrating lymphocyte (TIL) therapy is a novel cancer immunotherapy that patient’s own immune cells, present in the tumor microenvironment, are harvested, expanded in vitro and reinfused back to the patient to fight cancer. TILs after expansion with IL-2 showed 50 to 100 times higher efficacy compared to LAK cells in lung and hepatic metastatic tumors, although only with pretreatment with cyclophosphamide [11]. The first example of TIL therapy was in 1988, when Rosenberg et al. used autologous TILs along with the administration of IL -2 and pretreatment with cyclophosphamide to treat 20 patients with metastatic melanoma[12]. This study showed that the administration of TILs resulted in positive response in 9 out of 15 patients who were not treated with IL -2 before and in 2 of 5 who were treated with IL -2 in the past, compared to IL-2 and cyclophosphamide alone. A significant milestone in TIL therapy happened in 2002 when it was shown that the use of nonmyeloablative chemotherapy regimen as (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint conditioning to achieve lymphodepletion before TIL infusion in patients with metastatic melanoma causes tumor regression[13]. Neoantigen-reactive TIL therapy leverages the subset of tumor -infiltrating lymphocytes that are specific to the neoantigens expressed by tumor cells. Selection of neoantigen-reactive TIL pools can be accomplished using immunoassays to evaluate reactivity against candidate neoantigens nominated by next-generation sequencing of patient tumors[14]. Neoantigen-reactive TIL therapy has shown the potential to cure metastatic cancer. Zacharakis et al. presented the case of a patient with metastatic breast cancer who was treated with TILs specific to neoantigens created by nonsynonymous mutations in four genes, SLC3A2, KIAA0368, CADPS2 and CTSB, and lead to complete regression that persisted for more than 22 months[15]. Tran et al. published a case report of a patient with metastatic colorectal cancer. The patient was treated with T cells that target a neoantigen created by KRAS G12D[16]. All seven metastatic sites in the lungs regressed but one of them relapsed after 9 months of therapy because of the loss of chromosome 6 that encoded HLA-C*08:02 protein, which is necessary for the T cells to recognize the KRAS G12D tumor antigen. In 2024, the FDA approved the first TIL therapy for patients with melanoma who have been treated before with PD1 inhibitor[17]. Despite the success of neoantigen-reactive TIL therapy, it is not yet the norm. Zacharakis et al. examined the efficacy of the therapy in conjunction with pembrolizumab in patients with metastatic breast cancer[18]. TILs were isolated from 42 patients but only 6 of them were found eligible to be treated. From these patients, only 1 had a complete response which persisted over 5.5 years, 2 partial responses at 6 and 10 months, and tumor regression in 3 patients. This may be attributed to the immunosuppressive properties of tumor microenvironment in solid tumors, for example the poor infiltration of lymphocytes, the expression checkpoint inhibitors, poor (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint neoantigen presentation and chronic T cell receptor (TCR) activation [3]. To counteract these effects, there have been intensive effort s to override the intrinsic checkpoints that limit T cell activation pathways. Cytokine-inducible SH2-containing (CISH) protein is a member of suppressors of cytokine signaling (SOCS) family and a novel intracellular checkpoint[19]. It has been shown that all eight members of the SOCS family, which are SOCS1 -7 and CISH, play a key role in immune regulation[19], [20], [21], [22], [23] . Delconte et al. observed that CISH protein is a negative regulator of IL -15 in natural killer cells, and its deletion increases sensitivity to IL -15[24]. The deletion of CISH also increased production of IFN -γ and cytotoxicity against tumors, because of its interaction with JAK1 protein in the JAK -STAT signaling pathway[19], [22], [23], [25], [26], [27]. Palmer et al. also demonstrated that CISH is responsible for the ubiquitination and degradation of PLC -γ1 in T cells after TCR activation [19]. In vivo experiments showed that knocking out CISH protein increased the activation of T cells by enhancing the translocation of NFAT and NF κB in the nucleus, while it also increased the phosphorylation of ERK which contributes to higher T cell activity [28]. These results led to the development of genetically engineered TILs using CRISPR/Cas9 CISH knockout, resulting in a human phase I/II clinical trial (NCT04426669) to treat metastatic gastrointestinal cancer[4]. Here, we present a modeling study where we examine whether CISH KO can theoretically increase the activation and function of T cells, and whether there are other genetic modifications that could be used instead of , or in combination with CISH KO to achieve maximal effect on T cell activation. A mechanistic model has been previously developed using ordinary differential equations (ODEs) that captures the effects of NFAT with respect to important biomarkers for T cell activation, such as IL-2, TNF-α, CTLA-4 and Fas ligand (FasL), as described in Shin et al.[5]. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Based on their model, we expanded it to include PD1, CISH, and PLCγ1 proteins, as well as more protein interactions, to capture the effect of CISH KO in T cell function. IL-2 and TNF-α are both responsible for the activation and expansion of T cells [29], [30], [31], [32] . In addition, TNF -α plays a major role in recruiting other immune cells in the tumor site [33], [34]. On the other hand, PD1 acts as an immune checkpoint to prevent the overactivation of T cells and contributes to T cell exhaustion, while inhibiting the activation of TCR, PI3K and Ras signaling[35], [36], [37] . Furthermore, CTLA-4 has an inhibitory action on T cells by competing with CD28 for binding to the B7 protein on antigen presenting cells (APCs) and it inhibits the action of PP2A[38], [39], [40]. The role of FasL is to maintain immune homeostasis by inducing activation-induced cell death (AICD) in T cells[41]. Our results predict that CISH KO will be highly effective in elevating the expression of both activation and inhibitory markers. Moreover, our results indicate that multiplex gene-editing of CISH, PDCD1, FAS, CTLA4 may improve the efficacy of TILs, and also suggest that inhibition of GSK3β protein, which increases proliferation of T cells, could further improve T cell function in this context [42], [43], [44]. Model The formulation of the CISH model developed here is based on the paper of Shin et al.[5], originally developed to study activation -induced cell death (AICD) in T cells. In our study, the Shin et al. model was modified to integrate PLCγ1, CISH and PD1 and their interactions with other proteins (Figure 1). Enzyme-catalyzed reactions and indirect cascade effects were modeled using Michaelis-Menten reaction kinetics. For example, calcineurin activation, which is triggered by calcium release in the cytosol through the endoplasmic reticulum after the phosphor ylation of (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint PLCγ1, was simplified using a Michaelis-Menten equation (Supplementary Table A.1, reaction 6) to reduce model complexity. Similarly, the activation of TAK1, which is induced by TNF-α via tumor necrosis factor receptor (TNFR) was reduced to a single Michaelis -Menten equation (Supplementary Table A.1, reaction 43). Moreover, the Michaelis -Menten equations were decomposed into a two-step reaction, comprising a reversible binding step (𝐸 + 𝑆 ↔ 𝐸𝑆) and an irreversible catalytic step ( 𝐸𝑆 → 𝐸 + 𝑃), to enhance model stability during sensitivity analysis and maintain mass conservation. mRNA production rates were modeled as zeroth -order reactions (i.e. occurred at a constant rate), while degradation was modeled as a first -order reaction, with rates proportional to the concentration of the species being degraded. To account for the role of NFAT in protein transcription, Hill functions were used to describe the target gene mRNA expression dependence on NFAT. Mass -action kinetics were used to model association and dissociation events (Supplementary Table A.1). The strength of the model of Shin et al. lies in its detailed representation of the feedback loops that are involved in the NFAT expression[5]. In addition, our model incorporates the NFAT- mediated feedback loop for CISH expression, wherein NFAT activation enhances CISH transcription, and subsequently it leads to the ubiquitination of PLC γ1. Additionally, the model includes the regulation of PD1 expression via the NFAT translocation to the nucleus and the inhibitory effect of PD1 protein to the TCR, Ras and PI3K signaling. Our model consists of 160 free parameters and 74 ODEs and uses the same initial conditions as the original Shin et al. model. For the rest of the species, they were assumed to be 10 nM for the free intracellular proteins, 0 nM for the protein complexes, 0.1 nM for the surface proteins (e.g. CTLA-4) and mRNAs. The model was implemented in the Julia programming language[45]. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint

ODE solvers The DifferentialEquations.jl package in Julia was utilized in this study, since it provides a wide range of ODE solvers [46]. To solve the ODE system, the QNDF or the Rosenbrock23

were used [46], [47], [48] . The Rosenbrock23 solver was paired with the Tsit5 solver using an auto-switching algorithm to enable the switching between these two solvers, based on the stiffness of the ODE system[49]. The QNDF solver was used during the parameter estimation step, while the Rosenbrock23/Tsit5 solvers were used during the global sensitivity analysis and steady- state simulations because of their ability to take larger time steps, thus reducing the computational cost[50]. To further improve computational efficiency during global sensitivity analysis and steady-state simulations, numerical differentiation with central differencing was implemented. For parameter estimation, automatic differentiation was employed, providing higher accuracy but at a greater computational cost[51]. Parameter estimation Published in vitro experimental data were used to establish baseline parameter values [5], [39], [52], [53], [54], [55] . The lower and upper bounds of the parameter ranges are given in Supplementary Table A.3, while the estimated parameter values in Supplementary Table A.4. Because of large number of model parameters and the limited experimental data, many of the model parameters are non -identifiable[56]. To tackle this issue, Fides package for interior trust region reflective for boundary constrained optimization was used [57]. The objective function of the optimization problem is given by: (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint 𝜋(𝑦|ℎ, 𝜎) = 1 √2𝜎2 exp (− (𝑦 − ℎ)2 2𝜎2 ) (1) where 𝑦 is the experimental measurement, ℎ is the model observable and 𝜎 the standard deviation of the noise. A common problem in intracellular networks is the mismatch between the model outputs and the form of the experimental measurements, which in many cases are normalized or have different units. For this reason, a linear mapping between the observables and the experimental measurements was used[58]: 𝑦 = 𝑠𝑐𝑎𝑙𝑒 ∗ ℎ + 𝑜𝑓𝑓𝑠𝑒𝑡 (2) To estimate the gradient of the objective function with respect to the model parameters, forward automatic differentiation method was used. For computing the Hessian matrix, the Gauss- Newton method was used [59]. PEtab.jl library in the Julia programming language was used to implement parameter estimation[60]. The absolute and relative tolerances of the ODE solver were set to 10-3 and 10-6 respectively. Global Sensitivity Analysis To identify the sensitive parameters, the Morris method was used[61]. The Morris method is based on the calculation of elementary effects (𝐸𝐸), which are defined as: 𝐸𝐸𝑖 = [𝑌(𝑋1, 𝑋2, … , 𝑋𝑖 + 𝛥, … , 𝑋𝜅)− 𝑌(𝑋1, 𝑋2, … , 𝑋𝜅)] 𝛥 (3) (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint where 𝑌 is the output of the model, 𝑋𝑖, 𝑖 = 1, … , 𝑘 are the parameters of the model and 𝛥 is the step size, which takes values in { 1 𝑝−1 , … ,1 − 1 𝑝−1} and 𝑝 is the number of grid points in the parameter hyperspace. To sample the parameter hyperspace, 𝑟 trajectories of 𝑘 + 1 points are created, with 𝑘 elementary effects per trajectory, one for each parameter, and total 𝑟(𝑘 + 1) sample points. The trajectories are selected from a superset of trajectories to have the biggest spread in the parameter space[62]. The Morris method estimates the following sensitivity measures; 𝜇, 𝜇∗ and 𝜎: 𝜇𝑖 = 1 𝑟 ∑ 𝐸𝐸𝑖 𝑗 𝑟 𝑗=1 (4) 𝜇𝑖 ∗ = 1 𝑟 ∑ |𝐸𝐸𝑖 𝑗| 𝑟 𝑗=1 (5) 𝜎𝑖 2 = 1 𝑟 − 1 ∑(𝐸𝐸𝑖 𝑗 − 𝜇)2 (6) 𝑟 𝑗=1 Both, 𝜇 and 𝜇∗ are the means of the elementary effects and the absolute values of them respectively, and they are used to assess the overall influence of a parameter , while 𝜎2 is the variance of the distribution of the elementary effects and measures the non -linear and interactions with other parameters. Compared to 𝜇, the use of 𝜇∗prevents from potentially cancelling effect between positive and negative effects in the cases of non-monotonic or interactions with other parameters, resulting.

Model calibration (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint The model was calibrated using experimental time series data from multiple sources for CISH, PLCγ1, phosphorylated PLCγ1 (pPLCγ1), Akt1, phosphorylated MEK (pMEK), phosphorylated ERK (pERK), TNF-α, IL-2, mRNA coding for IL-2 (mIL-2), NFAT, PD1, CTLA- 4, and phosphorylated TCR (pTCR) [5], [39], [52], [53], [54], [55] . In all instances, wild -type T cells were stimulated with anti -CD28 and anti -CD3 antibodies. The CISH, PLC γ1 and phosphorylated PLCγ1 data were extracted from western blots using Fiji Image software, while the remaining data were extracted using WebPlotDigitizer[63], [64]. The initial stimulation of the T cell was modeled as a rectangular pulse function with a 2 -minute duration. Despite the large number of free parameters and the limited experimental data, Fides algorithm performed well in fitting the model to the experimental data. To estimate the parameters, 10000 starts from different initial points in the parameter space were used (Supplementary Figure A.1). NFAT, IL -2, and mIL -2 displayed longer peaks in protein expression following TCR stimulation (Figures 2A-C), compared to pAkt, pMEK, pERK, and pPLCγ1, which all exhibited a rapid increase followed by a fast reduction in intracellular concentrations (Figures 2D-F, Figure 2H). In contrast, CISH and TNF -α concentrations demonstrated a consistent increase over time (Figure 2I, 2K). However, it is important to note that the datasets for CISH, PLCγ1, and pPLCγ1 contained only four data points at early time points, potentially limiting the ability of the model to make accurate predictions at later times. PD1 concentration displays a much slower decay rate compared to the rest ( Figure 2L ). Additionally, CTLA -4 exhibits an initial peak after TCR activation and then a decrease that leads to a persistent expression (Figure 2J). The large number of model parameters and limited experimental data hinder the predictive power of the model due to identifiability and overfitting issues. In our case, parameter fitting is crucial for establishing baseline values for the intracellular biological processes that we consider (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint therapeutically unmodifiable, primarily the baseline generation and degradation rates of proteins and mRNAs, kgeni and kdegi , respectively (Supplementary Table A.1), which were held constant for the rest of the analysis. In the case of proteins that depend on the NFAT upregulation, CISH, RCAN, Carabin, CTLA -4 and PD1, and the expression of mRNA coding for NFAT, both the production rate of mRNAs and the proteins were considered as free parameters. Moreover, the Michaelis-Menten parameter, kmi, (Supplementary Table A.1) was kept constant based on the estimated values from the parameter fitting step. Using this approach, the number of free parameters was reduced from 160 to 71, which further reduced the computational cost of the sensitivity analysis. Furthermore, by limiting the number of free parameters during the rest of the analysis, we were able to obtain a better understanding of the important parts of the pathway. Critical parts of the pathway To identify if CISH is a good candidate for targeting and what other parts of the pathway can potentially be targeted independently or in conjunction with CISH KO, global sensitivity analysis (GSA) was performed using the Morris method. The outputs of Morris method are the 𝜇𝑖 ∗ (Equation 5 ), which measures the influence of each parameter, and 𝜎𝑖 2 (Equation 6 ) that corresponds to the variance of the output with respect to parameter 𝑖, which corresponds to a nonlinear or interacting behavior of the parameter. As an output for the sensitivity analysis, the steady -state concentrations of IL -2, TNF-α, PD1, FasL and CTLA-4 were used. It was assumed that T cells are, on average, in a continuously activated state in the tumor microenvironment [65], [66] . Compared to other global sensitivity analysis methods, for example Sobol, the Morris method is less accurate in terms of the ranking of the parameters, but because of the prohibitively high computational cost of the Sobol method (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint due to the number of free parameters and steady -state simulations, the Morris method was used for the purpose of this study. The expression level of IL-2 protein exhibits high sensitivity to the production velocity of CISH (v33) and the deactivation of calcineurin, v6a, with values of 𝜇𝑣33 ∗ and 𝜇𝑣6𝑎 ∗ equal to 3.36 and 3.14 respectively ( Figure 3A). Furthermore, the high variances of v33 and v6a parameters indicate a possible interaction with other parameters. The production rate of PD1, v35, and the dephosphorylation rate of pGSK3 β, v27, have a lower effect on IL -2 expression, but parameter v35 shows a high variance ( 𝜎𝑣35 2 = 142.9), similar to the production velocity of mNFAT, v6, which corresponds to the positive self-feedback loop of NFAT regulation and has a high variance, with a value 186.29. A summary of the parameters and their meaning is shown in Table 3.1. TNF-α is extremely sensitive to the deactivation of calcineurin, CISH production and dephosphorylation pGSK3β due to the high values of 𝜇𝑖 ∗ (Figure 3B), while the parameters also exhibit variances on the order of ~108. Furthermore, the parameters v27, v6a and v33 appear to influence FasL, PD1 and CTLA-4 (Figures 3C-E). According to the sensitivity analysis, PD1 and FasL are susceptible to changes in CISH production, dephosphorylation rate of pGSK3 β, and deactivation of calcineurin, but also these parameters probably interact with other parameters. On the other hand, CTLA-4 appears to be less sensitive to these parameters. In the case of CTLA -4, the dephosphorylation rate of v27 shows a high variance, which stems from interactions with other parameters or non-linear effects. Parameter v17, which influences FasL and CTLA-4, corresponds to the production of CTLA -4. An important observation is parameter v31, which corresponds to phosphorylation of IL -2 receptor (IL2R), has a high variance in FasL and PD1 respectively (variances equal to 9 .8*105 and 2.6*104), suggesting that the IL2R -PI3K-Rac1-Akt axis is potentially important. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Based on the global sensitivity analysis we confirm that theoretically CISH is an excellent target for genetic KO, via gene editing, as it is implicated in the expression of multiple distinct markers: CISH KO can increase the transcription of IL -2, TNF-α, PD1 and FasL. CTLA-4 is the most robust to parameter changes, compared to the other markers. Moreover, calcineurin appears to have an even more important role, since it directly acts on NFAT activation, while PDCD1 and GSΚ3Β gene-editing can also be used in conjunction with CISH knockout. Efficacy of CISH KO To analyze the long-term effects of CISH KO in T cells, 10000 cells were generated, under constant activation and with different parameter values, using sampling from a log -normal distribution with a mode equal to the fitted value of each parameter, 𝑋~ 𝐿𝑜𝑔𝑛𝑜𝑟𝑚𝑎𝑙(𝜇, 𝜎2) (7) where 𝜇 and 𝜎 are the mean and the standard deviation of ln(𝑋)~𝑁(𝜇, 𝜎2) distribution. The mode of the lognormal distribution is given by: 𝑚 = exp(𝜇 − 𝜎2)(8) It was assumed that 𝜎 is equal to 0.8, thus all the sampled parameters will have the same coefficient of variation (CV=0.947), which is calculated based on the following equation: 𝐶𝑉 = √exp(𝜎2)− 1 (9) (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint From the simulations, it observed that there is a large variability in the wild type T cells, with FasL exhibiting high expression in some cases and low expression in others. After knocking- out CISH, TNF-α and IL-2 increase, but also PD1, CTLA -4 and FasL, which is attributed to the rise of NFAT (Figure 4). This simultaneous rise in all markers can potentially explain the inability of CISH KO by itself to lead to tumor regression in mouse models of B16 melanoma[4], [28]. Therapy recommendations As noted in Palmer et al., CISH KO by itself was not sufficient to reduce tumor size in B16-melanoma bearing mice, but the combination of CISH KO and anti-PD1 block therapy had a high anti -tumoral efficacy [28]. Furthermore, Arthofer et al. reported an increase in tumor cell killing after inactivation of both CISH and PDCD1[67]. In addition, the global sensitivity analysis of the model shows the dephosphorylation rate of pGSK3 β to be important in some of these proteins. GSK3β is pivotal for deactivating NFAT via phosphorylation, since NFAT is inactivated when phosphorylated. To test this in the model, the production rates of CISH and PD1 mRNAs, v32 and v35 respectively, were set equal to zero to simulate simultaneous knockout of CISH and PDCD1 via multiplex gene editing[68]. In the case of simulated GSK3B knockout, both the initial conditions of GSK3β and pGSK3β were set to zero. Similarly to the previous section, 10000 cells for each case were simulated in a steady state. The synergistic inactivation of CISH and PDCD1 lead to a higher expression of IL -2, compared to CISH KO alone, while PDCD1 KO did not offer a significant improvement, as the IL -2 concentration distribution in this case is similar to the wild type cells (Figure 5A). When GSK3β is knocked-out in combination with CISH and PDCD1, the production of IL-2 is increased as evidenced by the higher IL -2 concentrations in the distribution (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint in this case. In addition to the increase of IL-2, the multiplex gene knockout of CISH, PDCD1 and GSK3B results in a better outcome for TNF -α, confirming the observation from the sensitivity analysis that TNF-α is the most sensitive to these changes (Figure 5B). However, FasL increases significantly in CISH KO and the combinations of CISH KO + PDCD1 KO or CISH KO + PDCD1 KO + GSK3B KO have a modest effect ( Figure 5C). CTLA-4 appears to increase significantly with the combination of CISH KO with PDCD1 KO and GSK3Β KO (Figure 5D). According to these results, CISH KO leads to the biggest increase in the biomarkers, compared to the individual gene knockouts, since it acts directly on PLC-γ1 and it has been shown that PLC-γ1 deficiency in T cell impairs their activation upon TCR stimulation [69]. CISH KO

in an increase in the transcription of IL-2 and TNF-α, with the latter being the most sensitive in this change. To improve therapeutic outcomes further, PDCD1 and GSK3B can also be knocked- out to enhance T cell activity. Because of the increase of FasL expression, which can lead to AICD when it binds to Fas receptor (FasR), a combination of CISH, PDCD1, GSK3B, CTLA4 and FAS knockout with multiplex gene-editing can be utilized to boost T cell reactivity, prevent exhaustion from PD1 and CTLA -4 and AICD from FasL. The benefits of targeting CTLA -4 and PD1 have been shown in numerous studies, including clinical trials[70], [71], [72], [73]. On the other hand, GSK3β is less popular as a potential target, despite its role in the inhibition of T cell proliferation as a result of the phosphorylation of NFAT [42], [43]. Previous studies demonstrated that GSK3β inhibition led to better activation of T cells by downregulating PD1 and LAG3 proteins, while increasing the expression of transcriptional regulator TBET[44], [74], [75], [76], [77], [78], [79] . These effects appear to balance the reduced T cell motility after the inhibition of GSK3β[80].

(which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Adoptive cell therapy , such as TIL therapy, is an emerging form of cancer therapy that leverages antitumoral lymphocytes to fight against cancer cells, resulting in tumor regression or even eradication[1], [7], [8]. Neoantigen TIL therapies manage to mitigate this barrier with higher cancer cell specificity, as has been reported in cases like metastatic bladder and breast cancers[81]. Despite the progress, there are many other limiting factors, such as the tumor microenvironment, which leads to suppressed proliferation and migration of T cells, low expression of neoantigens and the major histocompatibility complex (MHC) on cancer cells[3]. Traditional checkpoint inhibitors, such as anti -PD1 and anti -CTLA-4 antibodies, have shown beneficial therapeutic outcomes, but the high patient and tumor variability of PD1 and CTLA-4 expressions limit their efficacy [67], [70], [73], [82], [83], [84], [85], [86], [87] . CISH protein is an intracellular checkpoint that is responsible for the ubiquitination of PLC -γ1 and multiple studies have reported that CISH knockout enhances TIL reactivity against cancer cells. Palmer et al. found that CISH knockout alone does not lead to tumor regression, because of the increase in PD1, whereas a combination of CISH knockout and anti-PD1 therapy resulted in tumor regression in mice[28]. In this study, we developed an ODE signaling pathway model to capture the roles of CISH and NFAT in T cell activation. Our model was based on the model of Shin et al., and modified to capture the effect of CISH in the expression of functional output biomarkers[5], including IL-2, TNF-α, PD1, CTLA-4 and FasL. To estimate the parameters of the model, in vitro data from the literature have been used, where in all cases, T cells were stimulated with anti -CD28 and anti -CD3 antibodies to induce activation[5], [39], [53], [54], [55], [87]. There are significant differences in the transient behavior among the expression of proteins, for example NFAT and IL-2 exhibit a wider peak, compared to pAkt, pMEK, pERK, and pPLCγ1, which have a faster response followed by rapid decay. On the (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint other hand, CISH and TNF -α increase over time, but in the case of CISH, we only had four data points in early times, thus we cannot be confident that this is the actual behavior at longer times post-activation. Furthermore, this limitation applies to pPLC γ1 and PLC γ1, while CTLA -4 expression exhibits an initial peak and then settles to a steady state value. Despite the limited amount of data for some of the outputs, the parameter estimation mostly served the purpose of establishing the baseline values for t he protein degradation and generation rates, kdegi and kgeni respectively, as well as the Michaelis-Menten constants, kmi. To evaluate which parts of the pathway influence the outputs of the model the most, global sensitivity analysis was performed, using the Morris method[61]. The sensitivity analysis was used to find the most influential parameters with respect to the steady-state concentrations of IL-2, TNF- α, FasL, PD1 and CTLA-4, under constant stimulation of TCR. The estimated baseline parameters of kdegi, kgeni and kmi were kept constant during the global sensitivity analysis to gain a better understanding, since these parameters significantly influenced the outputs but do not contribute to improved insight into the biological mechanisms. From the analysis, it was conclud ed that CISH KO has a central role in the production of IL-2, TNF-α, CTLA-4, PD1 and FasL, and it potentially interacts with other parts of the pathway. An important observation is the influence of GSK3 β in expression of these proteins, indicating a possible additional target for multiplex gene -editing or drug targeting. GSK3β has been previously identified as a promising target for improving T cell activation, function and proliferation, despite the decrease in T cell migration[42], [43], [44], [80]. Furthermore, the steady-state simulations showed that CISH is the most crucial protein to target, but PD1 and GSK3β can also potentially be targeted to enhance TIL therapy effectiveness as well. In conclusion, CISH KO is effective in increasing both activation and inhibitory markers. To compensate for the increase in inhibitory gene expression, a multiplex gene-editing approach (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint is recommended where CISH is knocked-out in conjunction with knockout of CTLA-4, PD1, FAS and GSK3B to enhance T cell cytotoxicity and proliferation. The two most significant limitations of the current model are the limited experimental data that was used to calibrate the model and the focus on NFAT. A next step for improving the current model would be to collect more experimental time series data and include the NFκB transcription factor. ACKNOWLEDGMENTS B.R.W. acknowledges funding from Office of Discovery and Translation, NIH grants R21CA237789, R21AI163731, P01CA254849, P50CA136393, U54CA268069, R01AI146009, and Children’s Cancer Research Fund. B.S.M. acknowledges funding from NIH grants R01AI146009, R01AI161017, P01CA254849, P50CA136393, U24OD026641, U54CA232561, P30CA077598, U54CA268069, Children’s Cancer Research Fund, the Fanconi Anemia Research Fund, and the Randy Shaver Cancer and Community Fund. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Figures and Tables Figure 1. Model schematic for NFAT activation and downstream effects upon TCR activation. Black solid arrows and T -arrows indicate the activation and inhibition of proteins , respectively. Red arrows represent the expression of proteins that are dependent on NFAT transcription factor acting in the nucleus. Dashed black arrows indicate the translocation of proteins to the nucleus or in the extracellular space. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Figure 2. Model calibration based on 13 experimental outputs. The dots correspond to the experimental data obtained from references [5], [39], [52], [53], [54], [55] and the solid lines correspond to the model output. Equation (2) was used for the mapping between the model output and the experimental data. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Figure 3. Parameter sensitivity analysis using the Morris method for A) IL-2, B) TNF-α, C) FasL, D) PD1 and E) CTLA-4 steady state concentrations under constant activation of the TCR. The mean, 𝜇𝑖 ∗, is used to quantify the sensitivity of each output with respect to variations of a parameter, and the variance, 𝜎𝑖 2, for the nonlinear behavior or interactions with other parameters. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Table 1. Sensitive parameters based on Morris global sensitivity analysis method and their biological meaning. Parameter Meaning v33 Production of CISH protein v6a Deactivation of calcineurin v35 Production of mPD1 v36 Production of PD1 v27 Dephosphorylation of pGSK3β v6 Production of mNFAT v24 Dephosphorylation of ERK v20a Production of FasL v16 Production of mCTLA-4 v17 Production of CTLA-4 v31 Phosphorylation of inactive IL2R Figure 4. In silico steady state simulations for wild type and CISH KO T cells. CISH KO increases the expression of all markers. CISH KO increases the expression of IL -2, TNF-α, PD1, FasL and CTLA-4 simultaneously. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Figure 5. Steady-state concentration distributions of A) IL-2, B) TNF-α, C) FasL and D) CTLA- 4 for different therapy combinations. CISH KO + PD1 KO + GSK3 β KO leads to an increase in IL-2, TNF-α and FasL. CTLA-4 only requires CISH knockout to have higher expression, whereas the other outputs are further increased at steady -state by combining CISH KO with one or more other gene knockouts, such as PD1 and/or GSK3β. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint

[1] S. A. Rosenberg and N. P. Restifo, “Adoptive cell transfer as personalized immunotherapy for human cancer,” Science (1979), vol. 348, no. 6230, pp. 62–68, Apr. 2015, doi: 10.1126/SCIENCE.AAA4967/ASSET/25CD0751-5C0D-4AE1-AB3E- 242082A70CE0/ASSETS/GRAPHIC/348_62_F4.JPEG. [2] D. A. Morgan, F. W. Ruscetti, and R. Gallo, “Selective in vitro growth of T lymphocytes from normal human bone marrows,” Science, vol. 193, no. 4257, pp. 1007–1008, 1976, doi: 10.1126/SCIENCE.181845. [3] K. G. Anderson, I. M. Stromnes, and P. D. Greenberg, “Obstacles Posed by the Tumor Microenvironment to T cell Activity: A Case for Synergistic Therapies,” Cancer Cell, vol. 31, no. 3, pp. 311–325, Mar. 2017, doi: 10.1016/J.CCELL.2017.02.008/ASSET/0A6F3E9F- 354F-4864-B90D-DD8E7617B493/MAIN.ASSETS/GR2.JPG. [4] E. Lou et al., “Targeting the intracellular immune checkpoint CISH with CRISPR-Cas9- edited T cells in patients with metastatic colorectal cancer: a first-in-human, single- centre, phase 1 trial,” Lancet Oncol, vol. 26, no. 5, pp. 559–570, May 2025, doi: 10.1016/S1470-2045(25)00083-X. [5] S. Y. Shin, M. W. Kim, K. H. Cho, and L. K. Nguyen, “Coupled feedback regulation of nuclear factor of activated T-cells (NFAT) modulates activation-induced cell death of T cells,” Sci Rep, vol. 9, no. 1, p. 10637, Dec. 2019, doi: 10.1038/S41598-019-46592-Z. [6] C. Bonini and A. Mondino, “Adoptive T-cell therapy for cancer: The era of engineered T cells,” Eur J Immunol, vol. 45, no. 9, pp. 2457–2469, Sep. 2015, doi: 10.1002/EJI.201545552. [7] S. A. Rosenberg, J. J. Mulé, P. J. Spiess, C. M. Reichert, and S. L. Schwarz, “Regression of established pulmonary metastases and subcutaneous tumor mediated by the systemic administration of high-dose recombinant interleukin 2,” J Exp Med, vol. 161, no. 5, pp. 1169–1188, May 1985, doi: 10.1084/JEM.161.5.1169. [8] J. H. Donohue, M. Rosenstein, A. E. Chang, M. T. Lotze, R. J. Robb, and S. A. Rosenberg, “The systemic administration of purified interleukin 2 enhances the ability of sensitized murine lymphocytes to cure a disseminated syngeneic lymphoma.,” The Journal of Immunology, vol. 132, no. 4, pp. 2123–2128, Apr. 1984, doi: 10.4049/JIMMUNOL.132.4.2123. [9] S. A. Rosenberg et al., “Observations on the systemic administration of autologous lymphokine-activated killer cells and recombinant interleukin-2 to patients with metastatic cancer,” N Engl J Med, vol. 313, no. 23, pp. 1485–1492, Dec. 1985, doi: 10.1056/NEJM198512053132327. [10] S. A. Rosenberg et al., “Prospective randomized trial of high-dose interleukin-2 alone or in conjunction with lymphokine-activated killer cells for the treatment of patients with (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint advanced cancer,” J Natl Cancer Inst, vol. 85, no. 8, pp. 622–632, Apr. 1993, doi: 10.1093/JNCI/85.8.622. [11] S. A. Rosenberg, P. Spiess, and R. Lafreniere, “A New Approach to the Adoptive Immunotherapy of Cancer with Tumor-Infiltrating Lymphocytes,” Science (1979), vol. 233, no. 4770, pp. 1318–1321, 1986, doi: 10.1126/SCIENCE.3489291. [12] S. A. Rosenberg et al., “Use of tumor-infiltrating lymphocytes and interleukin-2 in the immunotherapy of patients with metastatic melanoma. A preliminary report,” N Engl J Med, vol. 319, no. 25, pp. 1676–1680, Dec. 1988, doi: 10.1056/NEJM198812223192527. [13] M. E. Dudley et al., “Cancer Regression and Autoimmunity in Patients After Clonal Repopulation with Antitumor Lymphocytes,” Science, vol. 298, no. 5594, p. 850, Oct. 2002, doi: 10.1126/SCIENCE.1076514. [14] P. F. Robbins et al., “Mining exomic sequencing data to identify mutated antigens recognized by adoptively transferred tumor-reactive T cells,” Nat Med, vol. 19, no. 6, pp. 747–752, Jun. 2013, doi: 10.1038/NM.3161,. [15] N. Zacharakis et al., “Immune recognition of somatic mutations leading to complete durable regression in metastatic breast cancer,” Nature Medicine 2018 24:6, vol. 24, no. 6, pp. 724–730, Jun. 2018, doi: 10.1038/s41591-018-0040-8. [16] E. Tran et al., “T-Cell Transfer Therapy Targeting Mutant KRAS in Cancer,” New England Journal of Medicine, vol. 375, no. 23, pp. 2255–2262, Dec. 2016, doi: 10.1056/NEJMOA1609279/SUPPL_FILE/NEJMOA1609279_DISCLOSURES.PDF. [17] J. Chesney et al., “Efficacy and safety of lifileucel, a one-time autologous tumor- infiltrating lymphocyte (TIL) cell therapy, in patients with advanced melanoma after progression on immune checkpoint inhibitors and targeted therapies: pooled analysis of consecutive cohorts of the C-144-01 study,” J Immunother Cancer, vol. 10, no. 12, p. e005755, Dec. 2022, doi: 10.1136/JITC-2022-005755. [18] N. Zacharakis et al., “Breast Cancers Are Immunogenic: Immunologic Analyses and a Phase II Pilot Clinical Trial Using Mutation-Reactive Autologous Lymphocytes,” J Clin Oncol, vol. 40, no. 16, pp. 1741–1754, Jun. 2022, doi: 10.1200/JCO.21.02170. [19] D. C. Palmer et al., “Cish actively silences TCR signaling in CD8+ T cells to maintain tumor tolerance,” Journal of Experimental Medicine, vol. 212, no. 12, pp. 2095–2113, Nov. 2015, doi: 10.1084/JEM.20150304. [20] D. J. Hilton et al., “Twenty proteins containing a C-terminal SOCS box form five structural classes,” Proc Natl Acad Sci U S A, vol. 95, no. 1, pp. 114–119, Jan. 1998, doi: 10.1073/PNAS.95.1.114. [21] R. Starr et al., “A family of cytokine-inducible inhibitors of signalling,” Nature 1997 387:6636, vol. 387, no. 6636, pp. 917–921, Jun. 1997, doi: 10.1038/43206. [22] T. A. Endo et al., “A new protein containing an SH2 domain that inhibits JAK kinases,” Nature 1997 387:6636, vol. 387, no. 6636, pp. 921–924, Jun. 1997, doi: 10.1038/43213. [23] T. Naka et al., “Structure and function of a new STAT-induced STAT inhibitor,” Nature 1997 387:6636, vol. 387, no. 6636, pp. 924–929, Jun. 1997, doi: 10.1038/43219. [24] R. B. Delconte et al., “CIS is a potent checkpoint in NK cell-mediated tumor immunity,” Nat Immunol, vol. 17, no. 7, pp. 816–824, Jun. 2016, doi: 10.1038/NI.3470. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint [25] R. Morris, N. J. Kershaw, and J. J. Babon, “The molecular details of cytokine signaling via the JAK/STAT pathway,” Protein Science, vol. 27, no. 12, pp. 1984–2009, Dec. 2018, doi: 10.1002/PRO.3519. [26] F. Verdier et al., “Proteasomes Regulate Erythropoietin Receptor and Signal Transducer and Activator of Transcription 5 (STAT5) Activation: POSSIBLE INVOLVEMENT OF THE UBIQUITINATED CIS PROTEIN,” Journal of Biological Chemistry, vol. 273, no. 43, pp. 28185–28190, Oct. 1998, doi: 10.1074/JBC.273.43.28185. [27] A. Matsumoto et al., “Suppression of STAT5 functions in liver, mammary glands, and T cells in cytokine-inducible SH2-containing protein 1 transgenic mice,” Mol Cell Biol, vol. 19, no. 9, pp. 6396–6407, Sep. 1999, doi: 10.1128/MCB.19.9.6396. [28] D. C. Palmer et al., “Internal checkpoint regulates T cell neoantigen reactivity and susceptibility to PD1 blockade,” Med (N Y), vol. 3, no. 10, pp. 682-704.e8, Oct. 2022, doi: 10.1016/J.MEDJ.2022.07.008. [29] A. K. Mehta, D. T. Gracias, and M. Croft, “TNF Activity and T cells,” Cytokine, vol. 101, p. 14, Jan. 2016, doi: 10.1016/J.CYTO.2016.08.003. [30] X. Chen, M. Bäumel, D. N. Männel, O. M. Z. Howard, and J. J. Oppenheim, “Interaction of TNF with TNF Receptor Type 2 Promotes Expansion and Function of Mouse CD4+CD25+ T Regulatory Cells,” The Journal of Immunology, vol. 179, no. 1, pp. 154–161, Jul. 2007, doi: 10.4049/JIMMUNOL.179.1.154. [31] V. Niederlova, O. Tsyklauri, M. Kovar, and O. Stepanek, “IL-2-driven CD8+ T cell phenotypes: implications for immunotherapy,” Trends Immunol, vol. 44, no. 11, pp. 890– 901, Nov. 2023, doi: 10.1016/J.IT.2023.09.003/ASSET/94B18407-EFBC-4066-B52A- 82ED670C9C48/MAIN.ASSETS/GR4.JPG. [32] S. H. Ross and D. A. Cantrell, “Signaling and Function of Interleukin-2 in T Lymphocytes,” Annu Rev Immunol, vol. 36, p. 411, Apr. 2018, doi: 10.1146/ANNUREV-IMMUNOL- 042617-053352. [33] L. B. Silva et al., “The Role of TNF-α as a Proinflammatory Cytokine in Pathological Processes,” Open Dent J, vol. 13, no. 1, pp. 332–338, Oct. 2019, doi: 10.2174/1874210601913010332. [34] B. L. Salomon, M. Leclerc, J. Tosello, E. Ronin, E. Piaggio, and J. L. Cohen, “Tumor necrosis factor α and regulatory T cells in oncoimmunology,” Front Immunol, vol. 9, no. MAR, p. 345465, Mar. 2018, doi: 10.3389/FIMMU.2018.00444/PDF. [35] K. S. Rallis, C. R. T. Hillyar, M. Sideris, and J. K. Davies, “T-cell-based Immunotherapies for Haematological Cancers, Part A: A SWOT Analysis of Immune Checkpoint Inhibitors (ICIs) and Bispecific T-Cell Engagers (BiTEs),” Anticancer Res, vol. 41, no. 3, pp. 1123–1141, Mar. 2021, doi: 10.21873/ANTICANRES.14870. [36] V. A. Boussiotis, “Molecular and Biochemical Aspects of the PD-1 Checkpoint Pathway,” N Engl J Med, vol. 375, no. 18, pp. 1767–1778, Nov. 2016, doi: 10.1056/NEJMRA1514296. [37] L. Chen and D. B. Flies, “Molecular mechanisms of T cell co-stimulation and co- inhibition,” Nat Rev Immunol, vol. 13, no. 4, pp. 227–242, Apr. 2013, doi: 10.1038/NRI3405. [38] X. Xu et al., “CTLA4 depletes T cell endogenous and trogocytosed B7 ligands via cis- endocytosis,” J Exp Med, vol. 220, no. 7, p. e20221391, Jul. 2023, doi: 10.1084/JEM.20221391. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint [39] H. M. Gibson et al., “Induction of the CTLA-4 gene in human lymphocytes is dependent on NFAT binding the proximal promoter,” J Immunol, vol. 179, no. 6, pp. 3831–3840, Sep. 2007, doi: 10.4049/JIMMUNOL.179.6.3831. [40] E. A. Wohlfert and R. B. Clark, “’Vive la Résistance!’--the PI3K-Akt pathway can determine target sensitivity to regulatory T cell suppression,” Trends Immunol, vol. 28, no. 4, pp. 154–160, Apr. 2007, doi: 10.1016/J.IT.2007.02.003. [41] E. Volpe, M. Sambucci, L. Battistini, and G. Borsellino, “Fas–Fas Ligand: Checkpoint of T Cell Functions in Multiple Sclerosis,” Front Immunol, vol. 7, no. SEP, p. 382, Sep. 2016, doi: 10.3389/FIMMU.2016.00382. [42] C. E. Rudd, “Upstream-downstream: CD28 co-signaling pathways and T cell function,” Immunity, vol. 4, no. 6, pp. 527–534, Jun. 1996, doi: 10.1016/S1074-7613(00)80479- 3/ASSET/0BF4FBBD-43D6-470D-BD13-99F95091B6D1/MAIN.ASSETS/GR4.JPG. [43] T. Ohteki et al., “Negative Regulation of T Cell Proliferation and Interleukin 2 Production by the Serine Threonine Kinase Gsk-3,” Journal of Experimental Medicine, vol. 192, no. 1, pp. 99–104, Jul. 2000, doi: 10.1084/JEM.192.1.99. [44] A. R. Crespo, S. G. Luna, B. Moës, A. Rodriguez, and C. E. Rudd, “The many faceted role of glycogen synthase kinase-3 (GSK-3) in T cells and cancer immunotherapy,” J Cancer Biol, vol. 5, no. 1, p. 11, Apr. 2024, doi: 10.46439/CANCERBIOLOGY.5.058. [45] J. Bezanson, A. Edelman, S. Karpinski, and V. B. Shah, “Julia: A Fresh Approach to Numerical Computing,” https://doi.org/10.1137/141000671, vol. 59, no. 1, pp. 65–98, Feb. 2017, doi: 10.1137/141000671. [46] C. Rackauckas and Q. Nie, “DifferentialEquations.jl – A Performant and Feature-Rich Ecosystem for Solving Differential Equations in Julia,” J Open Res Softw, vol. 5, no. 1, p. 15, May 2017, doi: 10.5334/JORS.151. [47] A. Sandu, J. G. Verwer, J. G. Blom, E. J. Spee, G. R. Carmichael, and F. A. Potra, “Benchmarking stiff ode solvers for atmospheric chemistry problems II: Rosenbrock solvers,” Atmos Environ, vol. 31, no. 20, pp. 3459–3472, Oct. 1997, doi: 10.1016/S1352- 2310(97)83212-8. [48] L. F. Shampine and M. W. Reichelt, “The MATLAB ODE Suite,” https://doi.org/10.1137/S1064827594276424, vol. 18, no. 1, pp. 1–22, Jul. 2006, doi: 10.1137/S1064827594276424. [49] C. Rackauckas and Q. Nie, “Confederated modular differential equation APIs for accelerated algorithm development and benchmarking,” Advances in Engineering Software, vol. 132, pp. 1–6, Jun. 2019, doi: 10.1016/J.ADVENGSOFT.2019.03.009. [50] C. Tsitouras, “Runge–Kutta pairs of order 5(4) satisfying only the first column simplifying assumption,” Computers & Mathematics with Applications, vol. 62, no. 2, pp. 770–775, Jul. 2011, doi: 10.1016/J.CAMWA.2011.06.002. [51] R. D. Neidinger, “Introduction to Automatic Differentiation and MATLAB Object-Oriented Programming,” https://doi.org/10.1137/080743627, vol. 52, no. 3, pp. 545–563, Aug. 2010, doi: 10.1137/080743627. [52] R. D. Reitsema et al., “Effect of age and sex on immune checkpoint expression and kinetics in human T cells,” Immunity and Ageing, vol. 17, no. 1, pp. 1–12, Dec. 2020, doi: 10.1186/S12979-020-00203-Y/FIGURES/7. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint [53] P. K. Chattopadhyay, J. Yu, and M. Roederer, “A live-cell assay to detect antigen-specific CD4+ T cells with diverse cytokine profiles,” Nature Medicine 2005 11:10, vol. 11, no. 10, pp. 1113–1117, Sep. 2005, doi: 10.1038/nm1293. [54] D. Lodygin et al., “A combination of fluorescent NFAT and H2B sensors uncovers dynamics of T cell activation in real time during CNS autoimmunity,” Nat Med, vol. 19, no. 6, pp. 784–790, Jun. 2013, doi: 10.1038/NM.3182. [55] M. E. Kalland, N. G. Oberprieler, T. Vang, K. Taskén, and K. M. Torgersen, “T cell-signaling network analysis reveals distinct differences between CD28 and CD2 costimulation responses in various subsets and in the MAPK pathway between resting and activated regulatory T cells,” J Immunol, vol. 187, no. 10, pp. 5233–5245, Nov. 2011, doi: 10.4049/JIMMUNOL.1101804. [56] R. N. Gutenkunst, J. J. Waterfall, F. P. Casey, K. S. Brown, C. R. Myers, and J. P. Sethna, “Universally Sloppy Parameter Sensitivities in Systems Biology Models,” PLoS Comput Biol, vol. 3, no. 10, p. e189, 2007, doi: 10.1371/JOURNAL.PCBI.0030189. [57] F. Fröhlich and P. K. Sorger, “Fides: Reliable trust-region optimization for parameter estimation of ordinary differential equation models,” PLoS Comput Biol, vol. 18, no. 7, p. e1010322, Jul. 2022, doi: 10.1371/JOURNAL.PCBI.1010322. [58] L. Schmiester, Y. Schälte, F. Fröhlich, J. Hasenauer, and D. Weindl, “Efficient parameterization of large-scale dynamic models based on relative measurements,” Bioinformatics, vol. 36, no. 2, pp. 594–602, Jan. 2020, doi: 10.1093/BIOINFORMATICS/BTZ581. [59] Å. Björck, “Numerical Methods for Least Squares Problems,” Numerical Methods for Least Squares Problems, Jan. 1996, doi: 10.1137/1.9781611971484. [60] L. Schmiester et al., “PEtab—Interoperable specification of parameter estimation problems in systems biology,” PLoS Comput Biol, vol. 17, no. 1, p. e1008646, Jan. 2021, doi: 10.1371/JOURNAL.PCBI.1008646. [61] M. D. Morris, “Factorial sampling plans for preliminary computational experiments,” Technometrics, vol. 33, no. 2, pp. 161–174, 1991, doi: 10.1080/00401706.1991.10484804. [62] A. Saltelli et al., “Global sensitivity analysis: The primer,” Global Sensitivity Analysis: The Primer, pp. 1–292, Jan. 2008, doi: 10.1002/9780470725184. [63] J. Schindelin et al., “Fiji: an open-source platform for biological-image analysis,” Nature

2012 9:7, vol. 9, no. 7, pp. 676–682, Jun. 2012, doi: 10.1038/nmeth.2019. [64] “WebPlotDigitizer.” Accessed: Feb. 24, 2025. [Online]. Available: https://automeris.io/docs/#cite-webplotdigitizer [65] N. E. Scharping et al., “Mitochondrial stress induced by continuous stimulation under hypoxia rapidly drives T cell exhaustion,” Nat Immunol, vol. 22, no. 2, pp. 205–215, Feb. 2021, doi: 10.1038/S41590-020-00834-9,. [66] D. R. Heintzman, E. L. Fisher, and J. C. Rathmell, “Microenvironmental influences on T cell immunity in cancer and inflammation,” Cellular & Molecular Immunology 2022 19:3, vol. 19, no. 3, pp. 316–326, Jan. 2022, doi: 10.1038/s41423-021-00833-2. [67] E. Arthofer et al., “Genetic editing of CISH enhances T cell effector programs independently of immune checkpoint cell surface ligand expression,” bioRxiv, p. 2021.08.17.456714, Aug. 2021, doi: 10.1101/2021.08.17.456714. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint [68] B. R. Webber et al., “Highly efficient multiplex human T cell engineering without double- strand breaks using Cas9 base editors,” Nature Communications 2019 10:1, vol. 10, no. 1, pp. 1–10, Nov. 2019, doi: 10.1038/s41467-019-13007-6. [69] G. Fu et al., “Phospholipase Cγ1 is essential for T cell development, activation, and tolerance,” J Exp Med, vol. 207, no. 2, p. 309, Feb. 2010, doi: 10.1084/JEM.20090880. [70] S. Mok, D. A. Çobanoğlu, H. Liu, J. J. Mancuso, and J. P. Allison, “Post-immunotherapy CTLA-4 Ig treatment improves antitumor efficacy,” Proc Natl Acad Sci U S A, vol. 121, no. 27, p. e2404661121, Jul. 2024, doi: 10.1073/PNAS.2404661121/SUPPL_FILE/PNAS.2404661121.SAPP.PDF. [71] A. Parvez et al., “PD-1 and PD-L1: architects of immune symphony and immunotherapy breakthroughs in cancer treatment,” Front Immunol, vol. 14, p. 1296341, Dec. 2023, doi: 10.3389/FIMMU.2023.1296341/PDF. [72] Y. Jiang, M. Chen, H. Nie, and Y. Yuan, “PD-1 and PD-L1 in cancer immunotherapy: clinical implications and future considerations,” Hum Vaccin Immunother, vol. 15, no. 5, p. 1111, May 2019, doi: 10.1080/21645515.2019.1571892. [73] P. De Silva, M. Aiello, C. Gu-Trantien, E. Migliori, K. Willard-Gallo, and C. Solinas, “Targeting CTLA-4 in cancer: Is it the ideal companion for PD-1 blockade immunotherapy combinations?,” Int J Cancer, vol. 149, no. 1, pp. 31–41, Jul. 2021, doi: 10.1002/IJC.33415. [74] A. Taylor, J. A. Harker, K. Chanthong, P. G. Stevenson, E. I. Zuniga, and C. E. Rudd, “Glycogen Synthase Kinase 3 Inactivation Drives T-bet-Mediated Downregulation of Co- receptor PD-1 to Enhance CD8(+) Cytolytic T Cell Responses,” Immunity, vol. 44, no. 2, pp. 274–286, Feb. 2016, doi: 10.1016/J.IMMUNI.2016.01.018. [75] A. Taylor and C. E. Rudd, “Glycogen Synthase Kinase 3 Inactivation Compensates for the Lack of CD28 in the Priming of CD8+ Cytotoxic T-Cells: Implications for anti-PD-1 Immunotherapy,” Front Immunol, vol. 8, no. DEC, Dec. 2017, doi: 10.3389/FIMMU.2017.01653. [76] J. Krueger, C. E. Rudd, and A. Taylor, “Glycogen synthase 3 (GSK-3) regulation of PD-1 expression and and its therapeutic implications,” Semin Immunol, vol. 42, Apr. 2019, doi: 10.1016/J.SMIM.2019.101295. [77] A. Taylor, D. Rothstein, and C. E. Rudd, “Small-Molecule Inhibition of PD-1 Transcription Is an Effective Alternative to Antibody Blockade in Cancer Therapy,” Cancer Res, vol. 78, no. 3, pp. 706–717, Feb. 2018, doi: 10.1158/0008-5472.CAN-17-0491. [78] A. Taylor, J. A. Harker, K. Chanthong, P. G. Stevenson, E. I. Zuniga, and C. E. Rudd, “Glycogen Synthase Kinase 3 Inactivation Drives T-bet-Mediated Downregulation of Co- receptor PD-1 to Enhance CD8(+) Cytolytic T Cell Responses,” Immunity, vol. 44, no. 2, pp. 274–286, Feb. 2016, doi: 10.1016/J.IMMUNI.2016.01.018. [79] C. E. Rudd, K. Chanthong, and A. Taylor, “Small Molecule Inhibition of GSK-3 Specifically Inhibits the Transcription of Inhibitory Co-receptor LAG-3 for Enhanced Anti-tumor Immunity,” Cell Rep, vol. 30, no. 7, pp. 2075-2082.e4, Feb. 2020, doi: 10.1016/J.CELREP.2020.01.076. [80] A. Taylor and C. E. Rudd, “Glycogen synthase kinase 3 (GSK-3) controls T-cell motility and interactions with antigen presenting cells,” BMC Res Notes, vol. 13, no. 1, p. 163, Mar. 2020, doi: 10.1186/S13104-020-04971-0. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint [81] Z. Zhang et al., “Neoantigen: A New Breakthrough in Tumor Immunotherapy,” Front Immunol, vol. 12, p. 672356, Apr. 2021, doi: 10.3389/FIMMU.2021.672356/BIBTEX. [82] N. Krishnamurthy et al., “High CTLA-4 transcriptomic expression correlates with high expression of other checkpoints and with immunotherapy outcome,” Ther Adv Med Oncol, vol. 16, p. 17588359231220510, Jan. 2024, doi: 10.1177/17588359231220510. [83] Y. Gonzalez-Montes et al., “Genetic variants of CTLA4 are associated with clinical outcome of patients with multiple myeloma,” Front Immunol, vol. 14, p. 1158105, 2023, doi: 10.3389/FIMMU.2023.1158105/FULL. [84] A. A. Davis and V. G. Patel, “The role of PD-L1 expression as a predictive biomarker: an analysis of all US Food and Drug Administration (FDA) approvals of immune checkpoint inhibitors,” J Immunother Cancer, vol. 7, no. 1, p. 278, Dec. 2019, doi: 10.1186/S40425- 019-0768-9. [85] K. Pawelczyk et al., “Role of PD-L1 Expression in Non-Small Cell Lung Cancer and Their Prognostic Significance according to Clinicopathological Factors and Diagnostic Markers,” Int J Mol Sci, vol. 20, no. 4, Feb. 2019, doi: 10.3390/IJMS20040824. [86] J. Lee et al., “Outstanding clinical efficacy of PD-1/PD-L1 inhibitors for pulmonary pleomorphic carcinoma,” Eur J Cancer, vol. 132, pp. 150–158, Jun. 2020, doi: 10.1016/J.EJCA.2020.03.029. [87] A. Martín-Ruiz et al., “Effects of anti-PD-1 immunotherapy on tumor regression: insights from a patient-derived xenograft model,” Sci Rep, vol. 10, no. 1, Dec. 2020, doi: 10.1038/S41598-020-63796-W. [88] B. Jones, D. Stekel, J. Rowe, and C. Fernando, “Is there a liquid state machine in the bacterium Escherichia coli?,” Proceedings of the 2007 IEEE Symposium on Artificial Life, CI-ALife 2007, pp. 187–191, 2007, doi: 10.1109/ALIFE.2007.367795. [89] H. Chen, K. Shiroguchi, H. Ge, and X. S. Xie, “Genome-wide study of mRNA degradation and transcript elongation in Escherichia coli,” Mol Syst Biol, vol. 11, no. 1, p. 781, Jan. 2015, doi: 10.15252/MSB.20145794. [90] S. Kühner et al., “Proteome organization in a genome-reduced bacterium,” Science (1979), vol. 326, no. 5957, pp. 1235–1240, Nov. 2009, doi: 10.1126/SCIENCE.1176343/SUPPL_FILE/KUEHNER.SOM.PDF. [91] J.-M. Anaya, Y. Shoenfeld, A. Rojas-Villarraga, R. A. Levy, and R. Cervera, “Autoimmunity,” Autoimmunity: From Bench to Bedside, pp. 237–378, Jul. 2013, Accessed: Jun. 11, 2025. [Online]. Available: https://www.ncbi.nlm.nih.gov/books/NBK459447/ [92] D. Davidia et al., “Global characterization of in vivo enzyme catalytic rates and their correspondence to in vitro kcat measurements,” Proc Natl Acad Sci U S A, vol. 113, no. 12, pp. 3401–3406, Mar. 2016, doi: 10.1073/PNAS.1514240113/SUPPL_FILE/PNAS.1514240113.SD01.XLSX. [93] Y. Wang, C. L. Liu, J. D. Storey, R. J. Tibshirani, D. Herschlag, and P. O. Brown, “Precision and functional specificity in mRNA decay,” Proc Natl Acad Sci U S A, vol. 99, no. 9, pp. 5860–5865, Apr. 2002, doi: 10.1073/PNAS.092538799,. [94] S. B. Cambridge, F. Gnad, C. Nguyen, J. L. Bermejo, M. Krüger, and M. Mann, “Systems- wide proteomic analysis in mammalian cells reveals conserved, functional protein (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint turnover,” J Proteome Res, vol. 10, no. 12, pp. 5275–5284, Dec. 2011, doi: 10.1021/PR101183K/SUPPL_FILE/PR101183K_SI_002.XLS. Supplementary Information To estimate the baseline generation rates for mRNAs, we assumed that the average mRNA length is 1000 nucleotides, the maximum elongation rate 70 nt/s and the average number of RNA polymerases[88], [89], [90] in a T-cell is 300. Moreover, it was assumed that the diameter of a T- cell is 10 μm, thus its volume is equal to 4 3 𝜋𝑟3 = 523.6 μm3[91]. Figure A. 1. Waterfall plot for 10000 parameter estimations with different initial parameter values. The ranking is from the best on the left to the worst on the right on the x -axis. The y -axis corresponds to the negative log-likelihood, calculated based on (1). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Table A. 1. Table of reaction rates for CISH model. No Reaction Reaction Rates re_1 TCR → pTCR kc1*E1*TCR/(km1+TCR) re_2 pTCR + CTLA4 → TCR + CTLA4 kc2*(CTLA4)*pTCR/(km2+pTCR) re_4 pTCR → iTCR v3*pTCR/(km3+pTCR) re_5 iTCR → TCR v4*iTCR/(km4+iTCR) re_6 CN → aCN kc5a*pPLCg*CN/(km5a+CN) re_7 aCN → CN v6a*aCN/(km6a+aCN) re_8 aCN+RCAN  aCN_RCAN kon*aCN*RCAN – kon*KD4*aCN_RCAN re_9 aCN+pRCAN  aCN_pRCAN kon*aCN*pRCAN – kon*KD5*aCN_pRCAN re_10 aCN+Carabin  aCN_Carabin kon*aCN*Carabin – kon*KD6*aCN_Carabin re_11 pNFAT + aCN  pNFAT_aCN kon* pNFAT*aCN – kon*KD7*pNFAT_aCN re_12 pNFAT_aCN → NFAT + aCN kp2*pNFAT_aCN re_13 pNFAT + aCN_pRCAN  pNFAT_aCN_pRCAN kon* pNFAT *aCN_pRCAN – kon*KD8* pNFAT_aCN_pRCAN re_14 pNFAT_aCN_pRCAN → NFAT + aCN_pRCAN kp3*pNFAT_aCN_pRCAN re_15 NFAT → pNFAT v5*NFAT/(km5 + NFAT) re_16 NFAT + GSK3  NFAT_GSK3 kon*pNFAT*GSK3 – kon*KD9*NFAT_GSK3 re_17 NFAT_GSK3 → pNFAT + GSK3 kp4*NFAT_GSK3 re_18  → mNFAT kgen1 re_19  → mNFAT v6*NFAT/(km6+NFAT) re_20 mNFAT →  kdeg1*mNFAT re_21  → NFAT v7*mNFAT/(km7+mNFAT) re_22 NFAT →  kdeg2*NFAT re_23 pNFAT →  kdeg3*pNFAT re_24  → mRCAN kgen2 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint re_25  → mRCAN v8*NFAT^2/(km8^2+NFAT^2) re_26 mRCAN →  kdeg4*mRCAN re_27  → RCAN v9*mRCAN/(km9+mRCAN) re_28 RCAN →  kdeg5*RCAN re_29 RCAN + pTAK1  RCAN_pTAK1 kon*RCAN*pTAK1 – kon*KD10* RCAN_pTAK1 re_30 RCAN_pTAK1 → pRCAN + pTAK1 kp5* RCAN_pTAK1 re_31 pRCAN → RCAN v10*pRCAN/(km10+pRCAN) re_32  → mCarabin kgen3 re_33  → mCarabin v11*NFAT^2/(km11^2+NFAT^2) re_34 mCarabin →  kdeg6*mCarabin re_35  → Carabin v12*mCarabin/(km12+mCarabin) re_36 Carabin →  kdeg7*Carabin re_37  → mTNFa kgen4 re_38  → mTNFa v13*NFAT^2/(km13^2+NFAT^2) re_39 mTNFa →  kdeg8*mTNFa re_40  → TNFa v14*mTNFa/(km14+mTNFa) re_41 TNFa →  kdeg9*TNFa re_42 TAK1 → pTAK1 v15*pTAK1/(km15+pTAK1) re_43 pTAK1 → TAK1 kc38*TAK1*TNFa/(km38 + TAK1) re_44  → mCTLA4 kgen5 re_45  → mCTLA4 v16*NFAT^2/(km16^2+NFAT^2) re_46 mCTLA4 →  kdeg10*mCTLA4 re_47  → CTLA4 v17*mCTLA4/(km17+mCTLA4) re_48 CTLA4 →  kdeg11*CTLA4 re_49  → mIL2 kgen6 re_50  → mIL2 v18a*(1+v18b*pERK/(km18b+pERK))*NFAT^2/(km1 8^2+NFAT^2) re_51 mIL2 →  kdeg12*mIL2 re_52  → IL2 v19*mIL2/(km19+mIL2) re_53 IL2 →  kdeg13*IL2 re_54  → mFasL kgen7 re_55  → mFasL v20a*(1+v20b*pERK/(km20b+pERK))*NFAT^5/(km2 0^5+NFAT^5) re_56 mFasL →  kdeg14*mFasL re_57  → FasL v21*mFasL/(km21+mFasL) re_58 FasL →  kdeg15*FasL re_59 Ras → aRas kc39*Ras*pPLCg/(km39 + Ras) re_60 aRas → Ras v22*aRas/(km22+aRas) re_61 Carabin + aRas  aRas_Carabin kon*aRas*Carabin – KD13*kon*aRas_Carabin (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint re_62 MEK + aRas  aRas_MEK kon*aRas*MEK – KD14*kon*aRas_MEK re_63 aRas_MEK → pMEK + aRas kp8*aRas_MEK re_64 MEK + Rac1GTP  Rac1GTP_MEK kon*Rac1GTP *MEK – KD15*kon*Rac1GTP_MEK re_65 Rac1GTP_MEK → pMEK + Rac1GTP kp9*Rac1GTP_MEK re_66 pMEK → MEK v23*pMEK/(km23+pMEK) re_67 ERK + pMEK  pMEK_ERK kon* ERK*pMEK – KD16*kon*pMEK_ERK re_68 pMEK_ERK → pMEK + pERK kp10*pMEK_ERK re_69 pERK → ERK v24*pERK/(km24+pERK) re_70 PI3K + pIL2R  PI3K_pIL2R kon*PI3K*pIL2R – KD17*kon*PI3K_pIL2R re_71 PI3K_pIL2R → aPI3K + pIL2R kp11*PI3K_pIL2R re_72 PI3K + pTCR  PI3K_pTCR kon*PI3K * pTCR – KD18*kon*PI3K_pTCR re_73 PI3K_pTCR → aPI3K + pTCR kp12*PI3K_pTCR re_74 aPI3K → PI3K v25*aPI3K/(km25+aPI3K) re_75 PP2A + CTLA4  PP2A_CTLA4 kon*PP2A*CTLA4 – KD19*kon*PP2A_CTLA4 re_76 PP2A_CTLA4 → aPP2A + CTLA4 kp13*PP2A_CTLA4 re_77 aPP2A → PP2A v26*aPP2A/(km26+aPP2A) re_78 Akt + aPI3K  Akt_aPI3K kon*Akt*aPI3K – KD20*kon*Akt_aPI3K re_79 Akt_aPI3K → pAkt + aPI3K kp14*Akt_aPI3K re_80 pAkt + aPP2A  pAkt_aPP2A kon*pAkt*aPP2A – KD21*kon*pAkt_aPP2A re_81 pAkt_aPP2A → Akt + aPP2A kp15*pAkt_aPP2A re_82 GSK3 + pAkt  GSK3_pAkt kon* pAkt *GSK3 – KD22*kon*GSK3_pAkt re_83 GSK3_pAkt → pGSK3 + pAkt kp16*GSK3_pAkt re_84 pGSK3 → GSK3 v27*pGSK3/(km27+pGSK3) re_85 Rac1GDP + aPI3K  Rac1GDP_aPI3K kon*Rac1GDP*aPI3K – KD23*kon* Rac1GDP_aPI3K (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint re_86 Rac1GDP_aPI3K → Rac1GTP + aPI3K kp17*Rac1GDP_aPI3K re_87 Rac1GTP → Rac1GDP v28*Rac1GTP/(km28+Rac1GTP) re_88 IL2R + IL2  IL2R_IL2 kon*IL2R*IL2 – KD24*kon* IL2R_IL2 re_89 IL2R_IL2 → pIL2R + IL2 kp18* IL2R_IL2 re_90 pIL2R → IL2R v29*pIL2R/(km29+pIL2R) re_91 pIL2R → iIL2R v30*pIL2R/(km30+pIL2R) re_92 iIL2R → pIL2R v31*iIL2R/(km31+iIL2R) re_93  → mCISH v32*NFAT^2/(km32^2+NFAT^2) re_94  → mCISH kgen8 re_95 mCISH →  kdeg16*mCISH re_96  → CISH v33*mCISH /(km33 + mCISH) re_97 CISH →  kdeg17*CISH re_98  → mPLCγ kgen9 re_99 mPLCγ →  kdeg18*mPLCγ re_100  → PLCγ v34*mPLCγ /(km34 + mPLCγ) re_101 PLCγ →  kdeg19*PLCγ re_102 PLCγ + pTCR  PLCγ_pTCR kon*PLCγ*pTCR – KD25*kon*PLCγ_pTCR re_103 PLCγ_pTCR → aPLCγ + pTCR kp19*PLCγ_pTCR re_104 aPLCγ + CISH  aPLCγ_CISH kon*aPLCγ*CISH – KD26*kon*aPLCγ_CISH re_105 aPLCγ_CISH → uPLCγ + CISH kp20*aPLCγ_CISH re_106 uPLCγ →  kdeg20*uPLCγ re_107  → mPD1 v35*NFAT^2/(km35^2 + NFAT^2) re_108 mPD1 →  kdeg21* mPD1 re_109  → PD1 v36*mPD1*(km36 + mPD1) re_110 PD1 →  kdeg22* PD1 re_111 pTCR + PD1 → TCR + PD1 kc37*(PD1)*pTCR/(km37+pTCR) re_112 aRas + PD1  aRas_PD1 kon*aRas*PD1 – KD28*kon*aRas_PD1 re_113 aRas_PD1 → Ras + PD1 kp22*aRas_PD1 re_114 aPI3K + PD1  aPI3K_PD1 kon*aPI3K*PD1 – KD29*kon*aPI3K_PD1 re_115 aPI3K _PD1 → PI3K + PD1 kp23*aPI3K_PD1 It is assumed that the receptors TCR, CTLA -4, PD1 and IL2R have different binding motifs for different proteins, so the reactions with blue color do not affect the mass balance for the receptors. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Table A.2. Table of the initial conditions (nM) for each species. TCR 10.0 pTCR 0.0 CTLA4 0.1 iTCR 0.0 CN 10.0 aCN 0.0 RCAN 0.1 aCN_RCAN 0.0 pRCAN 0.0 aCN_pRCAN 0.0 Carabin 0.1 aCN_Carabin 0.0 pNFAT 9.906 NFAT 0.001283 pNFAT_aCN 0.0 pNFAT_aCN_pRCAN 0.0 GSK3 4.507 NFAT_GSK3 0.0 mNFAT 0.1 mRCAN 0.1 pTAK1 0.0 RCAN_pTAK1 0.0 mCarabin 0.1 mTNFa 0.1 TNFa 0.12 TAK1 9.239 mCTLA4 0.1 mIL2 0.1 mFasL 0.1 IL2 0.1 FasL 0.1 Ras 10.0 aRas 0.0 aRas_Carabin 0.0 MEK 10.0 aRas_MEK 0.0 pMEK 0.0 Rac1GTP 0.1078 Rac1GTP_MEK 0.0 ERK 10.0 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint pMEK_ERK 0.0 pERK 0.0 PI3K 10.0 pIL2R 9.442 PI3K_pIL2R 0.0 PI3K_pTCR 0.0 aPI3K 0.0 PP2A 10.0 PP2A_CTLA4 0.0 aPP2A 0.0 Akt 8.853 Akt_aPI3K 0.0 pAkt 1.147 pAkt_aPP2A 0 GSK3_pAkt 0.0 pGSK3 5.493 Rac1GDP 9.892 Rac1GDP_aPI3K 0.0 IL2R 0.5258 IL2R_IL2 0.0 iIL2R 0.03187 mCISH 0.1 CISH 10.0 mPLCg 0.1 PLCg 10.0 PLCg_pTCR 0.0 pPLCg 0.0 pPLCg_CISH 0.0 PLCg_CISH 0.0 uPLCg 0.0 mPD1 0.1 PD1 0.1 aRas_PD1 0.0 aPI3K_PD1 0.0 E1 1.0 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint Table A.3. Upper and lower bounds of model parameters. Parameter Type Lower Limit Upper Limit Units Reference Description kon 10-3 1.0 nM-1min- 1 Assumed Association rate KDi 1.0 106 nM Assumed Binding affinity kci 0.0 6*104 min-1 [92] Catalytic rate in Michaelis-Menten equations kpi 0.0 6*104 min-1 [92] Catalytic rate vi (for catalytic reaction) 10-4 105 nM*min-1 [5], [92] Maximum velocity in Michaelis-Menten equations for unknown enzyme vi (for expression rate based on mRNA concentration) 10-4 103 nM*min-1 Assumed Maximum velocity in Michaelis-Menten equations for protein production based on mRNA expression kmi 10.0 105 nM [5], [8] Michaelis-Menten constant kgeni 10-6 0.4 nM*min-1 Estimated based on [1], [2], [3], [4] Baseline generation rate of mRNAs kdegi (mRNA) 0.008 0.23 min-1 [93] Baseline degradation rate of mRNAs (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint kdegi (proteins) 0.00008 0.0012 min-1 [94] Baseline degradation rate of proteins Table A.4. Fitted parameters. Parameter Values Units kon 0.744759 nM-1min-1 KD4 223700.1 nM KD5 168118.4 nM KD6 45002.22 nM KD7 35960.29 nM KD8 78363.78 nM KD9 722630.7 nM KD10 415214.8 nM KD13 453032.9 nM KD14 961656 nM KD15 293.516 nM KD16 906154.5 nM KD17 170870.2 nM KD18 677861.2 nM KD19 20799.91 nM KD20 451422.9 nM KD21 109834.2 nM KD22 243654.7 nM KD23 981692.4 nM KD24 9437.089 nM KD25 779131.4 nM KD26 247916.6 nM KD28 142418.6 nM KD29 426916.3 nM kp2 52917.64 min-1 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint kp3 47468.17 min-1 kp4 17341.5 min-1 kp5 46471.26 min-1 kp8 53.5484 min-1 kp9 58494.47 min-1 kp10 22373.75 min-1 kp11 37764.75 min-1 kp12 14643.93 min-1 kp13 47265.41 min-1 kp14 29050.98 min-1 kp15 45695.24 min-1 kp16 49236.73 min-1 kp17 130.4833 min-1 kp18 59510.23 min-1 kp19 28986.8 min-1 kp20 15644.21 min-1 kp22 53384.22 min-1 kp23 15354.08 min-1 kc5a 57388.77 min-1 kc1 26783.54 min-1 kc2 25227.83 min-1 kc37 10702.25 min-1 kc38 45070.62 min-1 kc39 59.07676 min-1 v6a 1940.749 nM*min-1 v3 42764.89 nM*min-1 v4 67904.71 nM*min-1 v5 16506.41 nM*min-1 v6 764.9546 nM*min-1 v7 841.3915 nM*min-1 v8 51.22405 nM*min-1 v9 243.4736 nM*min-1 v10 71250.36 nM*min-1 v11 951.8989 nM*min-1 v12 994.5142 nM*min-1 v13 887.4255 nM*min-1 v14 322.8619 nM*min-1 v15 64820.01 nM*min-1 v16 921.2605 nM*min-1 v17 367.4828 nM*min-1 v18a 66755.12 nM*min-1 v18b 90244.25 nM*min-1 v19 95.03896 nM*min-1 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint v20a 99478.15 nM*min-1 v20b 66403.2 nM*min-1 v21 754.1809 nM*min-1 v22 93567.16 nM*min-1 v23 73934.94 nM*min-1 v24 30499.61 nM*min-1 v25 17172.81 nM*min-1 v26 11686.16 nM*min-1 v27 29234.4 nM*min-1 v28 57883.74 nM*min-1 v29 30073.09 nM*min-1 v30 75811.64 nM*min-1 v31 92779.58 nM*min-1 v32 0.135135 nM*min-1 v33 999.735 nM*min-1 v34 69.27322 nM*min-1 v35 941.9566 nM*min-1 v36 441.9933 nM*min-1 km6a 91278.7 nM km5a 22127.18 nM km1 54590.13 nM km2 58320.07 nM km3 72420.89 nM km4 29245.67 nM km5 98188.35 nM km6 4236.314 nM km7 26841.6 nM km8 38745.78 nM km9 46150.41 nM km10 77655.69 nM km11 6575.344 nM km12 5268.861 nM km13 23090.31 nM km14 69872.06 nM km15 42230.83 nM km16 7757.289 nM km17 81674.67 nM km18 44431.25 nM km18b 22899.24 nM km19 74325.13 nM km20 7061.907 nM km20b 62219.96 nM km21 12292.1 nM (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint km22 887.574 nM km23 50621.05 nM km24 65056.24 nM km25 89176.11 nM km26 95710.72 nM km27 99814.19 nM km28 39089.5 nM km29 71266.83 nM km30 37001.62 nM km31 9628.729 nM km32 99987.58 nM km33 2930.999 nM km34 67837.13 nM km35 24243.93 nM km36 62854.4 nM km37 64706.31 nM km38 61341.04 nM km39 99669.55 nM kgen1 0.09331 nM*min-1 kgen2 0.123895 nM*min-1 kgen3 0.010411 nM*min-1 kgen4 0.234104 nM*min-1 kgen5 0.120324 nM*min-1 kgen6 0.081743 nM*min-1 kgen7 0.121839 nM*min-1 kgen8 0.399069 nM*min-1 kgen9 0.373886 nM*min-1 kdeg1 0.054886 min-1 kdeg2 0.000438 min-1 kdeg3 9.35E-05 min-1 kdeg4 0.190494 min-1 kdeg5 0.000476 min-1 kdeg6 0.009132 min-1 kdeg7 8.59E-05 min-1 kdeg8 0.110141 min-1 kdeg9 9.63E-05 min-1 kdeg10 0.15187 min-1 kdeg11 0.001123 min-1 kdeg12 0.229291 min-1 kdeg13 0.000541 min-1 kdeg14 0.149548 min-1 kdeg15 0.000816 min-1 kdeg16 0.11542 min-1 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint kdeg17 0.001096 min-1 kdeg18 0.227438 min-1 kdeg19 0.000107 min-1 kdeg20 0.000353 min-1 kdeg21 0.021341 min-1 kdeg22 8.27E-05 min-1 log10_sigma1 -0.47502 log10_scale1 -0.18046 log10_offset1 0.18815 log10_sigma2 -2.25911 log10_scale2 -2.29188 log10_offset2 -4.37036 log10_sigma3 -1.63241 log10_scale3 -1.36121 log10_offset3 -0.25962 log10_sigma4 -2.3683 log10_scale4 0.373133 log10_offset4 -4.99998 log10_sigma5 -0.88629 log10_scale5 1.055226 log10_offset5 0.062407 log10_sigma6 -0.77074 log10_scale6 1.490558 log10_offset6 -0.13829 log10_sigma7 1.341567 log10_scale7 -1.00375 log10_offset7 1.494577 log10_sigma8 2.227199 log10_scale8 0.940156 log10_offset8 -4.98264 log10_sigma9 -2.9999 log10_scale9 -3.62915 log10_offset9 -4.69992 log10_sigma10 1.294673 log10_scale10 0.684183 log10_offset10 2.120215 log10_sigma11 -0.77374 log10_scale11 -2.5305 log10_offset11 -0.76581 log10_sigma12 1.049281 log10_scale12 -1.54426 log10_offset12 1.226346 log10_sigma13 -2.9999 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint log10_scale13 -0.98793 log10_offset13 -0.05328

[1] B. Jones, D. Stekel, J. Rowe, and C. Fernando, “Is there a liquid state machine in the bacterium Escherichia coli?,” Proceedings of the 2007 IEEE Symposium on Artificial Life, CI -ALife 2007, pp. 187–191, 2007, doi: 10.1109/ALIFE.2007.367795. [2] H. Chen, K. Shiroguchi, H. Ge, and X. S. Xie, “Genome -wide study of mRNA degradation and transcript elongation in Escherichia coli,” Mol Syst Biol , vol. 11, no. 1, p. 781, Jan. 2015, doi: 10.15252/MSB.20145794. [3] S. Kühner et al., “Proteome organization in a genome -reduced bacterium,” Science (1979), vol. 326, no. 5957, pp. 1235–1240, Nov. 2009, doi: 10.1126/SCIENCE.1176343. [4] J.-M. Anaya, Y. Shoenfeld, A. Rojas -Villarraga, R. A. Levy, and R. Cervera, “Autoimmunity,” Autoimmunity: From Bench to Bedside, pp. 237–378, Jul. 2013, Accessed: Jun. 11, 2025. [Online]. Available: https://www.ncbi.nlm.nih.gov/books/NBK459447/ [5] D. Davidia et al. , “Global characterization of in vivo enzyme catalytic rates and their correspondence to in vitro kcat measurements,” Proc Natl Acad Sci U S A , vol. 113, no. 12, pp. 3401–3406, Mar. 2016, doi: 10.1073/PNAS.1514240113. [6] Y. Wang, C. L. Liu, J. D. Storey, R. J. Tibshirani, D. Herschlag, and P. O. Brown, “Precision and functional specificity in mRNA decay,” Proc Natl Acad Sci U S A , vol. 99, no. 9, pp. 5860 –5865, Apr. 2002, doi: 10.1073/PNAS.092538799,. [7] S. B. Cambridge, F. Gnad, C. Nguyen, J. L. Bermejo, M. Krüger, and M. Mann, “Systems -wide proteomic analysis in mammalian cells reveals conserved, functional protein turnover,” J Proteome Res, vol. 10, no. 12, pp. 5275–5284, Dec. 2011, doi: 10.1021/PR101183K. [8] S. Y. Shin, M. W. Kim, K. H. Cho, and L. K. Nguyen, “Coupled feedback regulation of nuclear factor of activated T-cells (NFAT) modulates activation-induced cell death of T cells,” Sci Rep, vol. 9, no. 1, p. 10637, Dec. 2019, doi: 10.1038/S41598-019-46592-Z. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted September 10, 2025. ; https://doi.org/10.1101/2025.09.05.674558doi: bioRxiv preprint