################################################################################
# #
# STEPS model of Cerebellar Long-Term Depression #
# #
# Project supervision: Erik De Schutter (2009-2016) #
# #
# Scripts Authors: #
# #
# Original development (2008-2012): Gabriela Antunes #
# Modified for STEPS 1.3 (2013): Iain Hepburn and Cory Simon #
# Addition of RKIP, improved calcium dynamics, and other #
# modifications (2014-2016): Anant Jain, Iain Hepburn, Himanshu Gangal #
# #
# #
################################################################################
# Example command-line usage:
# python LTD_STEPS_2.1.py ./data LTDsim 50 15 7 55 1 AMPAR Ca ERK MEK PKC
########################################################################
import sys
import os
###################### Command line arguments ##########################
# The path to store data
path=sys.argv[1]
# ID of the data (e.g. 'LTDsim')
dataid=sys.argv[2]
# Ca pulse peak (uM)
ca_peak=sys.argv[3]
# Ca pulse width (s)
ca_width=sys.argv[4]
# Raf number
raf_number = sys.argv[5]
# PKC number
pkc_number = sys.argv[6]
# Index of the data (to be converted to int, e.g. '7').
# This will also control the random number seed
idx = sys.argv[7]
# Data to save. Species names to be given as a subset of:
# 'AMPAR', 'Ca', 'AA', 'PKC', 'ERK', 'MEK', 'PLA2', 'Raf', 'RKIP'
# Activites of those species will be recorded over time
data_record = sys.argv[8:]
if len(data_record) == 0:
print "WARNING: Will not record any data"
datapath = path+'/'+dataid+'_'+ca_peak+'_'+ca_width+'_'+raf_number+'_'+pkc_number
# Create storage locations if necessary.
try: os.mkdir(path)
except: pass
try: os.mkdir(datapath)
except: pass
try: os.mkdir(datapath +'/'+ str(idx))
except: pass
# Automatically plot some figures or not. Basically plots activity
# over time for any species recorded
plotfigs = True
########################################################################
import datetime
import math
import numpy
from numpy import mean
from pylab import *
import pickle
import steps.model as smodel
import steps.solver as swmdirect
import steps.geom as swm
import steps.rng as srng
# Simualation data recording time step (s)
DT = 0.05
# Final time (s) simulated
INT = 3600
# Number of runs used to calculate the mean
# To see individual runs this should be set to 1
NITER = 1
# Initial time that will be discharged (s).
# The initial 10 minutes of all simulations should be discharged to allow the model to reach equilibrium.
tinit = 600
# Avogadro constant
Na = 6.02214129e23
########################################################################
"""
# Ca2+ pulse.
# The duration and amplitude of the pulse should be checked directly in the Ca2+ concentration
# The changes in [Ca2+] generated by the pulses are sensitive to stochasticity
# max_molar is the molar concentration per second injection (M/s)
# start defines the time the pulse starts (s)
# duration defines the duration of the pulse(s)
# start and duration are taken from command line arguments
"""
def square_pulse_ica(t,max_molar=float(ca_peak)*1e-6, duration=float(ca_width), start=1200.0):
if t < start or t > start + duration: return 0
else: return max_molar
########################################################################
def gen_model():
mdl=smodel.Model()
Ca = smodel.Spec('Ca', mdl)
PMCA = smodel.Spec('PMCA', mdl)
Ca1PMCA = smodel.Spec('Ca1PMCA', mdl)
NCX = smodel.Spec('NCX', mdl)
Ca1NCX = smodel.Spec('Ca1NCX', mdl)
Ca2NCX = smodel.Spec('Ca2NCX', mdl)
SERCA = smodel.Spec('SERCA', mdl)
Ca1SERCA = smodel.Spec('Ca1SERCA', mdl)
Ca2SERCA = smodel.Spec('Ca2SERCA', mdl)
PV = smodel.Spec('PV', mdl)
MgPV = smodel.Spec('MgPV', mdl)
CaPV = smodel.Spec('CaPV', mdl)
Mg2PV = smodel.Spec('Mg2PV', mdl)
Ca2PV = smodel.Spec('Ca2PV', mdl)
CBf = smodel.Spec('CBf', mdl)
CaCBf = smodel.Spec('CaCBf', mdl)
Ca2CBf = smodel.Spec('Ca2CBf', mdl)
CBs = smodel.Spec('CBs', mdl)
CaCBs = smodel.Spec('CaCBs', mdl)
Ca2CBs = smodel.Spec('Ca2CBs', mdl)
PKC = smodel.Spec('PKC', mdl)
Ca1PKC = smodel.Spec('Ca1PKC', mdl)
Ca3PKC = smodel.Spec('Ca3PKC', mdl)
AAPKC = smodel.Spec('AAPKC', mdl)
AACa1PKC = smodel.Spec('AACa1PKC', mdl)
AACa3PKC = smodel.Spec('AACa3PKC', mdl)
PKCstar = smodel.Spec('PKCstar', mdl)
PKCstar2 = smodel.Spec('PKCstar2', mdl)
PKCstar4 = smodel.Spec('PKCstar4', mdl)
PKCstar3 = smodel.Spec('PKCstar3', mdl)
Raf = smodel.Spec('Raf', mdl)
Rafstar = smodel.Spec('Rafstar', mdl)
PP5 = smodel.Spec('PP5', mdl)
PP5Rafstar = smodel.Spec('PP5Rafstar', mdl)
MEK = smodel.Spec('MEK', mdl)
RafstarMEK = smodel.Spec('RafstarMEK', mdl)
MEKp = smodel.Spec('MEKp', mdl)
PP2A = smodel.Spec('PP2A', mdl)
PP2AMEKp = smodel.Spec('PP2AMEKp', mdl)
RafstarMEKp = smodel.Spec('RafstarMEKp', mdl)
MEKstar = smodel.Spec('MEKstar', mdl)
PP2AMEKstar = smodel.Spec('PP2AMEKstar', mdl)
ERK = smodel.Spec('ERK', mdl)
MEKstarERK = smodel.Spec('MEKstarERK', mdl)
ERKp = smodel.Spec('ERKp', mdl)
MKP = smodel.Spec('MKP', mdl)
MKPERKp = smodel.Spec('MKPERKp', mdl)
MEKstarERKp = smodel.Spec('MEKstarERKp', mdl)
ERKstar = smodel.Spec('ERKstar', mdl)
MKPERKstar = smodel.Spec('MKPERKstar', mdl)
PLA2 = smodel.Spec('PLA2', mdl)
Ca1PLA2 = smodel.Spec('Ca1PLA2', mdl)
Ca2PLA2 = smodel.Spec('Ca2PLA2', mdl)
PLA2memb = smodel.Spec('PLA2memb', mdl)
Ca1PLA2memb = smodel.Spec('Ca1PLA2memb', mdl)
PLA2star1 = smodel.Spec('PLA2star1', mdl)
ERKstarPLA2 = smodel.Spec('ERKstarPLA2', mdl)
PLA2star2 = smodel.Spec('PLA2star2', mdl)
Ca1PLA2star2 = smodel.Spec('Ca1PLA2star2', mdl)
Ca2PLA2star2 = smodel.Spec('Ca2PLA2star2', mdl)
PLA2star2memb = smodel.Spec('PLA2star2memb', mdl)
Ca1PLA2star2memb = smodel.Spec('Ca1PLA2star2memb', mdl)
Ca2PLA2star2memb = smodel.Spec('Ca2PLA2star2memb', mdl)
PLA2star1APC = smodel.Spec('PLA2star1APC', mdl)
ERKstarCa1PLA2 = smodel.Spec('ERKstarCa1PLA2', mdl)
ERKstarCa2PLA2 = smodel.Spec('ERKstarCa2PLA2', mdl)
PP2ACa1PLA2star2 = smodel.Spec('PP2ACa1PLA2star2', mdl)
PP2ACa2PLA2star2 = smodel.Spec('PP2ACa2PLA2star2', mdl)
PP1Ca1PLA2star2 = smodel.Spec('PP1Ca1PLA2star2', mdl)
PP1Ca2PLA2star2 = smodel.Spec('PP1Ca2PLA2star2', mdl)
PLA2star2membAPC = smodel.Spec('PLA2star2membAPC', mdl)
Ca1PLA2star2membAPC = smodel.Spec('Ca1PLA2star2membAPC', mdl)
Ca2PLA2star2membAPC = smodel.Spec('Ca2PLA2star2membAPC', mdl)
PP1 = smodel.Spec('PP1', mdl)
PP1PLA2star2 = smodel.Spec('PP1PLA2star2', mdl)
PP2APLA2star2 = smodel.Spec('PP2APLA2star2', mdl)
AA = smodel.Spec('AA', mdl)
AMPAR = smodel.Spec('AMPAR', mdl)
PKCstarAMPAR = smodel.Spec('PKCstarAMPAR', mdl)
PKCstar2AMPAR = smodel.Spec('PKCstar2AMPAR', mdl)
PKCstar4AMPAR = smodel.Spec('PKCstar4AMPAR', mdl)
PKCstar3AMPAR = smodel.Spec('PKCstar3AMPAR', mdl)
PKCstarGRIPAMPAR = smodel.Spec('PKCstarGRIPAMPAR', mdl)
PKCstar2GRIPAMPAR = smodel.Spec('PKCstar2GRIPAMPAR', mdl)
PKCstar4GRIPAMPAR = smodel.Spec('PKCstar4GRIPAMPAR', mdl)
PKCstar3GRIPAMPAR = smodel.Spec('PKCstar3GRIPAMPAR', mdl)
AMPAR_P = smodel.Spec('AMPAR_P', mdl)
PP2AAMPAR_P = smodel.Spec('PP2AAMPAR_P', mdl)
PP2AGRIPAMPAR_P = smodel.Spec('PP2AGRIPAMPAR_P', mdl)
AMPARextra = smodel.Spec('AMPARextra', mdl)
AMPARextra_P = smodel.Spec('AMPARextra_P', mdl)
PP2AAMPARextra_P = smodel.Spec('PP2AAMPARextra_P', mdl)
GRIP = smodel.Spec('GRIP', mdl)
GRIPAMPAR = smodel.Spec('GRIPAMPAR', mdl)
GRIPAMPAR_P = smodel.Spec('GRIPAMPAR_P', mdl)
AMPARdend = smodel.Spec('AMPARdend', mdl)
AMPARdend_P = smodel.Spec('AMPARdend_P', mdl)
PP2AAMPARdend_P = smodel.Spec('PP2AAMPARdend_P', mdl)
AMPARcyt = smodel.Spec('AMPARcyt', mdl)
AMPARcyt_P = smodel.Spec('AMPARcyt_P', mdl)
PP2AAMPARcyt_P = smodel.Spec('PP2AAMPARcyt_P', mdl)
RKIP = smodel.Spec('RKIP', mdl)
RafRKIP = smodel.Spec('RafRKIP', mdl)
RKIPstar = smodel.Spec('RKIPstar', mdl)
RP = smodel.Spec('RP', mdl)
RKIPstarRP = smodel.Spec ('RKIPstarRP', mdl)
RafRKIPPKCstar = smodel.Spec('RafRKIPPKCstar', mdl)
RafRKIPPKCstar2= smodel.Spec ('RafRKIPPKCstar2', mdl)
RafRKIPPKCstar3= smodel.Spec ('RafRKIPPKCstar3', mdl)
RafRKIPPKCstar4= smodel.Spec ('RafRKIPPKCstar4', mdl)
MEKRKIP= smodel.Spec('MEKRKIP', mdl)
MEKRKIPPKCstar = smodel.Spec('MEKRKIPPKCstar', mdl)
MEKRKIPPKCstar2= smodel.Spec ('MEKRKIPPKCstar2', mdl)
MEKRKIPPKCstar3= smodel.Spec ('MEKRKIPPKCstar3', mdl)
MEKRKIPPKCstar4= smodel.Spec ('MEKRKIPPKCstar4', mdl)
extra = smodel.Volsys('extra', mdl)
vsys = smodel.Volsys('vsys', mdl)
s = smodel.Surfsys('memb', mdl)
ERs = smodel.Surfsys('ERmemb', mdl)
cytER = smodel.Volsys('cytER', mdl)
########################################################################
########################################################################
# Reactions
# Ca influx (Ca = calcium)
CainfluxR = smodel.Reac('Cainflux', vsys, rhs=[Ca])
# reaction constant for zero-order reaction units M/s in STEPS 1.2 and above
CainfluxR.kcst = 0.0
# Ca + PMCA <-> Ca1PMCA -> PMCA
Reac1 = smodel.SReac('pump2_f', s, ilhs = [Ca], slhs = [PMCA], srhs = [Ca1PMCA])
Reac1.kcst = 2500e6
Reac2 = smodel.SReac('Reac2', s, slhs = [Ca1PMCA], irhs = [Ca], srhs = [PMCA])
Reac2.kcst = 200
Reac3 = smodel.SReac('Reac3', s, slhs = [Ca1PMCA], srhs = [PMCA])
Reac3.kcst = 50
# Ca + NCX <-> Ca1NCX + Ca <-> Ca2NCX -> NCX
Reac4 = smodel.SReac('Reac4', s, ilhs = [Ca,Ca], slhs = [NCX], srhs = [Ca2NCX])
Reac4.kcst = 93.827e12
Reac5 = smodel.SReac('Reac5', s, slhs = [Ca2NCX], irhs = [Ca,Ca], srhs = [NCX])
Reac5.kcst = 4000.0
Reac6 = smodel.SReac('Reac6', s, slhs = [Ca2NCX], srhs = [NCX])
Reac6.kcst = 1000
# Ca + SERCA <-> Ca1SERCA +Ca <-> Ca2SERCA -> SERCA
Reac7 = smodel.SReac('Reac7', ERs, olhs = [Ca,Ca], slhs = [SERCA], srhs = [Ca2SERCA])
Reac7.kcst = 3428.7485714376226e12
Reac8 = smodel.SReac('Reac8', ERs, slhs = [Ca2SERCA], orhs = [Ca,Ca], srhs = [SERCA])
Reac8.kcst = 199.95577085780272
Reac9 = smodel.SReac('Reac9', ERs, slhs = [Ca2SERCA], srhs = [SERCA], irhs = [Ca,Ca])
Reac9.kcst = 50
# Leak
Reac10 = smodel.Reac('Reac10', vsys, rhs = [Ca])
# reaction constant for zero-order reaction units M/s in STEPS 1.2 and above
Reac10.kcst = 150/(Na*0.08e-15)
# PV + Ca <-> CaPV + Ca <-> Ca2PV
Reac15 = smodel.Reac('Reac15', vsys, lhs = [PV, Ca], rhs = [CaPV])
Reac15.kcst = 107e6
Reac16 = smodel.Reac('Reac16', vsys, lhs = [CaPV], rhs = [PV, Ca])
Reac16.kcst = 0.95
Reac17 = smodel.Reac('Reac17', vsys, lhs = [CaPV, Ca], rhs = [Ca2PV])
Reac17.kcst = 107e6
Reac18 = smodel.Reac('Reac18', vsys, lhs = [Ca2PV], rhs = [Ca, CaPV])
Reac18.kcst = 0.95
# PV + Mg <-> MgPV <-> Mg2PV, apparent rate constants, [Mg] = 590uM
Reac19 = smodel.Reac('Reac19', vsys, lhs = [PV], rhs = [MgPV])
Reac19.kcst = 472
Reac20 = smodel.Reac('Reac20', vsys, lhs = [MgPV], rhs = [PV])
Reac20.kcst = 25
Reac21 = smodel.Reac('Reac21', vsys, lhs = [MgPV], rhs = [Mg2PV])
Reac21.kcst = 472
Reac22 = smodel.Reac('Reac22', vsys, lhs = [Mg2PV], rhs = [MgPV])
Reac22.kcst = 25
# CBs + Ca <-> CaCBs + Ca <-> Ca2CBs, CBs - high affinity site
Reac23 = smodel.Reac('Reac23', vsys, lhs = [CBs, Ca], rhs = [CaCBs])
Reac23.kcst = 5.5e6
Reac24 = smodel.Reac('Reac24', vsys, lhs = [CaCBs], rhs = [CBs, Ca])
Reac24.kcst = 2.6
Reac25 = smodel.Reac('Reac25', vsys, lhs = [CaCBs, Ca], rhs = [Ca2CBs])
Reac25.kcst = 5.5e6
Reac26 = smodel.Reac('Reac26', vsys, lhs = [Ca2CBs], rhs = [CaCBs, Ca])
Reac26.kcst = 2.6
# CBf + Ca <-> CaCBf + Ca <-> Ca2CBf, CBf - medium affinity site
Reac27 = smodel.Reac('Reac27', vsys, lhs = [CBf, Ca], rhs = [CaCBf])
Reac27.kcst = 43.5e6
Reac28 = smodel.Reac('Reac28', vsys, lhs = [CaCBf], rhs = [CBf, Ca])
Reac28.kcst = 35.8
Reac29 = smodel.Reac('Reac29', vsys, lhs = [CaCBf, Ca], rhs = [Ca2CBf])
Reac29.kcst = 43.5e6
Reac30 = smodel.Reac('Reac30', vsys, lhs = [Ca2CBf], rhs = [CaCBf, Ca])
Reac30.kcst = 35.8
# PKC + Ca <-> CaPKC + 2Ca <-> Ca3PKC <-> Ca3PKC* (PKCstar)
Reac31 = smodel.Reac('Reac31', vsys, lhs = [PKC, Ca], rhs = [Ca1PKC])
Reac31.kcst = 13.3e6
Reac32 = smodel.Reac('Reac32', vsys, lhs = [Ca1PKC], rhs = [PKC, Ca])
Reac32.kcst = 12
Reac33 = smodel.Reac('Reac33', vsys, lhs = [Ca1PKC, Ca, Ca], rhs = [Ca3PKC])
Reac33.kcst = 1.e12
Reac34 = smodel.Reac('Reac34', vsys, lhs = [Ca3PKC], rhs = [Ca1PKC, Ca, Ca])
Reac34.kcst = 12
Reac35 = smodel.SReac('Reac35', s, ilhs = [Ca3PKC], srhs = [PKCstar])
Reac35.kcst = 11.3
Reac36 = smodel.SReac('Reac36', s, slhs = [PKCstar], irhs = [Ca3PKC])
Reac36.kcst = 0.23
# PKC + AA <-> AAPKC <-> AAPKC* (PKCstar2)
Reac37 = smodel.Reac('Reac37', vsys, lhs = [PKC, AA], rhs = [AAPKC])
Reac37.kcst = 1e6
Reac38 = smodel.Reac('Reac38', vsys, lhs = [AAPKC], rhs = [PKC, AA])
Reac38.kcst = 10
Reac39 = smodel.SReac('Reac39', s, ilhs = [AAPKC], srhs = [PKCstar2])
Reac39.kcst = 0.017
Reac40 = smodel.SReac('Reac40', s, slhs = [PKCstar2], irhs = [AAPKC])
Reac40.kcst = 0.0055
# Ca1PKC + AA <-> AACa1PKC <-> AACa1PKC* (PKCstar3) + 2Ca <-> AACa3PKC* (PKCstar4)
Reac41 = smodel.Reac('Reac41', vsys, lhs = [Ca1PKC, AA], rhs = [AACa1PKC])
Reac41.kcst = 1e6
Reac42 = smodel.Reac('Reac42', vsys, lhs = [AACa1PKC], rhs = [Ca1PKC, AA])
Reac42.kcst = 10
Reac43 = smodel.SReac('Reac43', s, ilhs = [AACa1PKC], srhs = [PKCstar3])
Reac43.kcst = 0.017
Reac44 = smodel.SReac('Reac44', s, slhs = [PKCstar3], irhs = [AACa1PKC])
Reac44.kcst = 0.0055
Reac45 = smodel.SReac('Reac45', s, slhs = [PKCstar3], ilhs = [Ca, Ca], srhs = [PKCstar4])
Reac45.kcst = 1.0e12
Reac46 = smodel.SReac('Reac46', s, slhs = [PKCstar4], srhs = [PKCstar3], irhs = [Ca, Ca])
Reac46.kcst = 12
# AAPKC + Ca <-> AACa1PKC <-> AACa1PKC + 2Ca <-> AACa3PKC <-> AACa3PKC* (PKCstar4)
Reac47 = smodel.Reac('Reac47', vsys, lhs = [AAPKC, Ca], rhs = [AACa1PKC])
Reac47.kcst = 13.3e6
Reac48 = smodel.Reac('Reac48', vsys, lhs = [AACa1PKC], rhs = [AAPKC, Ca])
Reac48.kcst = 12
Reac49 = smodel.Reac('Reac49', vsys, lhs = [AACa1PKC, Ca, Ca], rhs = [AACa3PKC])
Reac49.kcst = 1.0e12
Reac50 = smodel.Reac('Reac50', vsys, lhs = [AACa3PKC], rhs = [AACa1PKC, Ca, Ca])
Reac50.kcst = 12
Reac51 = smodel.SReac('Reac51', s, ilhs = [AACa3PKC], srhs = [PKCstar4])
Reac51.kcst = 11.3
Reac52 = smodel.SReac('Reac52', s, slhs = [PKCstar4], irhs = [AACa3PKC])
Reac52.kcst = 0.23
# Ca3PKC + AA <-> AACa3PKC
Reac53 = smodel.Reac('Reac53', vsys, lhs = [Ca3PKC, AA], rhs = [AACa3PKC])
Reac53.kcst = 1e6
Reac54 = smodel.Reac('Reac54', vsys, lhs = [AACa3PKC], rhs = [Ca3PKC, AA])
Reac54.kcst = 10
# AAPKC* + Ca <-> AACa1PKC*
Reac55 = smodel.SReac('Reac55', s, slhs = [PKCstar2], ilhs = [Ca], srhs = [PKCstar3])
Reac55.kcst = 13.3e6
Reac56 = smodel.SReac('Reac56', s, slhs = [PKCstar3], srhs = [PKCstar2], irhs = [Ca])
Reac56.kcst = 12
# Ca3PKC* + AA <-> AACa3PKC*
Reac57 = smodel.SReac('Reac57', s, ilhs = [AA], slhs = [PKCstar], srhs = [PKCstar4])
Reac57.kcst = 1e6
Reac58 = smodel.SReac('Reac58', s, slhs = [PKCstar4], irhs = [AA], srhs = [PKCstar])
Reac58.kcst = 10
# Ca3PKC* + RafRKIP <-> Ca3PKC*.RafRKIP -> Ca3PKC* + RafRKIP
Reac59 = smodel.SReac('Reac59', s, ilhs = [RafRKIP], slhs = [PKCstar], srhs = [RafRKIPPKCstar])
Reac59.kcst = 0.625e6
Reac60 = smodel.SReac('Reac60', s, slhs = [RafRKIPPKCstar], irhs = [RafRKIP], srhs = [PKCstar])
Reac60.kcst = 0.00245
Reac61 = smodel.SReac('Reac61', s, slhs = [RafRKIPPKCstar], srhs = [PKCstar], irhs = [RKIPstar, Rafstar])
Reac61.kcst = 0.0315
# AAPKC* + RafRKIP <-> AAPKC*.RafRKIP -> AAPKC* + RafRKIP
Reac62 = smodel.SReac('Reac62', s, ilhs = [RafRKIP], slhs = [PKCstar2], srhs = [RafRKIPPKCstar2])
Reac62.kcst = 0.625e6
Reac63 = smodel.SReac('Reac63', s, slhs = [RafRKIPPKCstar2], irhs = [RafRKIP], srhs = [PKCstar2])
Reac63.kcst = 0.00245
Reac64 = smodel.SReac('Reac64', s, slhs = [RafRKIPPKCstar2], srhs = [PKCstar2], irhs = [RKIPstar, Rafstar])
Reac64.kcst = 0.0315
# AACa1PKC* + RafRKIP -> AACa1PKC*.RafRKIP -> AACa1PKC* + RafRKIP
Reac65 = smodel.SReac('Reac65', s, ilhs = [RafRKIP], slhs = [PKCstar3], srhs = [RafRKIPPKCstar3])
Reac65.kcst = 0.625e6
Reac66 = smodel.SReac('Reac66', s, slhs = [RafRKIPPKCstar3], irhs = [RafRKIP], srhs = [PKCstar3])
Reac66.kcst = 0.00245
Reac67 = smodel.SReac('Reac67', s, slhs = [RafRKIPPKCstar3], srhs = [PKCstar3], irhs = [RKIPstar, Rafstar])
Reac67.kcst = 0.0315
# AACa3PKC* + RafRKIP <-> AACa3PKC*.RafRKIP -> AACa3PKC* + RafRKIP
Reac68 = smodel.SReac('Reac68', s, ilhs = [RafRKIP], slhs = [PKCstar4], srhs = [RafRKIPPKCstar4])
Reac68.kcst = 0.625e6
Reac69 = smodel.SReac('Reac69', s, slhs = [RafRKIPPKCstar4], irhs = [RafRKIP], srhs = [PKCstar4])
Reac69.kcst = 0.00245
Reac70 = smodel.SReac('Reac70', s, slhs = [RafRKIPPKCstar4], srhs = [PKCstar4], irhs = [RKIPstar, Rafstar])
Reac70.kcst = 0.0315
# Raf + RKIP --k1-> RafRKIP --k2-> Raf + RKIP
Reac71 = smodel.Reac('Reac71',vsys, lhs = [Raf , RKIP], rhs = [RafRKIP])
Reac71.kcst = 0.53e6
Reac72 = smodel.Reac('Reac72', vsys, lhs = [RafRKIP], rhs =[Raf , RKIP])
Reac72.kcst = 0.072
# Ca3PKC* + MEKRKIP <-> Ca3PKC*.MEKRKIP -> Ca3PKC* + MEKRKIP
Reac73 = smodel.SReac('Reac73', s, ilhs = [MEKRKIP], slhs = [PKCstar], srhs = [MEKRKIPPKCstar])
Reac73.kcst = 0.625e6
Reac74 = smodel.SReac('Reac74', s, slhs = [MEKRKIPPKCstar], irhs = [MEKRKIP], srhs = [PKCstar])
Reac74.kcst = 0.00245
Reac75 = smodel.SReac('Reac75', s, slhs = [MEKRKIPPKCstar], srhs = [PKCstar], irhs = [RKIPstar, MEKp])
Reac75.kcst = 0.0315
# AAPKC* + MEKRKIP <-> AAPKC*.MEKRKIP -> AAPKC* + MEKRKIP
Reac76 = smodel.SReac('Reac76', s, ilhs = [MEKRKIP], slhs = [PKCstar2], srhs = [MEKRKIPPKCstar2])
Reac76.kcst = 0.625e6
Reac77 = smodel.SReac('Reac77', s, slhs = [MEKRKIPPKCstar2], irhs = [MEKRKIP], srhs = [PKCstar2])
Reac77.kcst = 0.00245
Reac78 = smodel.SReac('Reac78', s, slhs = [MEKRKIPPKCstar2], srhs = [PKCstar2], irhs = [MEKp, RKIPstar])
Reac78.kcst = 0.0315
# AACa1PKC* + MEKRKIP -> AACa1PKC*.MEKRKIP -> AACa1PKC* + MEKRKIP
Reac79 = smodel.SReac('Reac79', s, ilhs = [MEKRKIP], slhs = [PKCstar3], srhs = [MEKRKIPPKCstar3])
Reac79.kcst = 0.625e6
Reac80 = smodel.SReac('Reac80', s, slhs = [MEKRKIPPKCstar3], irhs = [MEKRKIP], srhs = [PKCstar3])
Reac80.kcst = 0.00245
Reac81 = smodel.SReac('Reac81', s, slhs = [MEKRKIPPKCstar3], srhs = [PKCstar3], irhs = [MEKp, RKIPstar])
Reac81.kcst = 0.0315
# AACa3PKC* + MEKRKIP <-> AACa3PKC*.MEKRKIP -> AACa3PKC* + MEKRKIP
Reac82 = smodel.SReac('Reac82', s, ilhs = [MEKRKIP], slhs = [PKCstar4], srhs = [MEKRKIPPKCstar4])
Reac82.kcst = 0.625e6
Reac83 = smodel.SReac('Reac83', s, slhs = [MEKRKIPPKCstar4], irhs = [MEKRKIP], srhs = [PKCstar4])
Reac83.kcst = 0.00245
Reac84 = smodel.SReac('Reac84', s, slhs = [MEKRKIPPKCstar4], srhs = [PKCstar4], irhs = [RKIPstar, MEKp])
Reac84.kcst = 0.0315
# MEK + RKIP --k1-> MEKRKIP --k2-> MEK + RKIP
Reac85 = smodel.Reac('Reac85',vsys, lhs = [MEK , RKIP], rhs = [MEKRKIP])
Reac85.kcst = 0.53e6
Reac86 = smodel.Reac('Reac86', vsys, lhs = [MEKRKIP], rhs =[MEK , RKIP])
Reac86.kcst = 0.072
# RKIPstar + RP --k9-> RKIPstarRP --k10-> RKIPstar + RP
Reac87 = smodel.Reac('Reac87' , vsys, lhs =[RKIPstar, RP], rhs =[RKIPstarRP])
Reac87.kcst = 0.92e6
Reac88 = smodel.Reac('Reac88', vsys, lhs =[RKIPstarRP], rhs=[RKIPstar, RP])
Reac88.kcst = 0.00122
# RKIPstarRP --k11-> RKIP + RP
Reac89 = smodel.Reac('Reac89' , vsys, lhs =[RKIPstarRP], rhs = [RKIP, RP])
Reac89.kcst = 0.87
# PP5 (PP) + Raf* <-> PP5.Raf* -> PP5 + Raf
Reac90 = smodel.Reac('Reac90', vsys, lhs = [PP5, Rafstar], rhs = [PP5Rafstar])
Reac90.kcst = .55e6
Reac91 = smodel.Reac('Reac91', vsys, lhs = [PP5Rafstar], rhs = [PP5, Rafstar])
Reac91.kcst = 2
Reac92 = smodel.Reac('Reac92', vsys, lhs = [PP5Rafstar], rhs = [PP5, Raf])
Reac92.kcst = 0.5
# Raf* + MEK <-> Raf*.MEK -> Raf* + MEKp <->Raf*.MEKp -> Raf* + MEK* (MEKstar)
Reac93 = smodel.Reac('Reac93', vsys, lhs = [Rafstar, MEK], rhs = [RafstarMEK])
Reac93.kcst = 0.65e6
Reac94 = smodel.Reac('Reac94', vsys, lhs = [RafstarMEK], rhs = [Rafstar, MEK])
Reac94.kcst = 0.065
Reac95 = smodel.Reac('Reac95', vsys, lhs = [RafstarMEK], rhs = [Rafstar, MEKp])
Reac95.kcst = 1.0
Reac96 = smodel.Reac('Reac96', vsys, lhs = [Rafstar, MEKp], rhs = [RafstarMEKp])
Reac96.kcst = 0.65e6
Reac97 = smodel.Reac('Reac97', vsys, lhs = [RafstarMEKp], rhs = [Rafstar, MEKp])
Reac97.kcst = 0.065
Reac98 = smodel.Reac('Reac98', vsys, lhs = [RafstarMEKp], rhs = [Rafstar, MEKstar])
Reac98.kcst = 1.0
# PP2A + MEK* <-> PP2A.MEK* -> PP2A + MEKp <-> PP2A.MEKp -> PP2A + MEK
Reac99 = smodel.Reac('Reac99', vsys, lhs = [PP2A, MEKstar], rhs = [PP2AMEKstar])
Reac99.kcst = 0.75e6
Reac100 = smodel.Reac('Reac100', vsys, lhs = [PP2AMEKstar], rhs = [PP2A, MEKstar])
Reac100.kcst = 2
Reac101 = smodel.Reac('Reac101', vsys, lhs = [PP2AMEKstar], rhs = [PP2A, MEKp])
Reac101.kcst = 0.5
Reac102 = smodel.Reac('Reac102', vsys, lhs = [PP2A, MEKp], rhs = [PP2AMEKp])
Reac102.kcst = 0.75e6
Reac103 = smodel.Reac('Reac103', vsys, lhs = [PP2AMEKp], rhs = [PP2A, MEKp])
Reac103.kcst = 2
Reac104 = smodel.Reac('Reac104', vsys, lhs = [PP2AMEKp], rhs = [PP2A, MEK])
Reac104.kcst = 0.5
# MEK* + ERK <-> MEK*.ERK -> MEK* + ERKp <-> MEK*.ERKp -> MEK* + ERK* (ERKstar)
Reac105 = smodel.Reac('Reac105', vsys, lhs = [MEKstar, ERK], rhs = [MEKstarERK])
Reac105.kcst = 16.2e6
Reac106 = smodel.Reac('Reac106', vsys, lhs = [MEKstarERK], rhs = [MEKstar, ERK])
Reac106.kcst = 0.6
Reac107 = smodel.Reac('Reac107', vsys, lhs = [MEKstarERK], rhs = [MEKstar, ERKp])
Reac107.kcst = 0.15
Reac108 = smodel.Reac('Reac108', vsys, lhs = [MEKstar, ERKp], rhs = [MEKstarERKp])
Reac108.kcst = 16.2e6
Reac109 = smodel.Reac('Reac109', vsys, lhs = [MEKstarERKp], rhs = [MEKstar, ERKp])
Reac109.kcst = 0.6
Reac110 = smodel.Reac('Reac110', vsys, lhs = [MEKstarERKp], rhs = [MEKstar, ERKstar])
Reac110.kcst = 0.3
# MKP + ERK* <-> MKP.ERK* -> MKP + ERKp <-> MKP.ERKp -> MKP + ERK
Reac111 = smodel.Reac('Reac111', vsys, lhs = [MKP, ERKstar], rhs = [MKPERKstar])
Reac111.kcst = 13e6
Reac112 = smodel.Reac('Reac112', vsys, lhs = [MKPERKstar], rhs = [MKP, ERKstar])
Reac112.kcst = 0.396
Reac113 = smodel.Reac('Reac113', vsys, lhs = [MKPERKstar], rhs = [MKP, ERKp])
Reac113.kcst = 0.099
Reac114 = smodel.Reac('Reac114', vsys, lhs = [MKP, ERKp], rhs = [MKPERKp])
Reac114.kcst = 28e6
Reac115 = smodel.Reac('Reac115', vsys, lhs = [MKPERKp], rhs = [MKP, ERKp])
Reac115.kcst = 0.56
Reac116 = smodel.Reac('Reac116', vsys, lhs = [MKPERKp], rhs = [MKP, ERK])
Reac116.kcst = 0.14
# Ca + PLA2 <-> Ca1PLA2 + Ca <-> Ca2PLA2 <-> Ca2PLA2* (PLA2star1)
Reac117 = smodel.Reac('Reac117', vsys, lhs = [PLA2, Ca], rhs = [Ca1PLA2])
Reac117.kcst = 1.93e6
Reac118 = smodel.Reac('Reac118', vsys, lhs = [Ca1PLA2], rhs = [PLA2, Ca])
Reac118.kcst = 108
Reac119 = smodel.Reac('Reac119', vsys, lhs = [Ca, Ca1PLA2], rhs = [Ca2PLA2])
Reac119.kcst = 10.8e6
Reac120 = smodel.Reac('Reac120', vsys, lhs = [Ca2PLA2], rhs = [Ca1PLA2, Ca])
Reac120.kcst = 108
Reac121 = smodel.SReac('Reac121', s, ilhs = [Ca2PLA2], srhs = [PLA2star1])
Reac121.kcst = 300
Reac122 = smodel.SReac('Reac122', s, slhs = [PLA2star1], irhs = [Ca2PLA2])
Reac122.kcst = 15
# Ca1PLA2 <-> Ca1PLA2memb
Reac123 = smodel.SReac('Reac123', s, ilhs = [Ca1PLA2], srhs = [Ca1PLA2memb])
Reac123.kcst = 30
Reac124 = smodel.SReac('Reac124', s, slhs = [Ca1PLA2memb], irhs = [Ca1PLA2])
Reac124.kcst = 15
# PLA2 <-> PLA2memb
Reac125 = smodel.SReac('Reac125', s, ilhs = [PLA2], srhs = [PLA2memb])
Reac125.kcst = 3
Reac126 = smodel.SReac('Reac126', s, slhs = [PLA2memb], irhs = [PLA2])
Reac126.kcst = 15
# PLA2memb + Ca <-> Ca1PLA2memb + Ca <-> Ca2PLA2* (PLA2star1)
Reac127 = smodel.SReac('Reac127', s, slhs = [PLA2memb], ilhs = [Ca], srhs = [Ca1PLA2memb])
Reac127.kcst = 1.93e6
Reac128 = smodel.SReac('Reac128', s, ilhs = [Ca1PLA2memb], srhs = [PLA2memb], irhs = [Ca])
Reac128.kcst = 0.41
Reac129 = smodel.SReac('Reac129', s, slhs = [Ca1PLA2memb], ilhs = [Ca], srhs = [PLA2star1])
Reac129.kcst = 10.8e6
Reac130 = smodel.SReac('Reac130', s, ilhs = [PLA2star1], srhs = [Ca1PLA2memb], irhs = [Ca])
Reac130.kcst = 2.5
# PLA2star1 <-> PLA2star1APC -> PLA2star1 + AA
Reac131 = smodel.SReac('Reac131', s, slhs = [PLA2star1], srhs = [PLA2star1APC])
Reac131.kcst = 43
Reac132 = smodel.SReac('Reac132', s, slhs = [PLA2star1APC], srhs = [PLA2star1])
Reac132.kcst = 600
Reac133 = smodel.SReac('Reac133', s, slhs = [PLA2star1APC], irhs = [AA], srhs = [PLA2star1])
Reac133.kcst = 450
# ERK* + PLA2 <-> ERK*.PLA2 -> ERK* + PLA2p (PLA2star2)
Reac134 = smodel.Reac('Reac134', vsys, lhs = [ERKstar, PLA2], rhs = [ERKstarPLA2])
Reac134.kcst = 4e6
Reac135 = smodel.Reac('Reac135', vsys, lhs = [ERKstarPLA2], rhs = [ERKstar, PLA2])
Reac135.kcst = 1
Reac136 = smodel.Reac('Reac136', vsys, lhs = [ERKstarPLA2], rhs = [ERKstar, PLA2star2])
Reac136.kcst = 14
# ERK* + Ca1PLA2 <-> ERK*.Ca1PLA2 -> ERK* + Ca1PLA2p (Ca1PLA2star2)
Reac137 = smodel.Reac('Reac137', vsys, lhs = [ERKstar, Ca1PLA2], rhs = [ERKstarCa1PLA2])
Reac137.kcst = 4e6
Reac138 = smodel.Reac('Reac138', vsys, lhs = [ERKstarCa1PLA2], rhs = [ERKstar, Ca1PLA2])
Reac138.kcst = 1
Reac139 = smodel.Reac('Reac139', vsys, lhs = [ERKstarCa1PLA2], rhs = [ERKstar, Ca1PLA2star2])
Reac139.kcst = 14
# ERK* + Ca2PLA2 <-> ERK*.Ca2PLA2 -> ERK* + Ca2PLA2p (Ca2PLA2star2)
Reac140 = smodel.Reac('Reac140', vsys, lhs = [ERKstar, Ca2PLA2], rhs = [ERKstarCa2PLA2])
Reac140.kcst = 4e6
Reac141 = smodel.Reac('Reac141', vsys, lhs = [ERKstarCa2PLA2], rhs = [ERKstar, Ca2PLA2])
Reac141.kcst = 1
Reac142 = smodel.Reac('Reac142', vsys, lhs = [ERKstarCa2PLA2], rhs = [ERKstar, Ca2PLA2star2])
Reac142.kcst = 14
# PP2A + PLA2p <-> PP2A.PLA2p -> PP2A + PLA2
Reac143 = smodel.Reac('Reac143', vsys, lhs = [PP2A, PLA2star2], rhs = [PP2APLA2star2])
Reac143.kcst = 1.4e6
Reac144 = smodel.Reac('Reac144', vsys, lhs = [PP2APLA2star2], rhs = [PP2A, PLA2star2])
Reac144.kcst = 1.5
Reac145 = smodel.Reac('Reac145', vsys, lhs = [PP2APLA2star2], rhs = [PLA2, PP2A])
Reac145.kcst = 2.5
# PP2A + Ca1PLA2p <-> PP2A.Ca1PLA2p -> PP2A + Ca1PLA2
Reac146 = smodel.Reac('Reac146', vsys, lhs = [PP2A, Ca1PLA2star2], rhs = [PP2ACa1PLA2star2])
Reac146.kcst = 1.4e6
Reac147 = smodel.Reac('Reac147', vsys, lhs = [PP2ACa1PLA2star2], rhs = [PP2A, Ca1PLA2star2])
Reac147.kcst = 1.5
Reac148 = smodel.Reac('Reac148', vsys, lhs = [PP2ACa1PLA2star2], rhs = [PP2A, Ca1PLA2])
Reac148.kcst = 2.5
# PP2A + Ca2PLA2p <-> PP2A.Ca2PLA2p -> PP2A + Ca2PLA2
Reac149 = smodel.Reac('Reac149', vsys, lhs = [PP2A, Ca2PLA2star2], rhs = [PP2ACa2PLA2star2])
Reac149.kcst = 1.4e6
Reac150 = smodel.Reac('Reac150', vsys, lhs = [PP2ACa2PLA2star2], rhs = [PP2A, Ca2PLA2star2])
Reac150.kcst = 1.5
Reac151 = smodel.Reac('Reac151', vsys, lhs = [PP2ACa2PLA2star2], rhs = [PP2A, Ca2PLA2])
Reac151.kcst = 2.5
# PP1 + PLA2p <-> PP1.PLA2p -> PP1 + PLA2
Reac152 = smodel.Reac('Reac152', vsys, lhs = [PP1, PLA2star2], rhs = [PP1PLA2star2])
Reac152.kcst = 1.4e6
Reac153 = smodel.Reac('Reac153', vsys, lhs = [PP1PLA2star2], rhs = [PP1, PLA2star2])
Reac153.kcst = 1.5
Reac154 = smodel.Reac('Reac154', vsys, lhs = [PP1PLA2star2], rhs = [PLA2, PP1])
Reac154.kcst = 2.5
# PP1 + Ca1PLA2p <-> PP1.Ca1PLA2p -> PP1 + Ca1PLA2
Reac155 = smodel.Reac('Reac155', vsys, lhs = [PP1, Ca1PLA2star2], rhs = [PP1Ca1PLA2star2])
Reac155.kcst = 1.4e6
Reac156 = smodel.Reac('Reac156', vsys, lhs = [PP1Ca1PLA2star2], rhs = [PP1, Ca1PLA2star2])
Reac156.kcst = 1.5
Reac157 = smodel.Reac('Reac157', vsys, lhs = [PP1Ca1PLA2star2], rhs = [PP1, Ca1PLA2])
Reac157.kcst = 2.5
# PP1 + Ca2PLA2p <-> PP1.Ca2PLA2p -> PP1 + Ca2PLA2
Reac158 = smodel.Reac('Reac158', vsys, lhs = [PP1, Ca2PLA2star2], rhs = [PP1Ca2PLA2star2])
Reac158.kcst = 1.4e6
Reac159 = smodel.Reac('Reac159', vsys, lhs = [PP1Ca2PLA2star2], rhs = [PP1, Ca2PLA2star2])
Reac159.kcst = 1.5
Reac160 = smodel.Reac('Reac160', vsys, lhs = [PP1Ca2PLA2star2], rhs = [PP1, Ca2PLA2])
Reac160.kcst = 2.5
# PLA2p <-> PLA2** (PLA2star2memb) <-> PLA2**APC -> PLA2** + AA
Reac161 = smodel.SReac('Reac161', s, ilhs = [PLA2star2], srhs = [PLA2star2memb])
Reac161.kcst = 50
Reac162 = smodel.SReac('Reac162', s, slhs = [PLA2star2memb], irhs = [PLA2star2])
Reac162.kcst = 15
Reac163 = smodel.SReac('Reac163', s, slhs = [PLA2star2memb], srhs = [PLA2star2membAPC])
Reac163.kcst = 43
Reac164 = smodel.SReac('Reac164', s, slhs = [PLA2star2membAPC], srhs = [PLA2star2memb])
Reac164.kcst = 600
Reac165 = smodel.SReac('Reac165', s, slhs = [PLA2star2membAPC], srhs = [PLA2star2memb], irhs = [AA])
Reac165.kcst = 3600
# Ca1PLA2p <-> Ca1PLA2** (Ca1PLA2star2memb) <-> Ca1PLA2**APC -> Ca1PLA2** + AA
Reac166 = smodel.SReac('Reac166', s, ilhs = [Ca1PLA2star2], srhs = [Ca1PLA2star2memb])
Reac166.kcst = 30
Reac167 = smodel.SReac('Reac167', s, slhs = [Ca1PLA2star2memb], irhs = [Ca1PLA2star2])
Reac167.kcst = 15
Reac168 = smodel.SReac('Reac168', s, slhs = [Ca1PLA2star2memb], srhs = [Ca1PLA2star2membAPC])
Reac168.kcst = 43
Reac169 = smodel.SReac('Reac169', s, slhs = [Ca1PLA2star2membAPC], srhs = [Ca1PLA2star2memb])
Reac169.kcst = 600
Reac170 = smodel.SReac('Reac170', s, slhs = [Ca1PLA2star2membAPC], srhs = [Ca1PLA2star2memb], irhs = [AA])
Reac170.kcst = 3600
# Ca2PLA2p <-> Ca2PLA2** (Ca2PLA2star2memb) <-> Ca2PLA2**APC -> Ca2PLA2** + AA
Reac171 = smodel.SReac('Reac171', s, ilhs = [Ca2PLA2star2], srhs = [Ca2PLA2star2memb])
Reac171.kcst = 300
Reac172 = smodel.SReac('Reac172', s, slhs = [Ca2PLA2star2memb], irhs = [Ca2PLA2star2])
Reac172.kcst = 15
Reac173 = smodel.SReac('Reac173', s, slhs = [Ca2PLA2star2memb], srhs = [Ca2PLA2star2membAPC])
Reac173.kcst = 43
Reac174 = smodel.SReac('Reac174', s, slhs = [Ca2PLA2star2membAPC], srhs = [Ca2PLA2star2memb])
Reac174.kcst = 600
Reac175 = smodel.SReac('Reac175', s, slhs = [Ca2PLA2star2membAPC], srhs = [Ca2PLA2star2memb], irhs = [AA])
Reac175.kcst = 3600
# Ca + PLA2p <-> Ca1PLA2p
Reac176 = smodel.Reac('Reac176', vsys, lhs = [PLA2star2, Ca], rhs = [Ca1PLA2star2])
Reac176.kcst = 1.93e6
Reac177 = smodel.Reac('Reac177', vsys, lhs = [Ca1PLA2star2], rhs = [PLA2star2, Ca])
Reac177.kcst = 108
# Ca + PLA2** <-> Ca1PLA2**
Reac178 = smodel.SReac('Reac178', s, slhs = [PLA2star2memb], ilhs = [Ca], srhs = [Ca1PLA2star2memb])
Reac178.kcst = 1.93e6
Reac179 = smodel.SReac('Reac179', s, ilhs = [Ca1PLA2star2memb], srhs = [PLA2star2memb], irhs = [Ca])
Reac179.kcst = 0.41
# Ca + Ca1PLA2p <-> Ca2PLA2p
Reac180 = smodel.Reac('Reac180', vsys, lhs = [Ca, Ca1PLA2star2], rhs = [Ca2PLA2star2])
Reac180.kcst = 10.8e6
Reac181 = smodel.Reac('Reac181', vsys, lhs = [Ca2PLA2star2], rhs = [Ca1PLA2star2, Ca])
Reac181.kcst = 108
# Ca + Ca1PLA2** <-> Ca2PLA2**
Reac182 = smodel.SReac('Reac182', s, slhs = [Ca1PLA2star2memb], ilhs = [Ca], srhs = [Ca2PLA2star2memb])
Reac182.kcst = 10.8e6
Reac183 = smodel.SReac('Reac183', s, ilhs = [Ca2PLA2star2memb], srhs = [Ca1PLA2star2memb], irhs = [Ca])
Reac183.kcst = 2.5
# AA -> 0 ,degradation
Reac184 = smodel.SReac('Reac184', s, ilhs = [AA])
Reac184.kcst = .4
# Ca3PKC* + AMPAR <-> Ca3PKC*.AMPAR -> Ca3PKC* + AMPARp
Reac185 = smodel.SReac('Reac185', s, slhs = [PKCstar, AMPAR], srhs = [PKCstarAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac185.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac186 = smodel.SReac('Reac186', s, slhs = [PKCstarAMPAR], srhs = [PKCstar, AMPAR])
Reac186.kcst = 0.8
Reac187 = smodel.SReac('Reac187', s, slhs = [PKCstarAMPAR], srhs = [PKCstar, AMPAR_P])
Reac187.kcst = 0.3
# AAPKC* + AMPAR <-> AAPKC*.AMPAR -> AAPKC* + AMPARp
Reac188 = smodel.SReac('Reac188', s, slhs = [PKCstar2, AMPAR], srhs = [PKCstar2AMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac188.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac189 = smodel.SReac('Reac189', s, slhs = [PKCstar2AMPAR], srhs = [PKCstar2, AMPAR])
Reac189.kcst = 0.8
Reac190 = smodel.SReac('Reac190', s, slhs = [PKCstar2AMPAR], srhs = [PKCstar2, AMPAR_P])
Reac190.kcst = 0.3
# AACa1PKC* + AMPAR <-> AACa1PKC*.AMPAR -> AACa1PKC* + AMPARp
Reac191 = smodel.SReac('Reac191', s, slhs = [PKCstar3, AMPAR], srhs = [PKCstar3AMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac191.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac192 = smodel.SReac('Reac192', s, slhs = [PKCstar3AMPAR], srhs = [PKCstar3, AMPAR])
Reac192.kcst = 0.8
Reac193 = smodel.SReac('Reac193', s, slhs = [PKCstar3AMPAR], srhs = [PKCstar3, AMPAR_P])
Reac193.kcst = 0.3
# AACa3PKC* + AMPAR <-> AACa3PKC*.AMPAR -> AACa3PKC* + AMPARp
Reac194 = smodel.SReac('Reac194', s, slhs = [PKCstar4, AMPAR], srhs = [PKCstar4AMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac194.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac195 = smodel.SReac('Reac195', s, slhs = [PKCstar4AMPAR], srhs = [PKCstar4, AMPAR])
Reac195.kcst = 0.8
Reac196 = smodel.SReac('Reac196', s, slhs = [PKCstar4AMPAR], srhs = [PKCstar4, AMPAR_P])
Reac196.kcst = 0.3
# PP2A + AMPARp <-> PP2A.AMPARp -> PP2A + AMPAR
Reac197 = smodel.SReac('Reac197', s, ilhs = [PP2A], slhs = [AMPAR_P], srhs = [PP2AAMPAR_P])
Reac197.kcst = 0.6e6
Reac198 = smodel.SReac('Reac198', s, slhs = [PP2AAMPAR_P], irhs = [PP2A], srhs = [AMPAR_P])
Reac198.kcst = 0.17
Reac199 = smodel.SReac('Reac199', s, slhs = [PP2AAMPAR_P], srhs = [AMPAR], irhs = [PP2A])
Reac199.kcst = 0.25
# Ca3PKC* + GRIP.AMPA <-> Ca3PKC*.GRIP.AMPA -> Ca3PKC* + GRIP.AMPAp
Reac200 = smodel.SReac('Reac200', s, slhs = [PKCstar, GRIPAMPAR], srhs = [PKCstarGRIPAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac200.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac201 = smodel.SReac('Reac201', s, slhs = [PKCstarGRIPAMPAR], srhs = [PKCstar, GRIPAMPAR])
Reac201.kcst = 0.8
Reac202 = smodel.SReac('Reac202', s, slhs = [PKCstarGRIPAMPAR], srhs = [PKCstar, GRIPAMPAR_P])
Reac202.kcst = 0.3
# AAPKC* + GRIP.AMPA <-> AAPKC*.GRIP.AMPA -> AAPKC* + GRIP.AMPAp
Reac203 = smodel.SReac('Reac203', s, slhs = [PKCstar2, GRIPAMPAR], srhs = [PKCstar2GRIPAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac203.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac204 = smodel.SReac('Reac204', s, slhs = [PKCstar2GRIPAMPAR], srhs = [PKCstar2, GRIPAMPAR])
Reac204.kcst = 0.8
Reac205 = smodel.SReac('Reac205', s, slhs = [PKCstar2GRIPAMPAR], srhs = [PKCstar2, GRIPAMPAR_P])
Reac205.kcst = 0.3
# AACa1PKC* + GRIP.AMPA <-> AACa1PKC*.GRIP.AMPA -> AACa1PKC* + GRIP.AMPAp
Reac206 = smodel.SReac('Reac206', s, slhs = [PKCstar3, GRIPAMPAR], srhs = [PKCstar3GRIPAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac206.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac207 = smodel.SReac('Reac207', s, slhs = [PKCstar3GRIPAMPAR], srhs = [PKCstar3, GRIPAMPAR])
Reac207.kcst = 0.8
Reac208 = smodel.SReac('Reac208', s, slhs = [PKCstar3GRIPAMPAR], srhs = [PKCstar3, GRIPAMPAR_P])
Reac208.kcst = 0.3
# AACa3PKC* + GRIP.AMPA <-> AACa3PKC*.GRIP.AMPA -> AACa3PKC* + GRIP.AMPAp
Reac209 = smodel.SReac('Reac209', s, slhs = [PKCstar4, GRIPAMPAR], srhs = [PKCstar4GRIPAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac209.kcst = 0.4e6*(1.02e-12/1.87e-19)
Reac210 = smodel.SReac('Reac210', s, slhs = [PKCstar4GRIPAMPAR], srhs = [PKCstar4, GRIPAMPAR])
Reac210.kcst = 0.8
Reac211 = smodel.SReac('Reac211', s, slhs = [PKCstar4GRIPAMPAR], srhs = [PKCstar4, GRIPAMPAR_P])
Reac211.kcst = 0.3
# PP2A + GRIP.AMPAp <-> PP2A.GRIP.AMPAp -> PP2A + GRIP.AMPA
Reac212 = smodel.SReac('Reac212', s, ilhs = [PP2A], slhs = [GRIPAMPAR_P], srhs = [PP2AGRIPAMPAR_P])
Reac212.kcst = 0.6e6
Reac213 = smodel.SReac('Reac213', s, slhs = [PP2AGRIPAMPAR_P], irhs = [PP2A], srhs = [GRIPAMPAR_P])
Reac213.kcst = 0.17
Reac214 = smodel.SReac('Reac214', s, slhs = [PP2AGRIPAMPAR_P], srhs = [GRIPAMPAR], irhs = [PP2A])
Reac214.kcst = 0.25
# AMPAR + GRIP <-> GRIP.AMPA
Reac215 = smodel.SReac('Reac215', s, slhs = [AMPAR, GRIP], srhs = [GRIPAMPAR])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac215.kcst = 1e6*(1.02e-12/1.87e-19)
Reac216 = smodel.SReac('Reac216', s, slhs = [GRIPAMPAR], srhs = [GRIP, AMPAR])
Reac216.kcst = 7
# AMPARp + GRIP <-> GRIP.AMPAp (AMPAR = synaptic AMPAR)
Reac217 = smodel.SReac('Reac217', s, slhs = [AMPAR_P, GRIP], srhs = [GRIPAMPAR_P])
# units for 2D reaction based on square meters in STEPS 1.2 and above
# units for 2nd-order reaction are m^2/mol.s
Reac217.kcst = 1e6*(1.02e-12/1.87e-19)
Reac218 = smodel.SReac('Reac218', s, slhs = [GRIPAMPAR_P], srhs = [GRIP, AMPAR_P])
Reac218.kcst = 70
# AMPAR <-> AMPARextra (AMPARextra = extra-synaptic AMPAR)
Reac219 = smodel.SReac('Reac219', s, slhs = [AMPAR], srhs = [AMPARextra])
Reac219.kcst = 0.1
Reac220 = smodel.SReac('Reac220', s, slhs = [AMPARextra], srhs = [AMPAR])
Reac220.kcst = 0.02
# AMPARp <-> AMPARextra_P
Reac221 = smodel.SReac('Reac221', s, slhs = [AMPAR_P], srhs = [AMPARextra_P])
Reac221.kcst = 0.1
Reac222 = smodel.SReac('Reac222', s, slhs = [AMPARextra_P], srhs = [AMPAR_P])
Reac222.kcst = 0.02
# PP2A + AMPARextra_P <-> PP2A.AMPARextra_P -> PP2A + AMPARextra
Reac223 = smodel.SReac('Reac223', s, ilhs = [PP2A], slhs = [AMPARextra_P], srhs = [PP2AAMPARextra_P])
Reac223.kcst = 0.6e6
Reac224 = smodel.SReac('Reac224', s, slhs = [PP2AAMPARextra_P], irhs = [PP2A], srhs = [AMPARextra_P])
Reac224.kcst = 0.17
Reac225 = smodel.SReac('Reac225', s, slhs = [PP2AAMPARextra_P], srhs = [AMPARextra], irhs = [PP2A])
Reac225.kcst = 0.25
# AMPARextra <-> AMPARdend (AMPARdend = dendritic AMPAR)
Reac226 = smodel.SReac('Reac226', s, slhs = [AMPARextra], srhs = [AMPARdend])
Reac226.kcst = 0.02
Reac227 = smodel.SReac('Reac227', s, slhs = [AMPARdend], srhs = [AMPARextra])
Reac227.kcst = 0.00025
# AMPARextra_P <-> AMPARdend_P
Reac228 = smodel.SReac('Reac228', s, slhs = [AMPARextra_P], srhs = [AMPARdend_P])
Reac228.kcst = 0.02
Reac229 = smodel.SReac('Reac229', s, slhs = [AMPARdend_P], srhs = [AMPARextra_P])
Reac229.kcst = 0.00025
# PP2A + AMPARp_dend <-> PP2A.AMPARp_dend -> PP2A + AMPAR_dend
Reac230 = smodel.SReac('Reac230', s, ilhs = [PP2A], slhs = [AMPARdend_P], srhs = [PP2AAMPARdend_P])
Reac230.kcst = 0.6e6
Reac231 = smodel.SReac('Reac231', s, slhs = [PP2AAMPARdend_P], irhs = [PP2A], srhs = [AMPARdend_P])
Reac231.kcst = 0.17
Reac232 = smodel.SReac('Reac232', s, slhs = [PP2AAMPARdend_P], srhs = [AMPARdend], irhs = [PP2A])
Reac232.kcst = 0.25
# AMPARdend_P <-> AMPARcyt_P (AMPARcyt = cytosolic AMPAR)
Reac233 = smodel.SReac('Reac233', s, slhs = [AMPARdend_P], irhs = [AMPARcyt_P])
Reac233.kcst = 0.003
Reac234 = smodel.SReac('Reac234', s, ilhs = [AMPARcyt_P], srhs = [AMPARdend_P])
Reac234.kcst = 0.002
# PP2A + AMPARcyt_P <-> PP2A.AMPARcyt_P -> PP2A + AMPARcyt
Reac235 = smodel.Reac('Reac235', vsys, lhs = [PP2A, AMPARcyt_P], rhs = [PP2AAMPARcyt_P])
Reac235.kcst = 0.6e6
Reac236 = smodel.Reac('Reac236', vsys, lhs = [PP2AAMPARcyt_P], rhs = [PP2A, AMPARcyt_P])
Reac236.kcst = 0.17
Reac237 = smodel.Reac('Reac237', vsys, lhs = [PP2AAMPARcyt_P], rhs = [AMPARcyt, PP2A])
Reac237.kcst = 0.25
return mdl
########################################################################
# Geometric properties of the model: To change the size of the model,
# all volumetric compartments must be altered by the same ratio
# and all the areas must be scaled considering a spherical shape.
# Alterations in the size of the compartments will change automatically the population
# of species given in concentration (initial condition), but the species
# given in number of molecules must be altered manually by the same ratio
# to keep the balance among all the components of the model.
def gen_geom():
g = swm.Geom()
c = swm.Comp('vsys', g)
c.addVolsys('vsys')
c.vol = 0.08e-18
extra = swm.Comp('extra', g)
extra.addVolsys('extra')
extra.vol = 1.87e-22
s = swm.Patch('memb', g, icomp = c, ocomp = extra)
s.addSurfsys('memb')
s.area = 1.02e-12
cytER = swm.Comp('cytER', g)
cytER.addVolsys('cytER')
cytER.vol = 0.017e-18
ERs = swm.Patch('ERmemb', g, icomp = cytER, ocomp = c)
ERs.addSurfsys('ERmemb')
ERs.area = 0.32e-12
return g
######################################################################
def run_simulation():
rng=srng.create('mt19937',512)
m=gen_model()
g=gen_geom()
tpnts=numpy.arange(0.0,INT,DT)
ntpnts=tpnts.shape[0]
# initialize arrays for storage at each time point and at each run.
Ca = numpy.zeros((NITER, ntpnts))
PMCA = numpy.zeros((NITER, ntpnts))
Ca1PMCA = numpy.zeros((NITER, ntpnts))
NCX = numpy.zeros((NITER, ntpnts))
Ca1NCX = numpy.zeros((NITER, ntpnts))
Ca2NCX = numpy.zeros((NITER, ntpnts))
SERCA = numpy.zeros((NITER, ntpnts))
Ca1SERCA = numpy.zeros((NITER, ntpnts))
Ca2SERCA = numpy.zeros((NITER, ntpnts))
CBf = numpy.zeros((NITER, ntpnts))
CaCBf = numpy.zeros((NITER, ntpnts))
Ca2CBf = numpy.zeros((NITER, ntpnts))
CBs = numpy.zeros((NITER, ntpnts))
CaCBs = numpy.zeros((NITER, ntpnts))
Ca2CBs = numpy.zeros((NITER, ntpnts))
PV = numpy.zeros((NITER, ntpnts))
CaPV = numpy.zeros((NITER, ntpnts))
MgPV = numpy.zeros((NITER, ntpnts))
Ca2PV = numpy.zeros((NITER, ntpnts))
Mg2PV = numpy.zeros((NITER, ntpnts))
PP2A = numpy.zeros((NITER, ntpnts))
PKC = numpy.zeros((NITER, ntpnts))
Ca1PKC = numpy.zeros((NITER, ntpnts))
Ca3PKC = numpy.zeros((NITER, ntpnts))
AAPKC = numpy.zeros((NITER, ntpnts))
AACa1PKC = numpy.zeros((NITER, ntpnts))
AACa3PKC = numpy.zeros((NITER, ntpnts))
PKCstar = numpy.zeros((NITER, ntpnts))
PKCstar2 = numpy.zeros((NITER, ntpnts))
PKCstar4 = numpy.zeros((NITER, ntpnts))
PKCstar3 = numpy.zeros((NITER, ntpnts))
RKIPstarRP = numpy.zeros((NITER, ntpnts))
RafRKIPPKCstar = numpy.zeros((NITER, ntpnts))
RafRKIPPKCstar2 = numpy.zeros((NITER, ntpnts))
RafRKIPPKCstar3 = numpy.zeros((NITER, ntpnts))
RafRKIPPKCstar4 = numpy.zeros((NITER, ntpnts))
RafRKIP = numpy.zeros((NITER, ntpnts))
MEKRKIPPKCstar = numpy.zeros((NITER, ntpnts))
MEKRKIPPKCstar2 = numpy.zeros((NITER, ntpnts))
MEKRKIPPKCstar3 = numpy.zeros((NITER, ntpnts))
MEKRKIPPKCstar4 = numpy.zeros((NITER, ntpnts))
MEKRKIP = numpy.zeros((NITER, ntpnts))
Raf = numpy.zeros((NITER, ntpnts))
RKIPstar = numpy.zeros((NITER, ntpnts))
Rafstar = numpy.zeros((NITER, ntpnts))
PP5Rafstar = numpy.zeros((NITER, ntpnts))
PP5 = numpy.zeros((NITER, ntpnts))
MEK = numpy.zeros((NITER, ntpnts))
RafstarMEK = numpy.zeros((NITER, ntpnts))
MEKp = numpy.zeros((NITER, ntpnts))
PP2AMEKp = numpy.zeros((NITER, ntpnts))
RafstarMEKp = numpy.zeros((NITER, ntpnts))
MEKstar = numpy.zeros((NITER, ntpnts))
PP2AMEKstar = numpy.zeros((NITER, ntpnts))
ERK = numpy.zeros((NITER, ntpnts))
MEKstarERK = numpy.zeros((NITER, ntpnts))
ERKp = numpy.zeros((NITER, ntpnts))
MKP = numpy.zeros((NITER, ntpnts))
MKPERKp = numpy.zeros((NITER, ntpnts))
MEKstarERKp = numpy.zeros((NITER, ntpnts))
ERKstar = numpy.zeros((NITER, ntpnts))
MKPERKstar = numpy.zeros((NITER, ntpnts))
Ca1PLA2 = numpy.zeros((NITER, ntpnts))
Ca2PLA2 = numpy.zeros((NITER, ntpnts))
PLA2 = numpy.zeros((NITER, ntpnts))
PLA2memb = numpy.zeros((NITER, ntpnts))
Ca1PLA2memb = numpy.zeros((NITER, ntpnts))
PLA2star1 = numpy.zeros((NITER, ntpnts))
ERKstarPLA2 = numpy.zeros((NITER, ntpnts))
PLA2star2 = numpy.zeros((NITER, ntpnts))
AA = numpy.zeros((NITER, ntpnts))
Ca1PLA2star2 = numpy.zeros((NITER, ntpnts))
Ca2PLA2star2 = numpy.zeros((NITER, ntpnts))
PLA2star2memb = numpy.zeros((NITER, ntpnts))
Ca1PLA2star2memb = numpy.zeros((NITER, ntpnts))
Ca2PLA2star2memb = numpy.zeros((NITER, ntpnts))
PLA2star1APC = numpy.zeros((NITER, ntpnts))
ERKstarCa1PLA2 = numpy.zeros((NITER, ntpnts))
ERKstarCa2PLA2 = numpy.zeros((NITER, ntpnts))
PP2ACa1PLA2star2 = numpy.zeros((NITER, ntpnts))
PP2ACa2PLA2star2 = numpy.zeros((NITER, ntpnts))
PP1Ca1PLA2star2 = numpy.zeros((NITER, ntpnts))
PP1Ca2PLA2star2 = numpy.zeros((NITER, ntpnts))
PLA2star2membAPC = numpy.zeros((NITER, ntpnts))
Ca1PLA2star2membAPC = numpy.zeros((NITER, ntpnts))
Ca2PLA2star2membAPC = numpy.zeros((NITER, ntpnts))
PP1PLA2star2 = numpy.zeros((NITER, ntpnts))
PP1 = numpy.zeros((NITER, ntpnts))
PP2APLA2star2 = numpy.zeros((NITER, ntpnts))
PKCstarAMPAR = numpy.zeros((NITER, ntpnts))
PKCstar2AMPAR = numpy.zeros((NITER, ntpnts))
PKCstar4AMPAR = numpy.zeros((NITER, ntpnts))
PKCstar3AMPAR = numpy.zeros((NITER, ntpnts))
PKCstarGRIPAMPAR = numpy.zeros((NITER, ntpnts))
PKCstar2GRIPAMPAR = numpy.zeros((NITER, ntpnts))
PKCstar4GRIPAMPAR = numpy.zeros((NITER, ntpnts))
PKCstar3GRIPAMPAR = numpy.zeros((NITER, ntpnts))
PP2AGRIPAMPAR_P = numpy.zeros((NITER, ntpnts))
AMPAR = numpy.zeros((NITER, ntpnts))
AMPAR_P = numpy.zeros((NITER, ntpnts))
PP2AAMPAR_P = numpy.zeros((NITER, ntpnts))
AMPARextra = numpy.zeros((NITER, ntpnts))
AMPARextra_P = numpy.zeros((NITER, ntpnts))
PP2AAMPARextra_P = numpy.zeros((NITER, ntpnts))
GRIP = numpy.zeros((NITER, ntpnts))
GRIPAMPAR_P = numpy.zeros((NITER, ntpnts))
GRIPAMPAR = numpy.zeros((NITER, ntpnts))
AMPARdend = numpy.zeros((NITER, ntpnts))
AMPARdend_P = numpy.zeros((NITER, ntpnts))
PP2AAMPARdend_P = numpy.zeros((NITER, ntpnts))
AMPARcyt = numpy.zeros((NITER, ntpnts))
AMPARcyt_P = numpy.zeros((NITER, ntpnts))
PP2AAMPARcyt_P = numpy.zeros((NITER, ntpnts))
sim = swmdirect.Wmdirect(m, g, rng)
rng.initialize(int(idx)*100)
for j in arange(NITER):
print "Run number {0}".format(j)
sim.reset()
# initial conditions
sim.setCompConc('vsys', 'Ca', 0.045e-6)
sim.setCompConc('vsys', 'CBf', 37.775e-6)
sim.setCompConc('vsys', 'CaCBf', 2.1e-6)
sim.setCompConc('vsys', 'Ca2CBf', 0.125e-6)
sim.setCompConc('vsys', 'CBs', 36.25e-6)
sim.setCompConc('vsys', 'CaCBs', 3.4e-6)
sim.setCompConc('vsys', 'Ca2CBs', 0.125e-6)
sim.setCompConc('vsys', 'PV', 1.15e-6)
sim.setCompConc('vsys', 'CaPV', 8.4e-6)
sim.setCompConc('vsys', 'MgPV', 30.45e-6)
sim.setCompCount('vsys', 'PKC', int(pkc_number))
sim.setCompCount('vsys', 'Raf',int(raf_number))
sim.setCompConc('vsys', 'RKIP', 1.0e-6)
sim.setCompConc('vsys', 'MEK',1.5e-6)
sim.setCompConc('vsys', 'RP' ,3e-6)
sim.setCompConc('vsys', 'ERK',1.0e-6)
sim.setCompConc('vsys', 'MKP',0.26e-6)
sim.setCompConc('vsys', 'PP5', 1.0e-6)
sim.setCompCount('vsys', 'PP2A', 35)
sim.setCompConc('vsys', 'PLA2',0.4e-6)
sim.setCompCount('vsys', 'PP1', 30)
sim.setPatchCount('memb', 'PMCA', 10)
sim.setPatchCount('memb', 'NCX', 3)
sim.setPatchCount('ERmemb', 'SERCA', 5)
sim.setCompConc('cytER', 'Ca', 150e-6)
sim.setCompClamped('cytER', 'Ca', True) # clamped means the conc won't change as simulation runs.
sim.setPatchCount('memb', 'AMPAR', 3)
sim.setPatchCount('memb', 'AMPARextra', 16)
sim.setPatchCount('memb', 'GRIP', 22)
sim.setPatchCount('memb', 'GRIPAMPAR', 72)
sim.setPatchCount('memb', 'AMPARdend', 1600)
for i in range(ntpnts):
sim.run(tpnts[i])
if not i%1000: print "idx:", idx, " percent done", (tpnts[i]/INT*100)
Ca[j,i] = sim.getCompConc('vsys', 'Ca')
CBf[j,i] = sim.getCompConc('vsys', 'CBf')
CaCBf[j,i] = sim.getCompConc('vsys', 'CaCBf')
Ca2CBf[j,i] = sim.getCompConc('vsys', 'Ca2CBf')
CaCBs[j,i] = sim.getCompConc('vsys', 'CaCBs')
Ca2CBs[j,i] = sim.getCompConc('vsys', 'Ca2CBs')
PV[j,i] = sim.getCompConc('vsys', 'PV')
CaPV[j,i] = sim.getCompConc('vsys', 'CaPV')
MgPV[j,i] = sim.getCompConc('vsys', 'MgPV')
Ca2PV[j,i] = sim.getCompConc('vsys', 'Ca2PV')
Mg2PV[j,i] = sim.getCompConc('vsys', 'Mg2PV')
PP2A[j,i] = sim.getCompCount('vsys', 'PP2A')
PP5[j,i] = sim.getCompCount('vsys', 'PP5')
PKC[j,i] = sim.getCompCount('vsys', 'PKC')
Ca3PKC[j,i] = sim.getCompCount('vsys', 'Ca3PKC')
Ca1PKC[j,i] = sim.getCompCount('vsys', 'Ca1PKC')
AAPKC[j,i] = sim.getCompCount('vsys', 'AAPKC')
AACa1PKC[j,i] = sim.getCompCount('vsys', 'AACa1PKC')
AACa3PKC[j,i] = sim.getCompCount('vsys', 'AACa3PKC')
PKCstar[j,i] = sim.getPatchCount('memb', 'PKCstar')
PKCstar2[j,i] = sim.getPatchCount('memb', 'PKCstar2')
PKCstar3[j,i] = sim.getPatchCount('memb', 'PKCstar3')
PKCstar4[j,i] = sim.getPatchCount('memb', 'PKCstar4')
RafRKIP[j,i] = sim.getCompConc('vsys', 'RafRKIP')
RafRKIPPKCstar[j,i] = sim.getPatchCount('memb', 'RafRKIPPKCstar')
RafRKIPPKCstar2[j,i] = sim.getPatchCount('memb', 'RafRKIPPKCstar2')
RafRKIPPKCstar3[j,i] = sim.getPatchCount('memb', 'RafRKIPPKCstar3')
RafRKIPPKCstar4[j,i] = sim.getPatchCount('memb', 'RafRKIPPKCstar4')
MEKRKIP[j,i] = sim.getCompConc('vsys', 'MEKRKIP')
MEKRKIPPKCstar[j,i] = sim.getPatchCount('memb', 'MEKRKIPPKCstar')
MEKRKIPPKCstar2[j,i] = sim.getPatchCount('memb', 'MEKRKIPPKCstar2')
MEKRKIPPKCstar3[j,i] = sim.getPatchCount('memb', 'MEKRKIPPKCstar3')
MEKRKIPPKCstar4[j,i] = sim.getPatchCount('memb', 'MEKRKIPPKCstar4')
RKIPstar[j,i] = sim.getCompCount('vsys', 'RKIPstar')
RKIPstarRP[j,i] = sim.getCompCount('vsys', 'RKIPstarRP')
Raf[j,i] = sim.getCompCount('vsys', 'Raf')
Rafstar[j,i] = sim.getCompCount('vsys', 'Rafstar')
PP5Rafstar[j,i] = sim.getCompCount('vsys', 'PP5Rafstar')
MEK[j,i] = sim.getCompCount('vsys', 'MEK')
RafstarMEK[j,i] = sim.getCompCount('vsys', 'RafstarMEK')
MEKp[j,i] = sim.getCompCount('vsys', 'MEKp')
PP2AMEKp[j,i] = sim.getCompCount('vsys', 'PP2AMEKp')
RafstarMEKp[j,i] = sim.getCompCount('vsys', 'RafstarMEKp')
MEKstar[j,i] = sim.getCompCount('vsys', 'MEKstar')
ERK[j,i] = sim.getCompCount('vsys', 'ERK')
MEKstarERK[j,i] = sim.getCompCount('vsys', 'MEKstarERK')
ERKp[j,i] = sim.getCompCount('vsys', 'ERKp')
MKP[j,i] = sim.getCompCount('vsys', 'MKP')
MKPERKp[j,i] = sim.getCompCount('vsys', 'MKPERKp')
MEKstarERKp[j,i] = sim.getCompCount('vsys', 'MEKstarERKp')
ERKstar[j,i] = sim.getCompCount('vsys', 'ERKstar')
MKPERKstar[j,i] = sim.getCompCount('vsys', 'MKPERKstar')
PP2A[j,i] = sim.getCompCount('vsys', 'PP2A')
PLA2[j,i] = sim.getCompCount('vsys', 'PLA2')
Ca1PLA2[j,i] = sim.getCompCount('vsys', 'Ca1PLA2')
Ca2PLA2[j,i] = sim.getCompCount('vsys', 'Ca2PLA2')
PLA2memb[j,i] = sim.getPatchCount('memb', 'PLA2memb')
Ca1PLA2memb[j,i] = sim.getPatchCount('memb', 'Ca1PLA2memb')
PLA2star1[j,i] = sim.getPatchCount('memb', 'PLA2star1')
ERKstarPLA2[j,i] = sim.getCompCount('vsys', 'ERKstarPLA2')
PLA2star2[j,i] = sim.getCompCount('vsys', 'PLA2star2')
Ca1PLA2star2[j,i] = sim.getCompCount('vsys', 'Ca1PLA2star2')
Ca2PLA2star2[j,i] = sim.getCompCount('vsys', 'Ca2PLA2star2')
PLA2star2memb[j,i] = sim.getPatchCount('memb', 'PLA2star2memb')
Ca1PLA2star2memb[j,i] = sim.getPatchCount('memb', 'Ca1PLA2star2memb')
Ca2PLA2star2memb[j,i] = sim.getPatchCount('memb', 'Ca2PLA2star2memb')
AA[j,i] = sim.getCompCount('vsys', 'AA')
PP2APLA2star2[j,i] = sim.getCompCount('vsys', 'PP2APLA2star2')
PP1[j,i] = sim.getCompCount('vsys','PP1')
PP1PLA2star2[j,i] = sim.getCompCount('vsys', 'PP1PLA2star2')
PLA2star1APC[j,i] = sim.getPatchCount('memb', 'PLA2star1APC')
ERKstarCa1PLA2[j,i] = sim.getCompConc('vsys', 'ERKstarCa1PLA2')
ERKstarCa2PLA2[j,i] = sim.getCompConc('vsys', 'ERKstarCa2PLA2')
PP2ACa1PLA2star2[j,i] = sim.getCompConc('vsys', 'PP2ACa1PLA2star2')
PP2ACa2PLA2star2[j,i] = sim.getCompConc('vsys', 'PP2ACa2PLA2star2')
PP1Ca1PLA2star2[j,i] = sim.getCompConc('vsys', 'PP1Ca1PLA2star2')
PP1Ca2PLA2star2[j,i] = sim.getCompConc('vsys', 'PP1Ca2PLA2star2')
PLA2star2membAPC[j,i] = sim.getPatchCount('memb', 'PLA2star2membAPC')
Ca1PLA2star2membAPC[j,i] = sim.getPatchCount('memb', 'Ca1PLA2star2membAPC')
Ca2PLA2star2membAPC[j,i] = sim.getPatchCount('memb', 'Ca2PLA2star2membAPC')
Ca1PMCA[j,i] = sim.getPatchCount('memb', 'Ca1PMCA')
PMCA[j,i] = sim.getPatchCount('memb', 'PMCA')
NCX[j,i] = sim.getPatchCount('memb', 'NCX')
Ca2NCX[j,i] = sim.getPatchCount('memb', 'Ca2NCX')
SERCA[j,i] = sim.getPatchCount('ERmemb', 'SERCA')
Ca2SERCA[j,i] = sim.getPatchCount('ERmemb', 'Ca2SERCA')
PP1PLA2star2[j,i] = sim.getCompConc('vsys', 'PP1PLA2star2')
AMPAR[j,i] = sim.getPatchCount('memb', 'AMPAR')
AMPAR_P[j,i] = sim.getPatchCount('memb', 'AMPAR_P')
PP2AAMPAR_P[j,i] = sim.getPatchCount('memb', 'PP2AAMPAR_P')
AMPARextra[j,i] = sim.getPatchCount('memb', 'AMPARextra')
AMPARextra_P[j,i] = sim.getPatchCount('memb', 'AMPARextra_P')
PP2AAMPARextra_P[j,i] = sim.getPatchCount('memb', 'PP2AAMPARextra_P')
GRIP[j,i] = sim.getPatchCount('memb', 'GRIP')
GRIPAMPAR[j,i] = sim.getPatchCount('memb', 'GRIPAMPAR')
GRIPAMPAR_P[j,i] = sim.getPatchCount('memb', 'GRIPAMPAR_P')
PKCstarAMPAR[j,i] = sim.getPatchCount('memb', 'PKCstarAMPAR')
PKCstar2AMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar2AMPAR')
PKCstar3AMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar3AMPAR')
PKCstar4AMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar4AMPAR')
AMPARdend[j,i] = sim.getPatchCount('memb', 'AMPARdend')
AMPARdend_P[j,i] = sim.getPatchCount('memb', 'AMPARdend_P')
PP2AAMPARdend_P[j,i] = sim.getPatchCount('memb', 'PP2AAMPARdend_P')
AMPARcyt[j,i] = sim.getCompCount('vsys','AMPARcyt')
AMPARcyt_P[j,i] = sim.getCompCount('vsys', 'AMPARcyt_P')
PP2AAMPARcyt_P[j,i] = sim.getCompCount('vsys', 'PP2AAMPARcyt_P')
PKCstarGRIPAMPAR[j,i] = sim.getPatchCount('memb', 'PKCstarGRIPAMPAR')
PKCstar2GRIPAMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar2GRIPAMPAR')
PKCstar4GRIPAMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar4GRIPAMPAR')
PKCstar3GRIPAMPAR[j,i] = sim.getPatchCount('memb', 'PKCstar3GRIPAMPAR')
PP2AGRIPAMPAR_P[j,i] = sim.getPatchCount('memb', 'PP2AGRIPAMPAR_P')
sim.setCompReacK('vsys', 'Cainflux', square_pulse_ica(tpnts[i])) # changes reaction constant to be ...
#############################################################################
if 'AMPAR' in data_record:
AMPAR_mean = mean(AMPAR, axis=0)
AMPAR_P_mean = mean(AMPAR_P, axis=0)
PP2AAMPAR_P_mean = mean(PP2AAMPAR_P, axis=0)
GRIPAMPAR_mean = mean(GRIPAMPAR, axis=0)
GRIPAMPAR_P_mean = mean(GRIPAMPAR_P, axis=0)
PKCstarAMPAR_mean = mean(PKCstarAMPAR, axis=0)
PKCstar2AMPAR_mean = mean(PKCstar2AMPAR, axis=0)
PKCstar3AMPAR_mean = mean(PKCstar3AMPAR, axis=0)
PKCstar4AMPAR_mean = mean(PKCstar4AMPAR, axis=0)
PKCstarGRIPAMPAR_mean = mean(PKCstarGRIPAMPAR, axis=0)
PKCstar2GRIPAMPAR_mean = mean(PKCstar2GRIPAMPAR, axis=0)
PKCstar3GRIPAMPAR_mean = mean(PKCstar3GRIPAMPAR, axis=0)
PKCstar4GRIPAMPAR_mean = mean(PKCstar4GRIPAMPAR, axis=0)
PP2AGRIPAMPAR_P_mean = mean(PP2AGRIPAMPAR_P, axis=0)
AMPAR_syn_mean = AMPAR_mean + AMPAR_P_mean + PP2AAMPAR_P_mean + \
GRIPAMPAR_mean + GRIPAMPAR_P_mean + PKCstarAMPAR_mean + \
PKCstar2AMPAR_mean + PKCstar3AMPAR_mean + PKCstar4AMPAR_mean +\
PKCstarGRIPAMPAR_mean + PKCstar2GRIPAMPAR_mean + PKCstar3GRIPAMPAR_mean +\
PKCstar4GRIPAMPAR_mean + PP2AGRIPAMPAR_P_mean;
ofile_ampar=open(datapath +'/'+ str(idx) +'/AMPAR_data.txt','w')
pickle.dump(AMPAR_syn_mean, ofile_ampar)
ofile_ampar.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), AMPAR_syn_mean, 'g')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('Synaptic AMPAR (#)')
savefig(datapath +'/'+ str(idx) +'/Syn_AMPAR.png')
close()
if 'Ca' in data_record:
Ca_mean = mean(Ca, axis=0)
ofile_ca=open(datapath +'/'+ str(idx) +'/Ca_data.txt','w')
pickle.dump(Ca_mean, ofile_ca)
ofile_ca.close()
if plotfigs:
plot(linspace(0,INT,INT/DT),Ca_mean[:]*1e6,'k')
xlabel('Time (s)')
ylabel('Calcium (uM)')
xlim(1190, 1260)
savefig(datapath +'/'+ str(idx) +'/Ca_mean.png')
close()
if 'AA' in data_record:
AA_mean = mean(AA, axis=0)
AAPKC_mean = mean(AAPKC, axis=0)
AACa1PKC_mean = mean(AACa1PKC, axis=0)
AACa3PKC_mean = mean(AACa3PKC, axis=0)
AAA_mean = AA_mean+AAPKC_mean+AACa1PKC_mean+AACa3PKC_mean
ofile_aa=open(datapath +'/'+ str(idx) +'/AA_data.txt','w')
pickle.dump(AAA_mean, ofile_aa)
ofile_aa.close()
if plotfigs:
plot(linspace(0,INT,INT/DT),AAA_mean[:],label = 'AA')
xlabel('Time (s)')
ylabel('Arachidonic acid #')
savefig(datapath +'/'+ str(idx) +'/AA_mean.png')
close()
if 'PKC' in data_record:
PKCstar_mean = mean(PKCstar, axis=0)
PKCstar2_mean = mean(PKCstar2, axis=0)
PKCstar3_mean = mean(PKCstar3, axis=0)
PKCstar4_mean = mean(PKCstar4, axis=0)
RafRKIPPKCstar_mean = mean(RafRKIPPKCstar, axis=0)
RafRKIPPKCstar2_mean = mean(RafRKIPPKCstar2, axis=0)
RafRKIPPKCstar3_mean = mean(RafRKIPPKCstar3, axis=0)
RafRKIPPKCstar4_mean = mean(RafRKIPPKCstar4, axis=0)
PKCstarAMPAR_mean = mean(PKCstarAMPAR, axis=0)
PKCstar2AMPAR_mean = mean(PKCstar2AMPAR, axis=0)
PKCstar3AMPAR_mean = mean(PKCstar3AMPAR, axis=0)
PKCstar4AMPAR_mean = mean(PKCstar4AMPAR, axis=0)
PKCstarGRIPAMPAR_mean = mean(PKCstarGRIPAMPAR, axis=0)
PKCstar2GRIPAMPAR_mean = mean(PKCstar2GRIPAMPAR, axis=0)
PKCstar3GRIPAMPAR_mean = mean(PKCstar3GRIPAMPAR, axis=0)
PKCstar4GRIPAMPAR_mean = mean(PKCstar4GRIPAMPAR, axis=0)
MEKRKIPPKCstar_mean = mean(MEKRKIPPKCstar, axis=0)
MEKRKIPPKCstar2_mean = mean(MEKRKIPPKCstar2, axis=0)
MEKRKIPPKCstar3_mean = mean(MEKRKIPPKCstar3, axis=0)
MEKRKIPPKCstar4_mean = mean(MEKRKIPPKCstar4, axis=0)
PKC_mean = PKCstar_mean + PKCstar2_mean + PKCstar3_mean + PKCstar4_mean + \
RafRKIPPKCstar_mean + RafRKIPPKCstar2_mean + RafRKIPPKCstar3_mean + \
RafRKIPPKCstar4_mean + PKCstarAMPAR_mean + PKCstar2AMPAR_mean + \
PKCstar3AMPAR_mean + PKCstar4AMPAR_mean + \
MEKRKIPPKCstar_mean + MEKRKIPPKCstar2_mean + MEKRKIPPKCstar3_mean + \
MEKRKIPPKCstar4_mean + PKCstarGRIPAMPAR_mean + PKCstar2GRIPAMPAR_mean + \
PKCstar3GRIPAMPAR_mean + PKCstar4GRIPAMPAR_mean
ofile_pkc=open(datapath +'/'+ str(idx) +'/PKC_data.txt','w')
pickle.dump(PKC_mean, ofile_pkc)
ofile_pkc.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), PKC_mean, 'r')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('PKC (#)')
savefig(datapath +'/'+ str(idx) +'/PKC.png')
close()
if 'ERK' in data_record:
ERKstar_mean = mean(ERKstar, axis=0)
MKPERKstar_mean = mean(MKPERKstar, axis=0)
ERKstarPLA2_mean = mean(ERKstarPLA2, axis=0)
ERKstarCa1PLA2_mean = mean(ERKstarCa1PLA2, axis=0)
ERKstarCa2PLA2_mean = mean(ERKstarCa2PLA2, axis=0)
ERKp_mean = mean(ERKp, axis=0)
MKPERKp_mean = mean(MKPERKp, axis=0)
MEKstarERKp_mean = mean(MEKstarERKp, axis=0)
ERKstar_total = ERKstar_mean + MKPERKstar_mean + ERKstarPLA2_mean + ERKstarCa1PLA2_mean + ERKstarCa2PLA2_mean
ERKp_total = ERKp_mean + MKPERKp_mean + MEKstarERKp_mean
ERK_total = ERKstar_mean + MKPERKstar_mean + ERKstarPLA2_mean + ERKstarCa1PLA2_mean + ERKstarCa2PLA2_mean + \
ERKp_mean + MKPERKp_mean + MEKstarERKp_mean
ofile_erktotal=open(datapath +'/'+ str(idx) +'/ERKtotal_data.txt','w')
pickle.dump(ERK_total, ofile_erktotal)
ofile_erktotal.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), ERK_total, 'k', label ='total')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('ERK (#)')
savefig(datapath +'/'+ str(idx) +'/ERK.png')
close()
if 'MEK' in data_record:
MEKp_mean= mean(MEKp, axis=0)
PP2AMEKp_mean= mean(PP2AMEKp, axis=0)
RafstarMEKp_mean= mean(RafstarMEKp, axis=0)
MEKstar_mean= mean(MEKstar, axis=0)
MEKstarERK_mean= mean(MEKstarERK, axis=0)
MEKstarERKp_mean= mean(MEKstarERKp, axis=0)
PP2AMEKstar_mean = mean(PP2AMEKstar, axis=0)
MEK_total = MEKp_mean + PP2AMEKp_mean + RafstarMEKp_mean + \
MEKstar_mean + MEKstarERK_mean + MEKstarERKp_mean + PP2AMEKstar_mean
MEKstar_total = MEKstar_mean + MEKstarERK_mean + MEKstarERKp_mean + PP2AMEKstar_mean
MEKp_total = MEKp_mean + PP2AMEKp_mean + RafstarMEKp_mean
ofile_mektotal=open(datapath +'/'+ str(idx) +'/MEKtotal_data.txt','w')
pickle.dump(MEK_total, ofile_mektotal)
ofile_mektotal.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), MEK_total, 'k', label='total')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('MEK (#)')
savefig(datapath +'/'+ str(idx) +'/MEK.png')
close()
if 'PLA2' in data_record:
PLA2star1_mean = mean(PLA2star1, axis=0)
PLA2star1APC_mean = mean(PLA2star1APC, axis=0)
PLA2star2_mean = mean(PLA2star2, axis=0)
PLA2star2memb_mean = mean(PLA2star2memb, axis=0)
PLA2star2membAPC_mean = mean(PLA2star2membAPC, axis=0)
Ca1PLA2star2_mean = mean(Ca1PLA2star2, axis=0)
Ca2PLA2star2_mean = mean(Ca2PLA2star2, axis=0)
Ca1PLA2star2memb_mean = mean(Ca1PLA2star2memb, axis=0)
Ca2PLA2star2memb_mean = mean(Ca2PLA2star2memb, axis=0)
Ca1PLA2star2membAPC_mean = mean(Ca1PLA2star2membAPC, axis=0)
Ca2PLA2star2membAPC_mean = mean(Ca2PLA2star2membAPC, axis=0)
PP1PLA2star2_mean = mean(PP1PLA2star2, axis=0)
PP1Ca1PLA2star2_mean = mean(PP1Ca1PLA2star2, axis=0)
PP1Ca2PLA2star2_mean = mean(PP1Ca2PLA2star2, axis=0)
PP2APLA2star2_mean = mean(PP1PLA2star2, axis=0)
PP2ACa1PLA2star2_mean = mean(PP2ACa1PLA2star2, axis=0)
PP2ACa2PLA2star2_mean = mean(PP2ACa2PLA2star2, axis=0)
PLA2_mean = PLA2star1_mean + PLA2star1APC_mean + PLA2star2_mean + \
PLA2star2memb_mean + PLA2star2membAPC_mean + \
Ca1PLA2star2_mean + Ca2PLA2star2_mean + Ca1PLA2star2memb_mean + \
Ca2PLA2star2memb_mean + Ca1PLA2star2membAPC_mean + Ca2PLA2star2membAPC_mean + \
PP1PLA2star2_mean + PP1Ca1PLA2star2_mean + PP1Ca2PLA2star2_mean + \
PP2APLA2star2_mean + PP2ACa1PLA2star2_mean + PP2ACa2PLA2star2_mean
ofile_pla2=open(datapath +'/'+ str(idx) +'/PLA2_data.txt','w')
pickle.dump(PLA2_mean, ofile_pla2)
ofile_pla2.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), PLA2_mean, 'y')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('cPLA2 (#)')
savefig(datapath +'/'+ str(idx) +'/PLA2.png')
close()
if 'Raf' in data_record:
Rafstar_mean= mean(Rafstar, axis=0)
PP5Rafstar_mean= mean(PP5Rafstar, axis=0)
RafstarMEK_mean= mean(RafstarMEK, axis=0)
RafstarMEKp_mean= mean(RafstarMEKp, axis=0)
Raf_mean = Rafstar_mean + PP5Rafstar_mean + RafstarMEK_mean + RafstarMEKp_mean;
ofile_raf=open(datapath +'/'+ str(idx) +'/Raf_data.txt','w')
pickle.dump(Raf_mean, ofile_raf)
ofile_raf.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), Raf_mean, 'm')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('Raf (#)')
savefig(datapath +'/'+ str(idx) +'/Raf.png')
close()
if 'RKIP' in data_record:
RKIPstar_mean= mean(RKIPstar, axis=0)
RKIPstarRP_mean= mean(RKIPstarRP, axis=0)
RKIP_mean = RKIPstar_mean + RKIPstarRP_mean
ofile_rkip=open(datapath +'/'+ str(idx) +'/RKIP_data.txt','w')
pickle.dump(RKIP_mean, ofile_rkip)
ofile_rkip.close()
if plotfigs:
plot(linspace(0,INT,INT/DT), RKIP_mean, 'r')
ylim(ymin=0.0)
xlabel('Time (s)')
ylabel('RKIP (#)')
savefig(datapath +'/'+ str(idx) +'/RKIP.png')
close()
run_simulation()
