sys = BioSystem();

sys.AddConstant('k_tc',0.01);
sys.AddConstant('k_mdeg',0.003);
sys.AddConstant('k_tl',0.1);
sys.AddConstant('k_pdeg',0.03);
sys.AddConstant('k_1',2);
sys.AddConstant('k_2',1.5);
sys.AddConstant('k_mideg',0.003);
sys.AddConstant('k_cdeg',0.005);

dDNAdt = sys.AddCompositor('DNA',0);
dmiRNAdt = sys.AddCompositor('miRNA',10^11);
dRNAcomdt = sys.AddCompositor('RNAcom',0);
dmRNAdt = sys.AddCompositor('mRNA',0);
dGFPdt = sys.AddCompositor('GFP',0);

sys.AddPart(Part('translation',[dDNAdt dmiRNAdt dRNAcomdt dmRNAdt dGFPdt],[...
    Rate('0') ...
    Rate('k_2*RNAcom-k_1*miRNA*mRNA-k_mideg*miRNA') ...
    Rate('k_1*miRNA*mRNA-k_2*RNAcom-k_cdeg*RNAcom') ...    
    Rate('k_tc*DNA+k_2*RNAcom-k_1*miRNA*mRNA-k_mdeg*mRNA') ...
    Rate('k_tl*mRNA-k_pdeg*GFP')]));

[T,Y] = sys.run_pulses([...
    Pulse(0,'DNA',10), ...
    Pulse(12000,'DNA',0)]);

figure();
plot(T,Y(:,sys.CompositorIndex('DNA')),...
    T,Y(:,sys.CompositorIndex('mRNA')),...
    T,Y(:,sys.CompositorIndex('RNAcom')),...
    T,Y(:,sys.CompositorIndex('GFP')),...
    T,Y(:,sys.CompositorIndex('miRNA')));
legend('DNA','mRNA','RNAcom','GFP','miRNA');
xlabel('time');
ylabel('concentration');
