The flowchart shows how the populations of each of the 6 types of cells
effect each other. Because each individual cell has a chance of changing
into another type of cell at every unit of time, for the model to be steady
the number of cells that stop being type X must also equal the number of
cells that become type X.
For example, lets say that there are 5000 Ens cells at the beginning
of a unit of time. If during that unit of time 400 cells decide to become
Ins (the only type of cell an Ens can become) then 400 cells
of some other type have to decide to become Ens for the population
of Ens cells to remain at 5000. According to the flowchart those
new cells will come from the Ins population
The Deltas
So how does this flowchart become a system of equations? The chart shows
us a path of population changes, or the net change in population.
Example: The change in the anticipatory population (As) is equal
to the number of new anticipatory cells, minus the number of anticipatory
cells that die. For self, the new anticipatory cells is given by Is
times the probability of an Is cell becoming anticipatory (or k3),
and the number of anticipatory cells that die is given by k6 times the
number of anticipatory cells.
In an equation this statement looks like this:
As' = (k3 * Is)  (k6 * As)
Where Is, and As are current populations of initial and
anticipatory cells respectively, and As' is the new calculated population
of anticipatory cells.
Now we have a system of equations for the net change of each population:
Is' = (SI * N) + (k4 * Es)  (k5 * Is)  (k3 * Is)  (k2 * Is)
Es' = (k2 * Is)  (k4 * Es)
Ins' = ((1  SI) * N) + (k4 * Ens)  (k5 * Ins)  (k2 * Ins)
Ens' = (k2 * Ins)  (k4 * Ens)
As' = (k3 * Is)  (k6 * As)
Note that N is equal to k1*G in the graph.
We also need an equation to represent the steady size of the entire
system. This is accomplished by setting the number of cells entering the
system (being created by the stem cells) equal to the number of cells leaving
the system (by dying). N is the variable the represents the number
of cells produced per unit time.
N = (k5 * Is) + (k5 * Ins) + (k6 * As)
The N variable is also unknown, but it is different from other cell
populations because this value represents the number of cells that are
born and that die in the system, for each unit of time, when a steady state
is reached. So we can get a formula for the "change in N per
unit time" which is fixed to 0.
N' = (k5 * Is) + (k5 * Ins) + (k6 * As)  N
Thus we have
Is' = (SI * N) + (k4 * Es)  (k5 * Is)  (k3 * Is)  (k2 * Is)
Es' = (k2 * Is)  (k4 * Es)
Ins' = ((1  SI) * N) + (k4 * Ens)  (k5 * Ins)  (k2 * Ins)
Ens' = (k2 * Ins)  (k4 * Ens)
As' = (k3 * Is)  (k6 * As)
N' = (k5 * Is) + (k5 * Ins) + (k6 * As)  N
Although we how have 6 equations with 6 variables, we notice that these
equations are not linearly independent. In short, Is'
= Es' + Ins' + Ens' + As' + N'.
Another important input of the model is the total system size or Protecton
size (T). This is introduced by the following equation, which
effectively scales all the cell groups by the correct amount to have an
immune system of size T.
Is + Es + Ins + Ens + As = T
This is now a system of seven equations (six of which are linearly independent)
with six unknowns. For simplicy at this point, we will drop the Is'
equation, since it can be formed from linear combinations of other equations.
The k values and total system cell population are decided by the user
and the model represented by equation, once solved (see "The Steady, Stable
State"), can be used to obtain the steady state values of the populations
for each cell type. This can be performed quickly by a computer.
The web interface allows for the user to scan ranges of values for any
variable, creating a 2D graph. Furthermore it allows the user to scan 2
variables at the same time, to obtain a 3D graph. All the while, raw numbers
of the selected values are also returned.
A Steady, Stable State
For the model to have a steady state we would need to have population values
for the 6 types of cells that when a system is started at these levels,
they will stay at these levels as long as there is no outside interference.
The first requirement of a steady, stable state is that the net change
of each population equal zero. So our equation for the previous section
become...
0 = (k2 * Is)  (k4 * Es)
0 = ((1  SI) * N) + (k4 * Ens)  (k5 * Ins)  (k2 * Ins)
0 = (k2 * Ins)  (k4 * Ens)
0 = (k3 * Is)  (k6 * As)
0 = (k5 * Is) + (k5 * Ins) + (k6 * As)  N
0 = Is + Es + Ins + Ens + As  T
When the variables are separated from the known values we can form two
matrices one known and one unknown. Then we can use Maple (TM) to
solve for the unknown matrix. Like so...

k2*Is 
k4*Es 
0*Ins 
0*Ens 
0*As 
0*N 

= 

0 

0 
0 
(k2k5)*Ins 
k4*Ens 
0 
(1SI)*N 
0 
0 
0 
k2*Ins 
k4*Ens 
0 
0 
0 
k3*Is 
0 
0 
0 
k6*As 
0 
0 
k5*Is 
0 
k5*Ins 
0 
k6*As 
1*N 
0 
1*Is 
1*Es 
1*Ins 
1*Ens 
1*As 
0 
T 
Becomes...

k2 
k4 
0 
0 
0 
0 

* 

Is 

= 

0 

0 
0 
(k2k5) 
k4 
0 
(1SI) 
Es 
0 
0 
0 
k2 
k4 
0 
0 
Ins 
0 
k3 
0 
0 
0 
k6 
0 
Ens 
0 
k5 
0 
k5 
0 
k6 
1 
As 
0 
1 
1 
1 
1 
1 
0 
N 
T 
When Maple(TM) is given this system of equations and told to solve for
Is, Es, Ins, Ens, As & N with
T,
SI,
k2k6 being an input constant. We get these results:
Is = 
T * SI * k5 * k4 * k6 
D 

Ins = 
T (k3  k5 + SI * k5 + SI * k3)* k6 * k4 
D 

Ens = 
T * k2 * (k3  k5 + SI * k5 + SI * k3) * k6 
D 

Es = 
T * k2 * SI * k5 * k6 
D 

N = 
T * k6 * k5 * k4 * (k5  k3) 
D 

As = 
T * K * k5 * k4 * k3 
D 

Where D = (k5*k6*k2  k6*k4*k3  k6*k2*k3  K*k5*k4*k3 + K*k6*k4*k3 +
K*k6*k2*k3  k5*k4*k6).
So the system is solvable as long as D doesn't equal zero and all of
the populations are positive, for the given K values. This solution will
give us an instant calculation for the steady state, but how do we know
that it is a stable steady state? And what's the difference between
a steady state and a stable steady state? A stable steady state
is one that is resistant to long term change. A pendulum is a stable system
because, baring outside interference, the pendulum will always come to
the same state, straight down. A coin balanced on it's edge is an unstable
state because if an outside force knocks the coin over, the coin will never
be able to get up on it's own.
The Eigenvalues
Note: this section requires an understanding of linear algebra.
An eigenvalue is a number, y, that corresponds to a vector,
x,
and a matrix, A, such that A*x = y*x.
A homogeneous system of equations has nonzero solutions if and only
if the coefficient matrix (y*I  A) = 0 is singular; that is, if
and only if the determinant is zero.
The system of equations will equal the zero vector if Is = Es
= Ens = Ins = As = N = 0, but we are not interested
in an empty system. So we have to show that there is a non zero solution.
Because of the nature of this problem (we are applying this matrix
of k values an infinite number of times) the eigenvalues cannot be negative,
that would result in an alternating pair of values that doesn't converge.
Finding the determinant of this matrix is best done by computer, and
even then the solution isn't usable. Maple gives the determinant as a 5th
degree polynomial, the simplest coefficient of which is (k6 + k3 + 2*k2
+ 1 + 2*k5). Again not solvable by hand, but since the k values are all
between 0 and 1 we can run a simulation to see if negative eigenvalues
will be produced.
100 randomly picked sets of k values didn't produce a single negative
eigenvalue, so we can feel somewhat confident that the model is stable.
The Application
This is how the simulator works.

Get the immune system size (T), number of APCs per peptide, repertoire
size (R), fraction of new cells which are antiSelf (SI),
and the various rate constants (k2 to k5) from the user.

Start with all cell populations set to 0, except for the stem cells (N).

Use the solutions for the system of questions above, obtained from Maple(TM),
to quickly obtain the steady state values for each cell type, given the
user's inputs.
i.e. Is = (T*SI*k5*k4*k6)/D, and N = (T*k6*k5*k4*(k5k3))/D,
etc.

Display the predicted steady state values to the user.

For variables that are graphed, (1 variable gives a 2D graph, 2 variables
results in a 3D graph), when this option is selected, there are a minimum
of 100 samples taken per decade. Steps number 2 and 3 above are performed
for each new sample and then used to create the graph, which is then returned
to the user via the web page.