/*
 * A Mathematical Model for Chronic Myelogenous Leukemia (CML)
 * and T Cell Interaction
 * 
 * Model Status
 * 
 * This model runs in OpenCell to recreate the published results.
 * This particular version of the model recreates figure 6 from
 * the paper (where CML decreases over a 2 year period). The units
 * have been checked and they are consistent. The CellML model
 * also runs in COR however due to the timescale being in days
 * the model is not ideal for simulation in COR.
 * 
 * Model Structure
 * 
 * ABSTRACT: In this paper, we propose and analyse a mathematical
 * model for chronic myelogenous leukemia (CML), a cancer of the
 * blood. We model the interaction between naive T cells, effector
 * T cells, and CML cancer cells in the body, using a system of
 * ordinary differential equations which gives rates of change
 * of the three cell populations. One of the difficulties in modeling
 * CML is the scarcity of experimental data which can be used to
 * estimate parameters values. To compensate for the resulting
 * uncertainties, we use Latin hypercube sampling (LHS) on large
 * ranges of possible parameter values in our analysis. A major
 * goal of this work is the determination of parameters which play
 * a critical role in remission or clearance of the cancer in the
 * model. Our analysis examines 12 parameters, and identifies two
 * of these, the growth and death rates of CML, as critical to
 * the outcome of the system. Our results indicate that the most
 * promising research avenues for treatments of CML should be those
 * that affect these two significant parameters (CML growth and
 * death rates), while altering the other parameters should have
 * little effect on the outcome.
 * 
 * The complete original paper reference is cited below:
 * 
 * A mathematical model for chronic myelogenous leukemia (CML)
 * and T cell interaction, Helen Moore and Natasha K. Li, 2004,
 * Journal of Theoretical Biology, 11, 369-391. PubMed ID: 15038986
 * 
 * reaction diagram
 * 
 * [[Image file: moore_2004.png]]
 * 
 * Cell population diagram showing the population dynamics of each
 * cell type and how they interact with each other. These dynamics
 * and interactions are described by the mathematical model equations.
 */

import nsrunit;
unit conversion on;
unit day=86400 second^1;
unit first_order_rate_constant=1.1574074E-5 second^(-1);
unit cells = fundamental;
unit cells_per_microlitre=1E9 meter^(-3)*cells^1;
unit flux=1.1574074E4 meter^(-3)*second^(-1)*cells^1;
unit microlitre_per_cells_day=1.1574074E-14 meter^3*second^(-1)*cells^(-1);

math main {
	realDomain time day;
	time.min=0;
	extern time.max;
	extern time.delta;
	real Tn(time) cells_per_microlitre;
	when(time=time.min) Tn=1510.0;
	real sn flux;
	sn=0.37;
	real dn first_order_rate_constant;
	dn=0.23;
	real kn first_order_rate_constant;
	kn=0.062;
	real eta cells_per_microlitre;
	eta=720.0;
	real C(time) cells_per_microlitre;
	when(time=time.min) C=10000.0;
	real Te(time) cells_per_microlitre;
	when(time=time.min) Te=20.0;
	real alpha_n dimensionless;
	alpha_n=0.14;
	real alpha_e first_order_rate_constant;
	alpha_e=0.98;
	real de first_order_rate_constant;
	de=0.30;
	real gamma_e microlitre_per_cells_day;
	gamma_e=0.057;
	real Cmax cells_per_microlitre;
	Cmax=230000;
	real rc first_order_rate_constant;
	rc=0.0057;
	real dc first_order_rate_constant;
	dc=0.024;
	real gamma_c microlitre_per_cells_day;
	gamma_c=0.0034;

	// <component name="environment">

	// <component name="Tn">
	Tn:time=(sn-(dn*Tn+kn*Tn*(C/(C+eta))));

	// <component name="Te">
	Te:time=(alpha_n*kn*Tn*(C/(C+eta))+alpha_e*Te*(C/(C+eta))-(de*Te+gamma_e*C*Te));

	// <component name="C">
	C:time=(rc*C*ln(Cmax/C)-(dc*C+gamma_c*C*Te));

	// <component name="model_parameters">
}