Overview
To tackle the changing tumor antigens in vivo, our team adopt an antibody library strategy. Our initial idea, Plan A, is to place the entire antibody library on each CAR-NK92 cell, i.e. each CAR-NK92 contains all the necessary scFvs needed to correspond the tumor antigen. However, this idea is quickly dismissed by modeling predictions that if multiple scFvs are expressed simultaneously on a single cell, it can easily disrupt normal cells with low tumor antigens expression. Furthermore, research shows that transfection of multiple CARs on a single cell is challenging and expression of multiple CARs on a single cell is prone to ligand non-dependent tonic. Therefore, we choose Plan B: expressing multiple CARs on a large number of immune cells. Whether this plan supports our perspective involves a number of considerations: in terms of tumor-killing effects, we investigate the competing process between CAR-NK92 and tumor cells using the metacellular automata and Lotka-Volterra Predator-Prey model and predict the effectiveness of suicide genes being silenced subsequent to CAR recognition; in terms of safety, we calculate AP1903 dose and injection cycle for inactivating overdose of CAR-NK92, and simulate and predict the size of CAR-NK92 library through target antigen screening, in order to achieve high coverage rate. These considerations above are ultimately visualised in our Metacellular Automata Simulation.
Plan A
Modeling of number of CAR types
Introduction
There is a maximum carrying capacity of CARs on a single immune cell. We hope to employ a typical scfv to help estimate the maximum value of the cell membrane carrying capacity and to analyse the disadvantages of expressing too many CARs in a single cell by bioinformatic analysis.
Assumption
1.CAR-NK-92 is a regular sphere.
2.CAR length is much smaller than the cell radius , so the effect of CAR length is neglected.
Description
1.Maximum value of the cell membrane carrying capacity
Through electron micrographs, we measure the actual diameter of CAR-NK92 cells at about 3μm-5μm.1The surface area is calculated as about 113-320 μm2 without considering the cell deformation. We choose the commonly CD19 scFv(FMC63.3)protein sequence for homology modeling.2
Due to the flexibility and quantum nature of the molecule, it is impossible to have a constant value for its volume. According to the scale, 2×3×3 nm volume range can completely wrap around scFv.3
CAR surface aggregation can lead to misabsorption of VH and VL, resulting in tonic signaling.4
The fluid mosaic model proposed by Singer and Nicolson states that membrane proteins account for about 40%-50% of the membrane, and the occupied area of proteins is calculated to be about 56.5482 μm2.
Assuming that all the expressed proteins on the surface of NK-92 membrane are CARs, it can accommodate up to 6.98 × 105. In reality, it is impossible to express only CARs, because other functional proteins also need to be expressed. And there are disadvantages to expressing too many CARs in a single cell.
2.Disadvantages of expressing too many CARs in a single cell
Studies have shown that multi-targeted CAR-T shows positive efficacy in the treatment of multiple malignancies and becomes an effective way to prevent antigen evasion, but this type of CAR-T cells can be recognized by multiple target antigens, leading to an excessive increase in the level of NFAT activation and over-activation of immune cells, which can easily over-kill normal cells.5~8
Thus, we aimed to analyze the effect of the number of types of CARs expressed by individual immune cells on the rate of antigen evasion and accidental injury to normal cells through mathematical modeling.
2.1 Obtaining genes highly expressed in tumor cells
Taking colon cancer as an example, firstly we download transcription datas from TCGA database and use "Differential Expression Analysis" to obtain genes that are highly expressed in cancer cells compared to peritumoral tissue. The setting difference is 2.5-fold, P<0.05, and finally we obtain the volcano map and heat map.
To further filter genes whose expression products are membrane proteins, our team cooperate with ZJUintl IGEMer and writes a webscraping program in python.
Finally we obtain 54 genes.
2.2 Peritumoral Tissue killed Tendency
An important reason for normal cells to be accidentally injured is that the expression of some proteins is close to that of cancer cells. We use a function to process gene expression information to predict the tendency for normal cells low in expression of target antigens to be accidentally injured (Table 1).
2.3 Antigen Evasion Rate
Assuming only one type of CAR, the antigen evasion rate may be on a 9-point gradient from 0.1 to 0.9. Using the evasion rate of 0.1 as the threshold, the number of CAR types that reach the threshold limit is obtained for different initial evasion rates. In the hypothesis of an initial evasion rate of 0.6, the evasion rate drops to 0.7776 when the number of CAR types is 5. When the initial evasion rate is higher at 0.9, at least 22 CAR types are required to bring the evasion rate down below the threshold.
Results
Furthermore, research has shown that transfection of multiple CARs on a single cell is challenging and expression of multiple CARs on a single cell is prone to ligand non-dependent tonic signaling9. Therefore, we choose Plan B: expressing multiple CARs on a large number of immune cells.
References
1. Liu, Dongfang et al. “Chimeric antigen receptor (CAR)-modified natural killer cell-based immunotherapy and immunological synapse formation in cancer and HIV.” Protein & cell vol. 8,12 (2017): 861-877. doi:10.1007/s13238-017-0415-5.
2. Majzner, Robbie G et al. “Tuning the Antigen Density Requirement for CAR T-cell Activity.” Cancer discovery vol. 10,5 (2020): 702-723. doi:10.1158/2159-8290.CD-19-0945.
3. Asaadi Y, Jouneghani FF, Janani S, Rahbarizadeh F. A comprehensive comparison between camelid nanobodies and single chain variable fragments. Biomark Res.
4. Jayaraman, Jayapriya et al. “CAR-T design: Elements and their synergistic function.” EBioMedicine vol. 58 (2020): 102931. doi:10.1016/j.ebiom.2020.102931.
5. Gill, S., M.V.Maus, and D.L. Porter, Chimeric antigen receptor T cell therapy: 25years in themaking. Blood Reviews, 2016. 30(3): p. 157-167.
6. Abreu, T.R., et al., Currentchallenges and emerging opportunities of CAR-T cell therapies. Journal ofControlled Release, 2020. 319: p. 246-261.
7. Koichi Hirabayashi et al. NatureCancer, 2021, doi:10.1038/s43018-021-00244-2.
8. Han X, Wang Y, Wei J, Han W.Multi-antigen-targeted chimeric antigen receptor T cells for cancer therapy. JHematol Oncol. 2019;12(1):128.
9. Ajina A, Maher J. Strategies to Address Chimeric Antigen Receptor Tonic Signaling. Mol Cancer Ther. 2018 Sep;17(9):1795-1815. doi: 10.1158/1535-7163.MCT-17-1097. PMID: 30181329; PMCID: PMC6130819.
Appendix
1. Transcriptome data from TCGA database: first obtain the expression spectrum. The v22 version of gff3 (http://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_22/)and the v33 version of gff3 (http://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_33/gencode.v33.annotation.gff3.gz) are downloaded from GENCODE (https://www.gencodegenes.org/human/) gencode.v22.annotation.gff3.gz) , then extract the mapping information of GeneSymbol and ENSG_ID from them, and subsequently map ENSG_ID to GeneSymbol.
2. Limma (linear models for microarray data, DOI:10.1093/nar/gkv007) is a differential expression screening method based on a generalized linear model, and here we used the R package limma (version 3.40.6) for differential analysis to obtain different comparison groups versus controls differential genes between different comparison groups and controls.
please click here to see our Plan B !
Plan B
Mathematical modeling of CAR-NK92 proliferation
Overview
NK-92 cells are ideal chassis cells. We use the Lotka-Volterra Predator-Prey model to quantify the proliferative capacity of CAR-NK92 and its ability to kill tumor cells. And we get the ideal time point for injecting CAR-NK-92 library by modeling.
Introduction
Lotka-Volterra Predator-Prey model
Around 1925, the American mathematician Alfred J. Lotka (1880-1949) and the Italian mathematician Vito Volterra (1860-1940) proposed a model to describe the variation in the number of preys and predators. The two differential equations at the heart of this model, the Lotka-Volterra equation, give an indication of how the population size of predators and prey varies depending on the number and interaction of the two parties.
Description
Parameters
Assumption
The initial number of tumor cells is 30 and CAR-NK-92 is 5.
Process
However, in reality most tumor-killing systems do not exist periodic oscillations , but tend to some equilibrium point. Improvements are as follows:
x/N1 refers to the effect of environmental resources on x. The above two interaction model does not take into account environmental factors, so we add a logistic term.a/r y refers to the effect of the "predator" (CAR-NK-92) on the population size of tumour cells.
Similarly,
Plotting the (5) and (6), then there is an equilibrium point and the results are more realistic (Figure 2).
The meaning of the phase diagram is that if we take the initial value of any point, there is a trend in its evolutionary direction, and we can adjust the parameters, such as the initial number of CAR-NK-92, to make the trend point close to the y-axis (all tumor cells are killed) and provide guidance for subsequent experiments: what amount of CAR-NK-92 can completely eliminate tumor cells.
We adopt a more realistic scenario:
The initial number of tumor cells in a local tissue of 300 × 300 units size is 300, and the initial injection volume of CAR-NK-92 is N ml. There are five types of CAR-NK-92 recognized by this tumor surface antigen, five of each, and a total of 25 are able to expand and actually complete the killing process, and the rest will apoptosis due to the modulation of AP1903. The maximum environmental accommodation of tumor cells in local tissues is 5 × 104.
Results
In clinical situation, tumor cells are often difficult to be detected when they start to proliferate. As can be seen from Figure 4, when tumor cells first start to proliferate, it is difficult for CAR-NK92 to detect tumor cells and recruit around them due to the small number of tumor cells. When the tumor cells grow to a certain number, CAR-NK92 starts to initiate proliferation and killing mechanism.
Injected in the first 80 days, the killing effect is the best. The latest injection timing should be controlled on the 70th to 80th day, and the killing effect is basically the same when injected at this time point compared with the first day.
References
1. de Pillis LG, Radunskava AE, & Wiseman CL (2005) A validated mathematical model of cellmediated immune response to tumor growth. Cancer research 65(17):7950-7958.
2. Barsoum IB, Smallwood CA. Siemens DR, & Graham CH (2014) A mechanism of hypoxia-mediated escape from adaptive immunity in cancer cells. Cancer research 74(3):665-674.
3. Fotios Mpekris. Chrysovalantis Voutouri. JamesW. Baish.et al. Combining microenvironment normalization strategies toimprove cancer immunotherapy.Proceedings of the National Academy ofSciences Feb 2020,117 ⑺ 3728-3737.
Mathematical modeling of CAR-NK92 signaling pathway
Overview
Phosphorylation activation of the chimeric antigen receptor (CAR) is a key step in initiating suicide gene repression, which determines the ability of the CAR-NK92 cell library to make targeted compositional changes in different tumors. Since CARs contain well-characterized TCR signaling components, we draw on previous mathematical modeling of TCRs and improve the original mathematical model for structural differences between CARs and TCRs. Based on the improved model, the results interpret important characteristics of CAR activation, such as the sensitivity, specificity and transition time of activation.
Background
Original mathematical model
To explain the process of T-cell receptor phosphorylation from a theoretical and modeling perspective, McKeithan proposed a kinetic proofreading model for multi-step T-cell receptor phosphorylation in 1995 based on the experimental phenomenon of multi-locus phosphorylation.1 Although other models have been proposed over the years in the field of T-cell activation and have explained biological phenomena at different levels, there have been no major improvements or developments in McKeithan's kinetic calibration model.2This illustrates that the model can be used to analyze the specific recognition response of the TCR to exogenous molecules. Equations and their description are displayed in table 2.
Structural differences between CARs and TCRs
As previously described, CAR-NK92 uses a second-generation CAR as an antigen recognition receptor. The second-generation CAR differs structurally from the TCR in two major ways. First is the different number of phosphorylation sites: TCR has 10 tyrosine sites on its ITAM, meaning it has 20 steps of phosphorylation, while the second-generation CAR has only 3 tyrosine sites on its ITAM, meaning only 6 steps of phosphorylation.3Second, the extracellular and intracellular segments of the second-generation CAR are not independent of each other as in the TCR, which means that the signaling part cannot walk randomly between the kinase rich domains(KDR) and the phosphatase-rich environment(PRE), resulting in less rapid dephosphorylation of the CAR than in the TCR.4
Assumption and Description
As mentioned before due to structural differences, CARs dephosphorylate more slowly than TCRs. So the McKeithan's assumption of instantaneous dephosphorylation needs to be modified for CARs. We then set the unbound state Ri corresponding to the complex Ci after the dissociation of the exogenous molecule.Ri is converted to Ri-1 according to a dephosphorylation reaction in the opposite direction from Ci-1 to Ci. The exogenous molecule M can bind to any unbound state of Ri to form Ci. This leads to a reaction network of sequential phosphorylation from C0 to CN, followed by sequential dephosphorylation from RN to product, which is bridged by the reversible M-binding of Ri to form Ci.(see Figure 1)
For such a reaction network we still assume that for any i, the phosphorylation reaction from Ci to CN has the same rate α, the rate of dephosphorylation from Ri to Ri-1 is the same as β, and the rates of binding and dissociation for the reversible process of M binding Ri to form complex Ci are both kon and koff. So for each Ri and Ci in the reaction network, we can write its dynamical differential equation.(see Table 3)
CAR sensitivity and specificity with parameter β
It is assumed that the system is continuously activated by the exogenous signal and reaches and maintains a quasi-equilibrium state. The total number of receptors in the system and the concentration of exogenous molecules remain constant, and the number of receptors in the system in the bound and unbound states does not change. Letting the left side of dynamical differential equation be zero, a simple variable substitution gives:
In real-life biological systems, CAR-NK cells exhibit the same excellent specificity and sensitivity as T cells.5,2Is there then a reasonable interval for the dephosphorylation rate β that can combine these two important properties? Since our model does not assume a limited case for the value of β as the McKeithan model does, it can be used to analyze this question theoretically.
Activation characteristic time of the reaction system
Activation characteristic time is another important property of CAR. Here, we use the method of defining the transient activation time from the literature.8The activation characteristic time from moment 0, when CAR-NK92 cells start to be stimulated by exogenous molecules, to the final quasi-equilibrium state is defined as:
Results
Determination of the range of values of β
We define the proportion of CN to the total number of receptors to characterize sensitivity, with parameters taken in the physiological range respectively6:α=0.01s-1,kon=105s-1·mol-1·L,koff=0.1s-1,M=10-6mol/L,N=6.
In terms of sensitivity, which is in line with the result of McKeithan's model, the larger the value of β, the lower the proportion of CN in the total number of receptors; conversely, the higher the proportion. (see pink line in Figure 2)
In defining specificity, we draw on the experimental observations in the literature.7Two exogenous molecules with very similar morphological and biochemical properties are assumed, where molecule 1 (L1) is the true exogenous antigen and molecule 2 (L2) is not. The only difference is that the rate of dissociation of molecule 1 from the CAR(koff) is half that of molecule 2. The slower rate of dissociation will result in a higher proportion of CN being activated by molecule 1 than by molecule 2. The Ratio of L1 to L2 is defined as the specificity of the CAR cell for the two molecules. (see blue line in Figure 2)
Combining the results, the system has excellent sensitivity to molecule 1 when β is in the 7×10-3~9×10-3~s-1 interval with over 50% CARs activated, while the proportion of CN to the total number of receptors is more than an order of magnitude of that of molecule 2.
Calculation of the activation characteristic time
We adopt the idea and mathematical treatment of a completely random activation model of multiple phosphorylation sites in the literature.8In addition to the good specificity and sensitivity of the random activation strategy, the activation characteristic time is only in the order of tens to hundreds of seconds, a result consistent with biological experimental observations.7
Conclusion
Based on the original McKeithan kinetic calibration model, our CAR activation model is developed by relaxing various constraints, including the range of dephosphorylation rate β and the order of phosphorylation and dephosphorylation, to obtain a more consistent phosphorylation model for CAR. In terms of the results, the biological properties of CAR-NK92 cells in terms of sensitivity, specificity and time efficiency in the recognition of exogenous molecules are reflected more successfully and realistically. It has implications for predicting the activity of CAR-NK92 cell gene circuits in animal and clinical experiments.
References
1. McKeithan, T. W. Kinetic proofreading in T-cell receptor signal transduction. Proc. Natl. Acad. Sci. 92, 5042–5046 (1995).
2. George, A. J. T., Stark, J. & Chan, C. Understanding specificity and sensitivity of T-cell recognition. Trends Immunol. 26, 653–659 (2005).
3. Lindner, S. E., Johnson, S. M., Brown, C. E. & Wang, L. D. Chimeric antigen receptor signaling: Functional consequences and design implications. Sci. Adv. 6, eaaz3223 (2020).
4. Burroughs, N. J. & Van Der Merwe, P. A. Stochasticity and spatial heterogeneity in T-cell activation. Immunol. Rev. 216, 69–80 (2007).
5. Myers, J. A. & Miller, J. S. Exploring the NK cell platform for cancer immunotherapy. Nat. Rev. Clin. Oncol. 18, 85–100 (2021).
6. Stone, J. D., Chervin, A. S. & Kranz, D. M. T-cell receptor binding affinities and kinetics: impact on T-cell activity and specificity. Immunology 126, 165–176 (2009).
7. Daniels, M. A. et al. Thymic selection threshold defined by compartmentalization of Ras/MAPK signalling. Nature 444, 724–729 (2006).
8. Salazar, C. & Höfer, T. Versatile regulation of multisite protein phosphorylation by the order of phosphate processing and protein-protein interactions: Kinetic models of multisite phosphorylation. FEBS J. 274, 1046–1061 (2007).
Modeling of AP1903-induced apoptosis of CAR-NK92
Introduction
AP1903 is a novel small molecule dimer, which plays an important role in the process of suicide route induction and construction in this project.On the one hand, after AP1903 administration, a chimeric protein in cells changes from a static state to an activated state, which fuses the intracellular part of Fas receptor with the drug binding domain of FKBP, thus initiating Fas signal transduction and ultimately leading to cell apoptosis.1The iCas9 suicide gene we use in the suicide circuit has been shown to be clinically effective, as part of caspase-9,fused with the dimerization domain of FKBP-12,and the dimerization of FKBP-12 is closely related to AP1903.2On the other hand, AP1903 induces the CAR-activated suicide gene expression inhibition circuits.Therefore, AP1903 is involved in the induction of suicide gene circuit through the above two parts.
Description
We attempt to use AP1903 to limit the number of CAR-NK-92 cells with iCasp9 to a range that would make CAR-NK-92 safe and effective.Therefore, we construct a model to estimate the appropriate injection time and concentration of AP1903.(Table 1. for Parameters)
Concentration & Efficacy
It is known that the efficacy of AP1903 decreases with the decrease of drug concentration, so we establish the simplest linear model to represent this trend:
According to the enzyme kinetic theory, Michaelis-Menten equation can be used to describe the metabolism of free drugs.3Based on this, we optimize the linear model equation and obtained a relatively complete effect change equation:
According to this equation and relevant pharmacokinetic parameters of AP1903 in human body obtained by referring to literature, using MatLab software tools to assist,the pharmacodynamic and time-dependent images of AP1903 are plotted at intervals of 12h(Figure 1) and 24h(Figure 2).
Concentration & Time
In order to ensure the stable and efficient operation of the suicide circuit while ensuring safety, we try to keep the concentration of AP1903 within a relatively stable range.To achieve this goal, we first use the NSG mice as the initial model to inject AP1903 at equal intervals for many times, and try to use mathematical means to calculate the initial time of drug injection and the time interval of subsequent injection by using formulas:
As we know, after AP1903 with the initial concentration of α is introduced into the body, the concentration shows a decreasing trend over time:
The basic formula of AP1903 concentration decreasing with time is obtained:
Our aim is to limit the concentration of AP1903 to a certain range of multiple injections so that it can function well.Therefore, we assume that the concentration range is β~γ(β<γ), the initial injection time is T1, and the injection interval is T0. According to (1), we obtain:
According to (1) and (2), the concentration after the initial injection of AP1903 is:
After that, in the process of injecting AP1903 at equal intervals, the changes of AP1903 in each cycle are as follows:
According to (1), (2), (3) and (4), the initial injection time and interval injection time are respectively:
Results
According to literature review, the optimal concentration range of AP1903 for injection in NSG mice is 2.5-5mg /kg.2,4By substituting the data into (5) and (6), the initial injection time of AP1903 in the NSG mouse model is 0, and the interval injection time is 50.0106h (4 decimal places are reserved).This result is basically consistent with the data obtained by referring to relevant AP1903 literature,4 so the calculation formula is valid.
We attempt to further apply the formula to the prediction and guidance of clinical AP1903 injection(Figure 3). According to the equivalent dose conversion method, (5), (6) and model, we obtain that the initial injection time of AP1903 in human body is 3.8147h, and the interval injection time is 4.8060h (4 decimal places are reserved).
References
1. Iuliucci JD, Oliver SD, Morley S, Ward C, Ward J, Dalgarno D, Clackson T, Berger HJ. Intravenous safety and pharmacokinetics of a novel dimerizer drug, AP1903, in healthy volunteers. J Clin Pharmacol. 2001 Aug;41(8):870-9. doi: 10.1177/00912700122010771. PMID: 11504275.
2. Amatya C, Pegues MA, Lam N, Vanasse D, Geldres C, Choi S, Hewitt SM, Feldman SA, Kochenderfer JN. Development of CAR T Cells Expressing a Suicide Gene Plus a Chimeric Antigen Receptor Targeting Signaling Lymphocytic-Activation Molecule F7. Mol Ther. 2021 Feb 3;29(2):702-717. doi: 10.1016/j.ymthe.2020.10.008. Epub 2020 Oct 14. PMID: 33129371; PMCID: PMC7854354.
3. Zhang Q, Wang GJ. [Advances and related issues in the use of in vitro methods to predict metabolic clearance rate of new drugs]. Yao Xue Xue Bao. 2007 Oct;42(10):1023-8. Chinese. PMID: 18229605.
4. Fu W, Lei C, Wang C, Ma Z, Li T, Lin F, Mao R, Zhao J, Hu S. Synthetic libraries of immune cells displaying a diverse repertoire of chimaeric antigen receptors as a potent cancer immunotherapy. Nat Biomed Eng. 2022 Jun 6. doi: 10.1038/s41551-022-00895-1. Epub ahead of print. PMID: 35668107.
Mathematical modeling of CAR antibody library coverage
Introduction
The issue of CAR library security is certainly a key issue that cannot be bypassed. But, how large is a CAR library sufficient to cover all tumor antigens in vivo? Here, we provide a program simulation strategy for predicting the size of CAR library.
Assumption
1.As about 1010 ScFv known to be in the full CAR library1, we arrange and numbered all scFvs and their specific antigens from 1 to 1010 in the order of structural changes, which means that the narrower the gap between the ordinal number of a certain scFv and a certain antigen is, the higher the affinity is.
2.An one-dimensional array of 106 (or 107,108) ordinal numbers are randomly generated from the full library as the precast CAR library. Based on the results of the proliferation model, about 4~5 CAR-NK92 cells exerts target killing effects after CAR-NK92 library injection, so 5 ordinal numbers are randomly selected from the whole library species as the tumor antigen.
3.We subtract each value in the array from the 5 random ordinal numbers representing the antigen and take the minimum absolute value. We call this minimum absolute value the affinity gap. CAR libraries with the size of 106 have shown good coverage in the previous proliferation models. And after completing one hundred sets of pre-experiments for a CAR library of size 106, we find that the average value of affinity gap for each antigen for a CAR library of this size is around 5600. Considering the discrete values and large variance that appeared in the pre-experiments, we finally set the threshold value of the affinity gap at 10,000.
Description
In each set of simulated experiments, when the affinity gap of five antigens is inferior to the threshold value of the affinity gap, CAR library effectively covers all the tumor antigens in vivo in this set of experiment. When all sets of simulations are completed, the number of effectively covered groups divided by the total number of simulated groups is calculated as the coverage rate of this size CAR library. The flow chart of the program simulation experiment is presented in figure 1.
Results
After running the experimental simulation program, we find that the coverage rate can approach 100% when the library size reaches 107 and 108. However, the larger the CAR library capacity is, the higher the technical difficulty and cost of implementation. Furthermore, when the size reaches 107, there are about 20 effective scFvs per tumor antigen, which is sufficient to ensure that a sufficient number of CAR-NK92s are activated and perform tumor killing tasks. Therefore, a reasonable CAR bank size should be between 106 and 107. The results are presented in Table 1.
References
1. Sheets, M. D. et al. Efficient construction of a large nonimmune phage antibody library: The production of high-affinity human single-chain antibodies to protein antigens. Proc Natl Acad Sci USA 6 (1998) doi:10.1073/pnas.95.11.6157
Metacellular automata simulation
Overview
1.Simulate the process of tumor cell removal by CAR-NK92 library using metacellular automata and successfully resolved the tumor heterogeneity.
2.Figure out that NK92 with different CARs have different recognition ability to tumor cells:"Central Blossom" killing (one point of recognition and central point of breakthrough) and "Fireworks" killing ( multiple points of recognition and multiple points of breakthrough).
3.Find the "protective circle" of tumor cells and the "killing limit" mechanism of CAR-NK92 are important factors for the remnants of tumor, which provide guidance for subsequent experiments.
Introduction
·Synthetic biology
Synthetic biology is an engineering discipline. Part of the engineering is simulation and modeling so that system behavior can be determined prior to building the design.
To simulate the process of cell-to-cell confrontation, we NMU-China team introduced metacellular automata to visualize and simulate the dynamic process of proliferation and killing.
·Metacellular automata
A lattice dynamics model in which time, space and state are discrete, and spatial interactions and temporal causality are local, have the ability to simulate the spatio-temporal evolution of complex systems. In brief, metacellular automata is like a Go game. When you put a piece in a grid, it must have an effect on the grid around it. Metacellular automata is a similar idea, changing the state in some lattices and thus affecting the surrounding lattices.
Assumption
CAR-NK92 injected at any time point during the first 80 days is equally effective in killing tumor cells according to the Lotka-Volterra Predator-Prey model, so we assume that CAR-NK92 cells are injected on day 80 and all other parameters are the same.
Parameters
To make the simulation more realistic, we consider as many parameters as possible:
The number of cellular mutations is generally evaluated in terms of tumor mutation burden (TMB), which is generally expressed as the number of mutations per 1 Mb (1 megabase). Ten somatic mutations per Mb of coding DNA correspond to approximately 150 nonsynonymous mutations in expressed genes.1~2
However, only 10% of these nonsynonymous mutations can produce mutant peptides that bind with high affinity to immune cells.3 In turn, only 1% of the peptides that can bind with high affinity to immune cells can be recognized by immune cells in tumor patients and activate the killing mechanism to kill tumor cells.
This means that for 1000 nonsynonymous mutations, only 1 neoantigen may eventually be recognized by the immune cells. The probability of producing a neoantigen is about 1/1000.
To generate CAR response elements, a CAR inducible promoter containing 6 NFAT-REs in a minimal IL-2 promoter was placed upstream of the transcription factor Gal4-KRAB. Moreover, the inducible caspase-9 suicide gene that induces apoptosis upon specific binding with the small molecule dimerizer CID AP1903 was placed downstream of the expression under the control of the combined 5 SV40/UAS promoters(Figure 1).
Because the KRAB protein has been demonstrated to be capable of inhibiting all promoters within at least 3 kB 36 and because the forward construct would therefore inhibit other promoters, we generated the opposite construct, in which the two expression cassettes were cloned in a manner such that both promoters were at the opposite ends at a distance greater than 4 kb. The 2A peptide sequence was intercalated between iCASP9 and the GFP tag.(Figure 2)
All gene codes of CAR response elements were show in Figure 3. As we consider the mutation rate of stable genomes estimated to be 10-10/bp per cell generation,4the mutation rate of suicide gene is 2.4674×10-7.
mutation rate of suicide gene =(216/6+960+75+1308+442/5)×1010=2.4674×10-7
We add a random number function and appropriately adjusted the neighborhood rule to increase the randomness of proliferation direction. The random wandering model is also borrowed to make the boundary continuous by considering the cell movement factor, which makes the simulation more realistic.
To maintain the concentration of AP1903 at 2.5~5mg/kg that satisfies:
1.Don't affect CAR-NK92 that has already recognized the target antigen.
2.Keep the number of CAR-NK92s that temporarily fail to recognize the target antigen at a low number, rather than all apoptosis. After the emergence of neoantigen, the reserve can recognize them and start proliferation rapidly.
For more details, please click on Table 2.
Results
Conclusion
1.The proliferation rate of different mutants is different, which is affected by the density of other cells around at that time, the higher the density, the lower the resource utilization, the slower the proliferation.
2.CAR-NK92 that recognize different target antigens have different efficiency in binding to the target antigens, ranging from "Central Blossom" killing (one point of recognition and central breakthrough) to "Fireworks" killing (multiple points of recognition and multiple points of breakthrough).
3.The reasons why some tumor cells survive :
I.Protective circle: the mutants wrap up a small number of original tumor cells during the proliferation process. In the subsequent killing process, CAR0 that recognize the original tumor cells basically apoptotic due to the action of AP1903. After the protective circle of the mutant is cut open, CAR0 have withered away and original tumor cells have nothing to fear.
II. Killing limit: because CAR-NK92 cannot kill infinitely, we set the killing limit of a single CAR-NK92 cell to 20. It is possible that these surviving cells happen to encounter CAR-NK92 cells that have reached the kill limit during every wave of killing, and then escape by chance.
4.The safety of the gene lines is high. The suicide gene lines of various CAR-NK92 are not mutated during the whole simulation.
References
1. Hanahan D, Weinberg RA. Hallmarks of cancer: the next generation. Cell 2011;144:646-74.
2. Yarchoan M, Johnson BA, 3rd, Lutz ER, Laheru DA, Jaffee EM. Targeting neoantigens to augment antitumour immunity. Nature reviews Cancer 2017.
3. https://en.m.wikipedia.org/wiki/Game_theory
4. Schumacher TN, Schreiber RD. Neoantigens in cancer immunotherapy. Science 2015;348:69-74.