/*
 * A computational model on the modulation of mitogen-activated
 * protein kinase (MAPK) and Akt pathways in heregulin-induced
 * ErbB signalling
 * 
 * Model Status
 * 
 * This version of the model has been checked in COR and OpenCell
 * and it runs to replicate the published results. The units have
 * been checked and they are consistent.
 * 
 * Model Structure
 * 
 * ABSTRACT: ErbB tyrosine kinase receptors mediate mitogenic signal
 * cascade by binding a variety of ligands and recruiting the different
 * cassettes of adaptor proteins. In the present study, we examined
 * heregulin (HRG)-induced signal transduction of ErbB4 receptor
 * and found that the phosphatidylinositol 3'-kinase (PI3K)-Akt
 * pathway negatively regulated the extracellular signal-regulated
 * kinase (ERK) cascade by phosphorylating Raf-1 on Ser(259). As
 * the time-course kinetics of Akt and ERK activities seemed to
 * be transient and complex, we constructed a mathematical simulation
 * model for HRG-induced ErbB4 receptor signalling to explain the
 * dynamics of the regulation mechanism in this signal transduction
 * cascade. The model reflected well the experimental results observed
 * in HRG-induced ErbB4 cells and in other modes of growth hormone-induced
 * cell signalling that involve Raf-Akt cross-talk. The model suggested
 * that HRG signalling is regulated by protein phosphatase 2A as
 * well as Raf-Akt cross-talk, and protein phosphatase 2A modulates
 * the kinase activity in both the PI3K-Akt and MAPK (mitogen-activated
 * protein kinase) pathways.
 * 
 * The original paper reference is cited below:
 * 
 * A computational model on the modulation of mitogen-activated
 * protein kinase (MAPK) and Akt pathways in heregulin-induced
 * ErbB signalling, Mariko Hatakeyama, Shuhei Kimura, Takashi Naka,
 * Takuji Kawasaki, Noriko Yumoto, Mio Ichikawa, Jae-Hoon Kim,
 * Kazuki Saito, Mihoro Saeki, Mikako Shirouzu, Shigeyuki Yokoyama,
 * and Akihiko Konagaya, 2003, The Biochemical Journal, 373, 451-463.
 * PubMed ID: 12691603
 * 
 * model diagram
 * 
 * [[Image file: hatakeyama_networkdiagram_scaled_2003.png]]
 * 
 * A schematic diagram of the complete signalling network.
 */

import nsrunit;
unit conversion on;
unit s=1 second^1;
//Warning:  unit nm_ renamed from nm, as the latter is predefined in JSim with different fundamental units.
unit nm_=1E-6 meter^(-3)*mole^1;
unit flux=1E-6 meter^(-3)*second^(-1)*mole^1;
unit second_order_rate_constant=1E6 meter^3*second^(-1)*mole^(-1);
unit first_order_rate_constant=1 second^(-1);

