package simulation;
import java.util.Vector;
import miscellaneous.DirectoryManager;
import miscellaneous.FileUtility;
import miscellaneous.LSVList;
import reaction.Item;
import reaction.Pool;
import reaction.Reaction;
import recorders.Recorder;
import rtools.R;
public class SimplifiedModel {
public static String POL = "Pol Protein", TEMPLATE = "Template", rTEMPLATE = "rTemplate", RIB = "Ribosome", POLt = "Pol mRNA", Nt = "N mRNA", N = "N Protein";
public static double alphaLtNt = 0.1; // ratio of L transcript to N transcript
public static double ribRate = 50; // 18 nt/s; // 6 aa/s [Use nt/s here]
public static double ribSpacing = 238.5; // 238.5; //nt per ribosome
public static double ribDelay = ribSpacing / ribRate; // minimum time between ribs
public static double polRate = 15; // template chose to match expt data with kon = 1.0; 3.7 nt/s Iverson 1981; // nucleotides per second polumerates rate on template
public static double polSpacing = 500; // 170 nucleotides
public static double polDelay = polSpacing / polRate; // minimum time between pols
public static double lTemplate = 11000;
public static double lN = 1333;
public static double lPOL = 5000;
public static double tNt = (1.0 / 6.0) * lTemplate / polRate; // time delay for transcription
public static double tPOLt = lTemplate / polRate; // time delay for transcription
public static double tNp = lN / ribRate; // time delay for transcription
public static double tPOLp = lPOL / ribRate; // time delay for transcription
public static double initRibs = 5e5;
public static double initTemplates = 1;
public static double initPols = 50;
public static double kTranslation = 1 / (ribDelay * initRibs / 3);
public static double kReplication = 5e-12;
public static double tf = 30500; // 60 * 60 * 12;
public static double dt = 0.1;
public static double dtRecord = 1; // 60 * 5;
public static Pool cell = new Pool("Cell");
public static void main(String[] args)
{
runSimplifiedModel(true);
}
public Reaction getTranscriptionRxn()
{
Reaction ret = new Reaction();
// POL + TEMPLATE
ret.addReactant(cell, POL, 1, 1);
ret.addReactant(cell, TEMPLATE, 1, 1);
// ---> POL + TEMPLATE + Nt
ret.addProduct(cell, POL, (1 - alphaLtNt), tNt);
ret.addProduct(cell, POL, alphaLtNt, tPOLt);
ret.addProduct(cell, TEMPLATE, 1, polDelay);
ret.addProduct(cell, Nt, 1, tNt);
ret.addProduct(cell, POLt, alphaLtNt, tPOLt);
// K
ret.setk(1258 * 1000 * kReplication);
return ret;
}
public Reaction getNtRxn()
{
Reaction ret = new Reaction();
// POL + TEMPLATE
ret.addReactant(cell, RIB, 1, 1);
ret.addReactant(cell, Nt, 1, 1);
// ---> POL + TEMPLATE + Nt
ret.addProduct(cell, RIB, 1, tNp);
ret.addProduct(cell, Nt, 1, ribDelay);
ret.addProduct(cell, N, 1, tNp);
// K
ret.setk(kTranslation);
return ret;
}
public Reaction getPOLtRxn()
{
Reaction ret = new Reaction();
// POL + TEMPLATE
ret.addReactant(cell, RIB, 1, 1);
ret.addReactant(cell, POLt, 1, 1);
// ---> POL + TEMPLATE + Nt
ret.addProduct(cell, RIB, 1, tPOLp);
ret.addProduct(cell, POLt, 1, ribDelay);
ret.addProduct(cell, POL, 1, tPOLp);
// K
ret.setk(kTranslation);
return ret;
}
public Reaction getTemplateRxn()
{
Reaction ret = new Reaction();
// POL + TEMPLATE
ret.addReactant(cell, POL, 1, 1);
ret.addReactant(cell, TEMPLATE, 1, 1);
ret.addReactant(cell, N, 1258, 1);
// ---> POL + TEMPLATE + NmRNA
ret.addProduct(cell, POL, 1, tPOLt);
ret.addProduct(cell, TEMPLATE, 1, polDelay);
ret.addProduct(cell, rTEMPLATE, 1, tPOLt);
// K
ret.setk(kReplication);
return ret;
}
public Reaction getrTemplateRxn()
{
Reaction ret = new Reaction();
// POL + TEMPLATE
ret.addReactant(cell, POL, 1, 1);
ret.addReactant(cell, rTEMPLATE, 1, 1);
ret.addReactant(cell, N, 1258, 1);
// ---> POL + TEMPLATE + NmRNA
ret.addProduct(cell, POL, 1, tPOLt);
ret.addProduct(cell, rTEMPLATE, 1, polDelay);
ret.addProduct(cell, TEMPLATE, 1, tPOLt);
// K
ret.setk(kReplication);
return ret;
}
public static void runSimplifiedModel(boolean plot)
{
DirectoryManager.setHostDirectory("/Users/jaywarrick/Desktop");
SimplifiedModel mod = new SimplifiedModel();
Vector<Double> resultsTemplates = new Vector<Double>();
Vector<Double> resultsrTemplates = new Vector<Double>();
Vector<Double> resultsRIBs = new Vector<Double>();
Vector<Double> resultsNts = new Vector<Double>();
Vector<Double> resultsPOLts = new Vector<Double>();
Vector<Double> resultsNs = new Vector<Double>();
Vector<Double> resultsPols = new Vector<Double>();
Vector<Double> time = new Vector<Double>();
// Initialize Cell Amounts
cell.initialize(new Item(TEMPLATE, initTemplates, 0.0));
cell.initialize(new Item(RIB, initRibs, 0.0));
cell.initialize(new Item(POL, initPols, 0.0));
Vector<Reaction> rxns = new Vector<Reaction>();
// Add reactions
rxns.add(mod.getNtRxn());
rxns.add(mod.getPOLtRxn());
rxns.add(mod.getrTemplateRxn());
rxns.add(mod.getTemplateRxn());
rxns.add(mod.getTranscriptionRxn());
LSVList reactions = new LSVList();
reactions.add("Reactions");
for (Reaction r : rxns)
{
reactions.add(r.toString());
}
Recorder.log(reactions.toString(), "Tester");
double t = 0;
double count = -1.0;
double temp = 0;
while (t <= tf)
{
// Add matured items to the pool
cell.mature(t);
temp = Math.floor(t / dtRecord);
if(count < temp)
{
count = temp;
resultsTemplates.add(cell.getCount(TEMPLATE));
resultsrTemplates.add(cell.getCount(rTEMPLATE));
resultsRIBs.add(cell.getCount(RIB));
resultsNts.add(cell.getCount(Nt));
resultsPOLts.add(cell.getCount(POLt));
resultsNs.add(cell.getCount(N));
resultsPols.add(cell.getCount(POL));
time.add(t);
}
// calculate rates based on each reaction in the simulation
for (Reaction rxn : rxns)
{
rxn.calculateBaseRate();
rxn.applyStoichAndRate();
}
double actualDt = cell.simpleLinearIntegrate(dt, t, 1);
Recorder.log("" + t, SimplifiedModel.class);
t = t + actualDt;
}
if(plot)
{
R.eval("x <- 0;");
R.eval("graphics.off();");
plot("templates", resultsTemplates, time);
plot("rTemplates", resultsrTemplates, time);
plot("RIBs", resultsRIBs, time);
plot("Nts", resultsNts, time);
plot("POLts", resultsPOLts, time);
plot("Ns", resultsNs, time);
plot("POLs", resultsPols, time);
}
}
public static void plot(String name, Vector<Double> data, Vector<Double> time)
{
R.makeVector(name, data);
R.makeVector("time", time);
R.eval("write.table(" + name + ",'/Users/jaywarrick/Desktop/temp/" + name + ".txt',row.names=FALSE);");
String file = R.startPlot("pdf", 10, 8, 300, 12, "Helvetica", null);
R.eval("plot(time, " + name + ", type = 'l', col='red', log='y')");
R.endPlot();
try
{
FileUtility.openFileDefaultApplication(file);
}
catch (Exception e)
{
e.printStackTrace();
}
}
}