#########################################################################
## multiscaleOneCompt.py ---
## This program is part of 'MOOSE', the
## Messaging Object Oriented Simulation Environment.
## Copyright (C) 2014 Upinder S. Bhalla. and NCBS
## It is made available under the terms of the
## GNU General Public License version 2 or later.
## See the file COPYING.LIB for the full notice.
#########################################################################
from __future__ import print_function
import sys
import os
os.environ['NUMPTHREADS'] = '1'
import math
import numpy
import pylab
import matplotlib.pyplot as plt
import moose
import proto18
scriptDir = os.path.dirname( os.path.realpath( __file__ ) )
#EREST_ACT = -70e-3
def loadElec():
library = moose.Neutral( '/library' )
moose.setCwe( '/library' )
proto18.make_Ca()
proto18.make_Ca_conc()
proto18.make_Na()
proto18.make_K_DR()
proto18.make_K_A()
# Disable all the prototypes.
for x in moose.wildcardFind( "/library/##" ):
x.tick = -1
model = moose.Neutral( '/model' )
cellId = moose.loadModel(
os.path.join( scriptDir, 'soma.p')
, '/model/elec', "Neutral"
)
moose.setCwe( '/' )
return cellId
def loadChem():
chem = moose.Neutral( '/model/chem' )
modelId = moose.loadModel(
os.path.join( scriptDir, '..', 'genesis', 'chanPhosphByCaMKII.g' )
, '/model/chem', 'ee'
)
stoich = moose.Stoich( '/model/chem/kinetics/stoich' )
ksolve = moose.Ksolve( '/model/chem/kinetics/ksolve' )
stoich.compartment = moose.element( '/model/chem/kinetics' )
stoich.ksolve = ksolve
stoich.path = "/model/chem/##"
def makeModel():
loadElec()
loadChem()
makeAdaptors()
def makeAdaptors():
##################################################################
# set up adaptor for elec model Ca -> chem model Ca
# Here it is easy because we don't have to deal with different
# sizes of electrical and chemical compartments.
adaptCa = moose.Adaptor( '/model/chem/kinetics/adaptCa' )
chemCa = moose.element( '/model/chem/kinetics/Ca' )
elecCa = moose.element( '/model/elec/soma/Ca_conc' )
moose.connect( elecCa, 'concOut', adaptCa, 'input' )
moose.connect( adaptCa, 'output', chemCa, 'setConc' )
adaptCa.inputOffset = 0.0 #
adaptCa.outputOffset = 0.00008 # 80 nM offset in chem.
adaptCa.scale = 0.0008
# set up adaptor for chem model chan -> elec model chan.
adaptChan = moose.Adaptor( '/model/chem/kinetics/adaptChan' )
chemChan = moose.element( '/model/chem/kinetics/chan' )
elecChan = moose.element( '/model/elec/soma/K_A' )
# The Adaptor has to request the output conc of the chemical pool,
# since there isn't an output message to deliver this value.
moose.connect( adaptChan, 'requestOut', chemChan, 'getConc' )
moose.connect( adaptChan, 'output', elecChan, 'setGbar' )
adaptChan.inputOffset = 0.0 #
adaptChan.outputOffset = 0.0
adaptChan.scale = 1e-5 #
def addPlot( objpath, field, plot, tick ):
if moose.exists( objpath ):
tab = moose.Table( '/graphs/' + plot )
obj = moose.element( objpath )
moose.connect( tab, 'requestOut', obj, field )
tab.tick = tick
return tab
else:
print(("failed in addPlot(", objpath, field, plot, tick, ")"))
return 0
[docs]def main():
"""
This example builds a simple multiscale model involving
electrical and chemical signaling, but without spatial dimensions.
The electrical cell model is in a single compartment and has
voltage-gated channels, including a voltage-gated Ca channel for
Ca influx, and a K_A channel which is regulated by the chemical
pathways.
The chemical model has calcium activating Calmodulin which activates
CaM-Kinase II. The kinase phosphorylates the K_A channel to inactivate
it.
The net effect of the multiscale activity is a positive feedback
loop where activity increases Ca, which activates the kinase,
which reduces K_A, leading to increased excitability of the cell.
In this example this results
in a bistable neuron. In the resting state the cell does not fire,
but if it is activated by a current pulse the cell will continue to
fire even after the current is turned off. Application of an
inhibitory current restores the cell to its silent state.
Both the electrical and chemical models are loaded in from model
description files, and these files could be replaced if one wished
to define different models. However, there
are model-specific Adaptor objects needed to map activity between the
models of the two kinds. The Adaptors connect specific model entities
between the two models. Here one Adaptor connects the electrical
Ca_conc object to the chemical Ca pool. The other Adaptor connects
the chemical pool representing the K_A channel to its conductance
term in the electrical model.
"""
runtime = 4
elecDt = 50e-6
chemDt = 0.005
ePlotDt = 0.5e-3
cPlotDt = 0.0025
makeModel()
moose.setClock( 8, ePlotDt )
moose.setClock( 18, cPlotDt )
for i in range( 0, 10 ):
moose.setClock( i, elecDt )
for i in range( 10, 18 ):
moose.setClock( i, chemDt )
graphs = moose.Neutral( '/graphs' )
caplot = addPlot( '/model/elec/soma/Ca_conc', 'getCa', 'somaCa', 8 )
vmplot = addPlot( '/model/elec/soma', 'getVm', 'somaVm', 8 )
ikplot = addPlot( '/model/elec/soma/K_A', 'getIk', 'KAIk', 8 )
addPlot( '/model/chem/kinetics/chan', 'getConc', 'chan', 18 )
addPlot( '/model/chem/kinetics/Ca', 'getConc', 'Ca', 18 )
addPlot( '/model/chem/kinetics/CaM', 'getConc', 'CaM', 18 )
addPlot( '/model/chem/kinetics/Ca_CaM_CaMKII', 'getConc', 'enz', 18 )
hsolve = moose.HSolve( '/model/elec/hsolve' )
hsolve.dt = elecDt
hsolve.target = '/model/elec/soma'
moose.reinit()
moose.showfield( '/model/elec/soma' )
moose.element( '/model/elec/soma' ).inject = 0e-12
moose.start( runtime )
for i in moose.wildcardFind( '/model/elec/soma/#[ISA=ChanBase]' ):
print( "{} {} {}".format( i.name, i.Gbar, i.modulation ))
moose.element( '/model/elec/soma' ).inject = 1e-12
moose.start( runtime )
moose.element( '/model/elec/soma' ).inject = 0e-12
moose.start( runtime )
moose.element( '/model/elec/soma' ).inject = -1e-12
moose.start( runtime )
moose.element( '/model/elec/soma' ).inject = 0e-12
moose.start( runtime )
fig = plt.figure( figsize = (12,10) )
t = numpy.arange( 0, caplot.vector.size, 1 ) * caplot.dt
p1 = fig.add_subplot( 411 )
p1.plot( t, caplot.vector, label="Ca elec" )
p1.legend()
p2 = fig.add_subplot( 412 )
p2.plot( t, vmplot.vector, label="Vm" )
p2.legend()
p3 = fig.add_subplot( 413 )
p3.plot( t, ikplot.vector, label="Ik for K_A" )
p3.legend()
p4 = fig.add_subplot( 414 )
for x in moose.wildcardFind( '/graphs/#[FIELD(tick)=18]' ):
t = numpy.arange( 0, x.vector.size, 1 ) * x.dt
p4.plot( t, x.vector, label=x.name )
p4.legend()
plt.show()
quit()
if __name__ == '__main__':
main()