math main {
	realDomain time s;
	time.min=0;
	extern time.max;
	extern time.delta;
	real kf1 second_order_rate_constant;
	kf1=0.0012;
	real kb1 first_order_rate_constant;
	kb1=0.00076;
	real kf2 second_order_rate_constant;
	kf2=0.01;
	real kb2 first_order_rate_constant;
	kb2=0.1;
	real kf3 first_order_rate_constant;
	kf3=1;
	real kb3 first_order_rate_constant;
	kb3=0.01;
	real kf34 first_order_rate_constant;
	kf34=0.001;
	real kb34 first_order_rate_constant;
	kb34=0;
	real V4 flux;
	V4=62.5;
	real k4 nm_;
	k4=50;
	real kf5 second_order_rate_constant;
	kf5=0.1;
	real kb5 first_order_rate_constant;
	kb5=1;
	real kf6 first_order_rate_constant;
	kf6=20;
	real kb6 first_order_rate_constant;
	kb6=5;
	real kf7 first_order_rate_constant;
	kf7=60;
	real kb7 first_order_rate_constant;
	kb7=546;
	real kf8 first_order_rate_constant;
	kf8=2040;
	real kb8 second_order_rate_constant;
	kb8=15700;
	real kf9 first_order_rate_constant;
	kf9=40.8;
	real kb9 second_order_rate_constant;
	kb9=0;
	real V10 flux;
	V10=0.0154;
	real k10 nm_;
	k10=340;
	real kf23 second_order_rate_constant;
	kf23=0.1;
	real kb23 first_order_rate_constant;
	kb23=2;
	real kf24 first_order_rate_constant;
	kf24=9.85;
	real kb24 first_order_rate_constant;
	kb24=0.0985;
	real kf25 first_order_rate_constant;
	kf25=45.8;
	real kb25 second_order_rate_constant;
	kb25=0.047;
	real V26 flux;
	V26=2620;
	real k26 nm_;
	k26=3680;
	real R(time) nm_;
	when(time=time.min) R=80;
	real Shc(time) nm_;
	when(time=time.min) Shc=1000;
	real PI3K(time) nm_;
	when(time=time.min) PI3K=10;
	real HRG(time) nm_;
	when(time=time.min) HRG=10;
	real R_HRG(time) nm_;
	when(time=time.min) R_HRG=0;
	real R_HRG2(time) nm_;
	when(time=time.min) R_HRG2=0;
	real Internalisation(time) nm_;
	when(time=time.min) Internalisation=0;
	real RP(time) nm_;
	when(time=time.min) RP=0;
	real R_Shc(time) nm_;
	when(time=time.min) R_Shc=0;
	real R_ShP(time) nm_;
	when(time=time.min) R_ShP=0;
	real ShP(time) nm_;
	when(time=time.min) ShP=0;
	real R_ShGS(time) nm_;
	when(time=time.min) R_ShGS=0;
	real ShGS(time) nm_;
	when(time=time.min) ShGS=0;
	real GS(time) nm_;
	when(time=time.min) GS=10;
	real R_PI3K(time) nm_;
	when(time=time.min) R_PI3K=0;
	real R_PI3Kstar(time) nm_;
	when(time=time.min) R_PI3Kstar=0;
	real PI3Kstar(time) nm_;
	when(time=time.min) PI3Kstar=0;
	real two dimensionless;
	two=2;
	real RasGTP(time) nm_;
	when(time=time.min) RasGTP=0;
	real kf11 first_order_rate_constant;
	kf11=0.222;
	real k11 nm_;
	k11=0.181;
	real V12 flux;
	V12=0.289;
	real k12 nm_;
	k12=0.0571;
	real RasGDP(time) nm_;
	when(time=time.min) RasGDP=120;
	real Akt_PIPP(time) nm_;
	when(time=time.min) Akt_PIPP=0.0;
	real RAF_star(time) nm_;
	when(time=time.min) RAF_star=100;
	real kf13 first_order_rate_constant;
	kf13=1.53;
	real k13 nm_;
	k13=11.7;
	real kf14 first_order_rate_constant;
	kf14=0.00673;
	real k14 nm_;
	k14=8.07;
	real E nm_;
	E=7;
	real RAF(time) nm_;
	when(time=time.min) RAF=0;
	real MEKP(time) nm_;
	when(time=time.min) MEKP=0;
	real MEKPP(time) nm_;
	when(time=time.min) MEKPP=0;
	real kf27 first_order_rate_constant;
	kf27=16.9;
	real k27 nm_;
	k27=39.1;
	real V28 flux;
	V28=17000;
	real k28 nm_;
	k28=9.02;
	real kf29 second_order_rate_constant;
	kf29=507;
	real kb29 first_order_rate_constant;
	kb29=234;
	real V30 flux;
	V30=20000;
	real k30 nm_;
	k30=80000;
	real kf31 first_order_rate_constant;
	kf31=0.107;
	real Akt.k31 nm_;
	Akt.k31=4.35;
	real V32 flux;
	V32=20000;
	real k32 nm_;
	k32=80000;
	real kf33 first_order_rate_constant;
	kf33=0.211;
	real Akt.k33 nm_;
	Akt.k33=12;
	real Akt.k16 nm_;
	Akt.k16=2200;
	real Akt.k18 nm_;
	Akt.k18=60;
	real P(time) nm_;
	when(time=time.min) P=800;
	real PIP3(time) nm_;
	when(time=time.min) PIP3=0;
	real Akt(time) nm_;
	when(time=time.min) Akt=10;
	real Akt_PIP3(time) nm_;
	when(time=time.min) Akt_PIP3=0;
	real Akt_PIP(time) nm_;
	when(time=time.min) Akt_PIP=0;
	real Akt.PP2A nm_;
	Akt.PP2A=11.4;
	real Akt.one dimensionless;
	Akt.one=1;
	real MEK.PP2A nm_;
	MEK.PP2A=11.4;
	real MEK(time) nm_;
	when(time=time.min) MEK=120;
	real kf15 first_order_rate_constant;
	kf15=3.5;
	real k15 nm_;
	k15=317;
	real kf16 first_order_rate_constant;
	kf16=0.058;
	real MEK.k16 nm_;
	MEK.k16=2200;
	real kf17 first_order_rate_constant;
	kf17=2.9;
	real k17 nm_;
	k17=317;
	real kf18 first_order_rate_constant;
	kf18=0.058;
	real MEK.k18 nm_;
	MEK.k18=60;
	real MEK.k31 nm_;
	MEK.k31=4.35;
	real MEK.k33 nm_;
	MEK.k33=12;
	real MEK.one dimensionless;
	MEK.one=1;
	real MKP3 nm_;
	MKP3=2.4;
	real ERK(time) nm_;
	when(time=time.min) ERK=1000;
	real ERKP(time) nm_;
	when(time=time.min) ERKP=0;
	real ERKPP(time) nm_;
	when(time=time.min) ERKPP=0;
	real kf19 first_order_rate_constant;
	kf19=9.5;
	real k19 nm_;
	k19=146000;
	real kf20 first_order_rate_constant;
	kf20=0.3;
	real k20 nm_;
	k20=160;
	real kf21 first_order_rate_constant;
	kf21=16;
	real k21 nm_;
	k21=146000;
	real kf22 first_order_rate_constant;
	kf22=0.27;
	real k22 nm_;
	k22=60;
	real ERK.one dimensionless;
	ERK.one=1;

	// <component name="environment">

	// <component name="PI3K">
	R:time=((-1)*kf1*R*HRG+kb1*R_HRG);
	HRG:time=((-1)*kf1*R*HRG+kb1*R_HRG);
	R_HRG:time=(kf1*R*HRG-kb1*R_HRG-two*(kf2*R_HRG*R_HRG-kb2*R_HRG2));
	R_HRG2:time=(kf2*R_HRG*R_HRG-kb2*R_HRG2-(kf3*R_HRG2-kb3*RP)+V4*RP/(k4+RP));
	RP:time=(kf3*R_HRG2-kb3*RP-V4*RP/(k4+RP)-(kf5*RP*Shc-kb5*R_Shc)+(kf8*R_ShP-kb8*ShGS*RP)-(kf23*RP*PI3K-kb23*R_PI3K)+(kf25*R_PI3Kstar-kb25*RP*PI3Kstar)-(kf34*RP-kb34*Internalisation));
	Internalisation:time=(kf34*RP-kb34*Internalisation);
	R_Shc:time=(kf5*RP*Shc-kb5*R_Shc-(kf6*R_Shc-kb6*R_ShP));
	Shc:time=((-1)*(kf5*RP*Shc-kb5*R_Shc)+V10*ShP/(k10+ShP));
	R_ShP:time=(kf6*R_Shc-kb6*R_ShP-(kf7*R_ShP-kb7*R_ShGS));
	GS:time=((-1)*(kf7*R_ShP-kb7*R_ShGS)+(kf9*ShGS-kb9*GS*ShP));
	ShP:time=(kf9*ShGS-kb9*GS*ShP-V10*ShP/(k10+ShP));
	R_ShGS:time=(kf7*R_ShP-kb7*R_ShGS-(kf8*R_ShP-kb8*ShGS*RP));
	ShGS:time=(kf8*R_ShP-kb8*ShGS*RP-(kf9*ShGS-kb9*GS*ShP));
	R_PI3K:time=(kf23*RP*PI3K-kb23*R_PI3K-(kf24*R_PI3K-kb24*R_PI3Kstar));
	PI3K:time=((-1)*(kf23*RP*PI3K-kb23*R_PI3K)+V26*PI3Kstar/(k26+PI3Kstar));
	R_PI3Kstar:time=(kf24*R_PI3K-kb24*R_PI3Kstar-(kf25*R_PI3Kstar-kb25*RP*PI3Kstar));
	PI3Kstar:time=(kf25*R_PI3Kstar-kb25*RP*PI3Kstar-V26*PI3Kstar/(k26+PI3Kstar));

	// <component name="RasGDPtoRasGTP">
	RasGDP:time=((-1)*(kf11*ShGS*RasGDP/(k11+RasGDP))+V12*RasGTP/(k12+RasGTP));
	RasGTP:time=(kf11*ShGS*RasGDP/(k11+RasGDP)-V12*RasGTP/(k12+RasGTP));

	// <component name="RAF">
	RAF:time=(kf14*(Akt_PIPP+E)*RAF_star/(k14+RAF_star)-kf13*RasGTP*RAF/(k13+RAF));
	RAF_star:time=((-1)*(kf14*(Akt_PIPP+E)*RAF_star/(k14+RAF_star))+kf13*RasGTP*RAF/(k13+RAF));

	// <component name="Akt">
	P:time=(V28*PIP3/(k28+PIP3)-kf27*PI3Kstar*P/(k27+P));
	PIP3:time=((-1)*(V28*PIP3/(k28+PIP3))+kf27*PI3Kstar*P/(k27+P)-(kf29*PIP3*Akt-kb29*Akt_PIP3));
	Akt:time=((-1)*(kf29*PIP3*Akt-kb29*Akt_PIP3));
	Akt_PIP3:time=(kf29*PIP3*Akt-kb29*Akt_PIP3-V30*Akt_PIP3/(k30*(Akt.one+Akt_PIP/k32)+Akt_PIP3)+kf31*Akt.PP2A*Akt_PIP/(Akt.k31*(Akt.one+MEKP/Akt.k16+MEKPP/Akt.k18+Akt_PIPP/Akt.k33)+Akt_PIP));
	Akt_PIP:time=(V30*Akt_PIP3/(k30*(Akt.one+Akt_PIP/k32)+Akt_PIP3)-kf31*Akt.PP2A*Akt_PIP/(Akt.k31*(Akt.one+MEKP/Akt.k16+MEKPP/Akt.k18+Akt_PIPP/Akt.k33)+Akt_PIP)-V32*Akt_PIP/(k32*(Akt.one+Akt_PIP3/k30)+Akt_PIP)+kf33*Akt.PP2A*Akt_PIPP/(Akt.k33*(Akt.one+MEKP/Akt.k16+MEKPP/Akt.k18+Akt_PIP/Akt.k31)+Akt_PIPP));
	Akt_PIPP:time=(V32*Akt_PIP/(k32*(Akt.one+Akt_PIP3/k30)+Akt_PIP)-kf33*Akt.PP2A*Akt_PIPP/(Akt.k33*(Akt.one+MEKP/Akt.k16+MEKPP/Akt.k18+Akt_PIP/Akt.k31)+Akt_PIPP));

	// <component name="MEK">
	MEK:time=((-1)*(kf15*RAF_star*MEK/(k15*(MEK.one+MEKP/k17)+MEK))+kf16*MEK.PP2A*MEKP/(MEK.k16*(MEK.one+MEKPP/MEK.k18+Akt_PIP/MEK.k31+Akt_PIPP/MEK.k33)+MEKP));
	MEKP:time=(kf15*RAF_star*MEK/(k15*(MEK.one+MEKP/k17)+MEK)-kf16*MEK.PP2A*MEKP/(MEK.k16*(MEK.one+MEKPP/MEK.k18+Akt_PIP/MEK.k31+Akt_PIPP/MEK.k33)+MEKP)-kf17*RAF_star*MEKP/(k17*(MEK.one+MEK/k15)+MEKP)+kf18*MEK.PP2A*MEKPP/(MEK.k18*(MEK.one+MEKP/MEK.k16+Akt_PIP/MEK.k31+Akt_PIPP/MEK.k33)+MEKPP));
	MEKPP:time=(kf17*RAF_star*MEKP/(k17*(MEK.one+MEK/k15)+MEKP)-kf18*MEK.PP2A*MEKPP/(MEK.k18*(MEK.one+MEKP/MEK.k16+Akt_PIP/MEK.k31+Akt_PIPP/MEK.k33)+MEKPP));

	// <component name="ERK">
	ERK:time=((-1)*(kf19*MEKPP*ERK/(k19*(ERK.one+ERKP/k21)+ERK))+kf20*MKP3*ERKP/(k20*(ERK.one+ERKPP/k22)+ERKP));
	ERKPP:time=(kf21*MEKPP*ERKP/(k21*(ERK.one+ERK/k19)+ERKP)-kf22*MKP3*ERKPP/(k22*(ERK.one+ERKP/k20)+ERKPP));
	ERKP:time=(kf19*MEKPP*ERK/(k19*(ERK.one+ERKP/k21)+ERK)-kf20*MKP3*ERKP/(k20*(ERK.one+ERKPP/k22)+ERKP)-kf21*MEKPP*ERKP/(k21*(ERK.one+ERK/k19)+ERKP)+kf22*MKP3*ERKPP/(k22*(ERK.one+ERKP/k20)+ERKPP));
}