The flowchart shows how the populations of each of the 8 types of cells affect each other (9 types if one includes G, the thymus). In order to maintain a constant system size, the thymus produces an equal number of cells as those that go to Cell Death.

The Equations:

Example: The change in the Initial self population (Is) is equal to the number of new Initial cells from the thymus [ = + k1 * G * SI ], plus cells converting from Es to Is at a rate k4  [ = + Es * k4 ] , minus Is to Es at a rate k2 [ =  - Is * k2 ],  minus Is to As at a rate k3 [ = - Is * k3], minus the number of initial cells that die [ = - Is * k5 ].

In equation form, this statement looks like this:

Is' = (k1*G* SI) + (Es*k4) - (Is*k2) - (Is*k3) - (Is*k5);

Now we have a system of equations to calculate the cell population in each group. It is convenient to think of X' as the (current time step) variables that we calculate knowing X data (time step -1). The program executes the equations in the execution order (ExecOrder) as listed in the table. However, the most meaningful equations are those where the bulk of the cell are converted; namely, those listed as ExecOrder #6, which are listed first. Equations in the other ExecOrders are mostly supporting rules. These rules were seperated out, because some calculations were too complex to be included as a single rule within the bulk conversion equations (#6). Most of the detailed complexities are a result of constraints that are placed on certain variables. These supporting rules include such constraints as not to exceed the number of APCs, or to normalize outgoing rates from the auto-catalytic rate equations, or to simply prevent cell division if there isn't enough cells to displace for such a thing to occur, without exceeding the defined total system size.
 

ExecOrder	Variable		Equation
6)	Is'	=	Is + ((N)*SI) + (Es1*k4) - (Is*k2) - (Is*k3) - (Is*k5) + ((M)*SI);
6)	As'	=	As + (Is*k3) - P - Q;
6)	Es'	=	Es1 + (Is*k2) - (Es1*k4) + Q;
6)	Inse'	=	Inse + ((N)*(1-SI)*L) + (Ense1*k4) - (Inse*k2) - (Inse*k3) - (Inse*k5) + (M*(1-SI)*L);
6)	Anse'	=	Anse + (Inse*k3) - U - V;
6)	Ense'	=	Ense1 + (Inse*k2) - (Ense1*k4) + V;
6)	Insu'	=	Insu1 + ((N)*(1-SI)*(1-L)) - (Insu1*k2) - (Insu1*k5) + (Ensu1*k4) + (M*(1-SI)*(1-L));
6)	Ensu'	=	Ensu1 + (Insu1*k2) - (Ensu1*k4);
where: (most are internal variables, with the exception of T and N)
1)	M	=	Cells Added to the system to initially grow it. Equals 0 once the system reaches Total size (T).
1)	T'	=	(T + M) = (Is' + As' + Es' + Inse' + Anse' + Ense' + Insu' + Ensu')
2)	P	=	As*k6;
			Note: P is normalized with Q if ((P+Q) > As)
2)	Q	=	(As/(R*SI)) * ((Es+Ense*XRxn)/(R*SI)) * k7 * PIEs * (R*SI);
			Note: Q must be <= # APC, and normalized with P if ((P+Q) > As)
2)	U	=	Anse*k6;
			Note: U is normalized with V if ((U+V) > As)
2)	V	=	(Anse/(R*(1-SI)*L)) * ((Ense+Es*XRxn)/(R*(1-SI)*L)) * k7 * PIEnsu * (R*(1-SI)*L);
			Note: V must be <= #APC, and normalized with U if ((U+V) > As)
