Our project aims to regulate the activity of CAR-T cells in order to reduce the risk of inflammatory factor storm which do damage to patients. In order to further reveal how many drugs and how often we should give patients, we have modeled to solve these problems. In order to guide the amount of administration, we expect to be able to know how the amount of administration or the time after administration affect the effect of drug, in order to avoid the calculation trouble caused by the excessive definition of the variable, we divided the whole process into ①dox pharmacokinetics of Dox②Tet-On system's response and U24 protein expression and ③the impact of cell surface concentration of CAR and ligand on the immune intensity.
Dox is a tetracycline analogue that induces Tet-On system activation. Compared with tetracycline, it has better stability and lower toxicity. Therefore, we choose Dox as our inducer.
Observing the structure of the entire T cell activation signaling pathway and selecting the DNA-RNAP complex for sensitivity analysis, it can be easily seen that CAR is a key node and have biggest coefficient (assume kinetic order always close to 1), so regulate the amount of CAR is reasonable.
Figure 1:pathway of T-Cell activation[1]
As we all known,drugs are not all distributed in cells but are distributed between plasma and cells in human's body. We consulted the literature and obtained the transfer rate constant of drugs form the digestive system to plasma or form cytosol to the plasma as well as redistribution rate constant form cytosol to the plasma [2],[3]. After consider the drug's decomposition rate is consistent in the plasma and cell, we list reaction equation as follows
By solving the homogeneous ODEs(ordinary differential equations )obtained by the law of mass action, we finally get the function of the concentration of the drug in the cytosol about time and plot it as follows(assume dosage is 2mmol/kg). As can be seen from the analysis of the Dox -time relationship, Dox concentration-Time relationship follows a similar manner to a three-exponential function will reach a peak after about 8 hours, and then decrease quickly, the concentration of dox has dropped very low at about 4000 minutes. We also plot the situation of successive administration and find the drug concentration in the cytosol will stabilize when administered 2mmol per kg per 1000 minutes with certain parameters we use.
Figure 2: time dependent cytosol Dox concentration
Figure 3: situation of successive administration 2mmol Dox per kg per 1000 minutes
The change of dox over time is obtained. Next step is to consider the Response of the Tet-on system. Due to the complexity of the
above functions, we can not find the analytical solution so we assume that under normal conditions, the content of rtTA in the cells
is in steady state and the transcription rate is determined by the following formula:.
Among them, (1+dox[t]/kd) and (((dox[t])/kd))/((1+dox[t]/kd)) represent the ratio of time-dependent Dox-unbound /binding DNA, B and A
are the rate of leaked transcription and the rate of induced transcription, respectively, and kd is the apparent dissociation
constant. We implicitly include rtTA and some other constant in kinetic constants above according to literature and experiments.
Having obtained certain data [4]、[5]、[6]、[7]., we can get certain transcription rate and solve the differential equation with
considering the degradation of RNA which is listed below to obtain the relationship between RNA concentration and time .
Figure 4 : RNA concentration-Time relationship
After considering the translation rate constant by estimated from the kinetic constant and the ribosome concentration[7],[8],[9], we
try to solve the ODE about translation. However, due to the use of interpolation function, we can't solve the equation. At first, we
want to solve it by the perturbation method which mean determine the protein yield without degradation firstly and use the result to
estimate the degradation rate at each moment, but we found that it is difficult to construct a suitable perturbation function over a
long period of time because the value of the perturbation always become too colossal, thereby we having to give up this method.
Next, we remember that we can use the exponential decay law to estimate this result, so we estimate the yield-time relationship
according to the expression and plot it.
Qualitatively, the model is correct, but the drawing of the rate of production shows that such treatment overestimates the degradation. We use the average time, which is obtained by integrate the amount of protein without considering degradation with time and divide the result by the amount of protein, instead of the time to improve calculation result but it has not improved much. We give the relationship between degradation rate constant and u24 concentration-Time relationship, too. We find that the degradation rate constant affect maximum U24 concentration and retention time in great extent. It may possible to use variant degradation tag to control concentration of U24 and density of CAR in wide range.
Figure 5: u24 concentration-Time relationship
Figure 6: u24 concentration-Time relationship with consider degradation dependent on average time
(a) | (b) |
(c) | (d) |
(a)kde=0 (a)kde=0.003(a)kde=0.012(a)kde=0.0332
After surveying the equation and using the new micro-element ktrans rna[t] E^(-kdeu24(t0-t)dt instead of dy and integrate, we get an
integral upper bound function. After finding the numerical solution of it and plot it as follows, we find those attenuation result is
quite accuracy, RNA degrade slowly and lasts longer. Further calculation will be obtained after we get higher computing power.
In general, the subsequent refinement changes quantitatively, but has little effect on the peak shape, and is in good agreement with
the tendency of amount of dox.
F igure8: u24 concentration-Time with change integrate micro-element
Because we have some doubts about the accuracy of some constants use above, we hope to use Hill equation to get another new model. After obtain the protein yield-time relationship by experiment, we find that it roughly matches the Hill equation. In the form, we use the Hill equation to fit it which is often used. Since cytoplasm Dox concentration can be regarded as a constant during cell culture, according to the experimental data, expression
is fitted to obtain expression of concentration of u24, with fewer guess components than we obtain in above chapter. In this expression, A0 and A1 stand for the production rate of leak express or of induced express, ki stand for affinity constant between dox and DNA, h is Hill constant and B stand for the Cytotoxic of Dox. Data is obtained by our experiment about cell-culture with concentration gradient Dox.
Figure 9: u24 concentration-Time determined by Hill equation and experiment
The expression is
Obtained the Hill equation about Dox and u24, we can use Dox[t] which is obtained above instead constant Dox concentration and get a new graph about U24 and time.
We know that in immunology, surface receptors are generally transmit information to the downstream pathway through phosphorylation. Due to the phosphorylated region on the TCR is mainly the ITAM domain and there are three ITAM domains on the CAR, it is reasonable to use the concentration of CAR which three ITAMs are phosphorylated as the indicator of immune signal intensity.
We chose a model called "separation model" to describe the behavior of CAR during the immunization process. The main content of
this model is that large phosphatases will be "squeezed" out of the immune synapse to form a kinase-rich domain(KRD) during the
process of immunization process, T cell receptors that bind to pMHCs will undergo a series of phosphorylation and initiate downstream
pathways by recruiting molecules such as ZAP70. After reference to the dynamic proofreading model, we made a dynamic mode diagram of
CAR phosphorylation.[10],[11],[12],[13],[14]
Figure 10:Schematic of the model
According to the law of mass action, we list a series of ordinary differential equations(ODE) as follows to describe phosphorylation and use the density of the three-phosphorylated receptor(A3) as an indicator of immune intensity. Since the ordinary differential equations have no analytical solution, we adopt two approaches to solve this problem.
ODEs:
T^' [t]=-konT[t]*M[t]+koff(a1[t]+a2[t]+a3[t]+a0[t])
M^' [t]=-konT[t]*M[t]+koff(a1[t]+a2[t]+a3[t]+a0[t])
a0^' [t]=konT[t]*M[t]-koffa0[t]-pa0[t]+qa1[t]
a1^' [t]=pa0[t]-qa1[t]+qa2[t]-pa1[t]-koffa1[t]
a2^' [t]=pa1[t]-qa2[t]+negaa3[t]-posia2[t]-koffa2[t]
a3^' [t]=posia2[t]-negaa3[t]-koffa3[t]
T[0]=100,M[0]=300,a0[0]=0,a1[0]=0,a2[0]=0,a3[0]=0
Using matlab®, we get the numerical solution of the system of ODEs. It can be seen that the tertiary phosphorylated receptors gradually increased and eventually almost all of the receptors were phosphorylated. We also made a diagram about how the cell surface CAR density affects the density of the phosphorylated receptors.
Figure 11:A3-t
Figure 12:T0-M0-A3
We assume that all phosphorylated receptors are in steady state during the immune response which changes ODEs into a homogeneous equation, making it easier to solve. Using Mathematica®, we express A3 as a function of CAR receptor density T0 and its ligand density M0 and draw a three-dimensional graph of it. It can be seen that A3 is almost linearly positively correlated for T0 and M0, respectively. Taking the ratio of the density of A3 produced as a specificity index, we plot how T0 density affects CAR-T's specificity of cells have different ligand numbers and see when ratio of density of ligands between two cells is low, specificity depends on density of T0 linearly unless T0 lower than a certain number, which will make specificity vanish. When the ratio approach or exceed 0.3, specificity vanish totally. Using mathmatica®, we also plot how kinetic parameter affect our result.
Figure 13:Effects of CAR receptor density T0 and its ligand density M0 on A3
(a) | (b) |
(c) | (d) |
The diffusion of receptor between the kinase-rich domain and the phosphatase-rich environment is important. We use the Monte Carlo algorithm and the random walking model to simulate the diffusion of the CAR receptor on the cell surface and calculate the average retention time and equilibrium density ratio of those two regions. After add a diffusion-related term to the ODEs and chemical equation, we get a new Numerical Solution and find an interesting result——there is no steady state in this model, instead, production of A3m show a peak shape after the formation of immune synapses, and then returns to 0. We also used this model to do some simulations of KRD formation. New chemical equations we use are show below.
Figure 16: Result of simulation diffusion by Monte Carlo algorithm and the random walking mode
Figure 17:A3-t relationship with consider diffusion
Figure 18:a short movie simulates CAR on cell surface during T-cell activation,in which blue point stand for A0 and red point stand for A3
According to the literature, phosphorylated CAR will form cluster on the cell surface, and then zap70 will be recruited to initiate the downstream reaction. Here, we use ZAP70 concentration to more accurately reflect the immune intensity. To make calculate easier, we define the biggest cluster as a tetramer and assuming that each CAR recruits only one ZAP70, according to the aggregation- recruitment pattern and the low of mass action, and using the steady-state approximation, we solve the expressions of the different state complexes and determine the relationship of the amount of activated ZAP70. T0, M0 with certain various kinetic constants, the final result is too complicated to listed here[10]
Figure 19: schematic diagram of more accurately model(A) Schematic of the model(B) State space of an individual cluster. Each state is identified by the doublet of integers en;mT: the number of ligated receptors (n) and the number of those receptors that are ZAP70 positive (m). Note that m<=n.(C) Possible stochastic transitions, and associated rates, from a representative state (n;m)[15]
Based on the relationship between dox and A3, we can get a new function for the dose and actual effect and understand the
relationship between dose, administration time and immune intensity better. Combined with more clinical trials, this new function
will be more instructive. In combination with relationship between Dox and immune intensity, we can get guidance on the dose of dox
to the patient, too.
Using the modeling we modeled, we can learn that Dox should be 8 hours administered before the symptoms worsen. If long-acting
inhibition is required, it should be administered at a frequency of about 1 time/day (based on the dose of 5 mmol/kg). The
combination of models and more experimental studies in the clinic, we can better guide the application.