3)	Es1	=	Es + Effector cells that are now dividing that went through the As to Es (via k7) before.
3)	Ense1	=	Ense + Effector cells that are now dividing that went though the Anse to Ense (via k7) before.
4)	Insu1	=	Insu - Es1 - Ense1
5)	N	=	P + (Is*k5) + U + (Inse*k5) + (Insu1*k5) = cells that die = cells that are re-born = k1*G
Executed when Foreign AG is initially introduced, to convert the % of unengaged cells to the engaged type:
7)	Inse'	=	Insu'*L';
7)	Ense'	=	Ensu'*L';
8)	Insu'	=	Insu' - Inse';
8)	Ensu'	=	Ensu' - Ense';
PIE Calculations:
9)	PIE-s'	=	(1-e(-(As'/APC)/(SI*R))) * (1-e(-(Es'/APC)/(SI*R))) * APC;
9)	PIE-nse'	=	(1-e(-(Anse'/APC)/(L'*(1-SI)*R))) * (1-e(-(Ense'/APC)/(L'*(1-SI)*R))) * APC;
9)	PIE-nsu'	=	(1-e(-(Insu'/APC)/((1-L')*(1-SI)*R))) * (1-e(-(Ensu'/APC)/((1-L')*(1-SI)*R))) * APC;
Effector / peptide calculations:
10)	Es'/p	=	Es'/(R*SI);
10)	Ens'/p	=	Ense'/(R*(1-SI)*L');
10)	Ensu'/p	=	Ensu'/(R*(1-SI)*(1-L'));

Further model notes:

Autocatalytic k7 rate:
The autocatalytic rate k7, which convert As-type to Es-type cells (also true for nonself engaged cells, Anse to Ense) can use some explanation. First, the simple case without being autocatalytic or interacting with APC, the number of cells going through this pathway would simply be As*k7 = # of new Es cells through this pathway. Adding catalytic feedback from Es makes the equation: As*Es*k7. This equation implies that any self-specific peptide on As interacts with any self-specific peptide on Es. This isn't the case as they must be specific to the exact same peptide for them to interact. Thus, we must divide out the number of peptides for each As and Es which gives us the following for each peptide: (As/ (R* SI))* (Es/ (R* SI))* k7. Now, we include the probability that an As-APC-Es interaction will happen in the same physical location, namely a Priming Inductive Event (PIEs) would yield, (As/(R*SI)) * (Es/(R*SI)) * k7 * PIEs. For more details on PIE, see the math details in the Th-Genesis Minimal Model.  However, recall that this calculation is in terms of per peptide, so multiply by the number of self peptides (R*SI) gives the entire amount of cells converting from As to Es as: (As/(R*SI)) * (Es/(R*SI)) * k7 * PIEs * (R*SI).

k6 vs k7 competition on As & Anse:
Since k6 is large (default value of 0.1386) and k7 is autocatalyic, and thus effectively changes with respect to the concentration of the anticipatory and effecotor cells, there is a competition for As (likewise for Anse). We place a few constraints: The first constraint is that no more than the number of APC/p*(number of peptides) may go through the k7 pathway. And second, if the pull from the k6 and k7 pathway exceeds the number of Anticipatory cells, the two pathways are normalized. The Normalizing Equation for values A and B, with to a maximum of C is A' = C*A/(A+B) and B' = C*B/(A+B).

Instant conversion of Insu to Inse and Ensu to Ense cells:
Equations in the execution order of #7 & #8 in the table above result in instant appearance of a large number of Inse and Ense cells. This conversion is equal to the nonself antigenic load times the unengaged cells. This is performed because of the actuality that Insu and Ensu cells are no different than Inse and Ense cells. Their only difference is that Inse and Ense cells would react with the newly introduced foreign antigen. Thus these unengaged cells, would now engage, and thus are renamed accordingly. This also helps in keeping the math straight, as non-self engaged becomes symmetrical mathematically to self.

Steady State:
Steady state values can be observed by allowing the system to continue until the numbers within each cell groups converge on a number and stays there regardless of additional time steps. It is important to note that the system is generally allowed sufficient time after the initial system growth stage to become steady state, before the introduction of the foreign antigen. After addition of the antigen, the system again will adjust itself to find a new steady state. Generally, it does. However on occasions the system will oscillate slightly. This is due to the fact that cell divisions are timed to be exactly n-steps, and cease to divide when there is no Insu to displace. Eventually Insu type cells re-generate, and after some delay, effector cell divisions gets up to speed and take them down again. Thus causing a small oscillation.

The Application

1. Get the k and system values from the user.
2. Starts with all cell populations = 0.
3. Iterate through the Equations for SEED number of steps, while L=0

and gradually grow the system via an internal stem cells variable M = (input Total) / (input SEED steps),
4. Iterate through the Equations for (BEFORE - SEED) number of steps while L=0.
5. Iterate through the Equations for AFTER steps, with L= the input value for load.
6. The simulator runs to completion. The resulting data can be visualized through the online graphing, or using the raw data output.

Last Modified: November 1, 2002