Package PyDSTool :: Package Toolbox :: Module neuralcomp
[hide private]
[frames] | no frames]

Source Code for Module PyDSTool.Toolbox.neuralcomp

   1  """
 
   2      An example set of basic compartmental model ModelSpec classes
 
   3      for use in computational neuroscience modeling.
 
   4  
 
   5      Basic classes provided:
 
   6      channel, compartment
 
   7  
 
   8      For multi-compartment models:
 
   9      soma, dendr_compartment, synapse, neuron, network, exc_synapse, inh_synapse
 
  10  
 
  11  Rob Clewley, September 2005.
 
  12  """ 
  13  from __future__ import division 
  14  
 
  15  from PyDSTool import * 
  16  from copy import copy 
  17  
 
  18  # which generators will specs made from these classes be compatible with
 
  19  compatGens = findGenSubClasses('ODEsystem') 
  20  
 
21 -class compatODEComponent(Component):
22 compatibleGens=compatGens 23 targetLangs=targetLangs # (from common.py) -- all are compatible with ODEs
24
25 -class compatODELeafComponent(LeafComponent):
26 compatibleGens=compatGens 27 targetLangs=targetLangs # (from common.py) -- all are compatible with ODEs
28 29 # ----------------------------------------------------------------------------- 30 31 # basic classes -- not intended for direct use 32 # indicates that only ODE-based Generators are compatible
33 -class channel(compatODELeafComponent):
34 pass
35
36 -class compartment(compatODEComponent):
37 pass
38 39 compartment.compatibleSubcomponents=(channel,) 40 channel.compatibleContainers=(compartment,) 41
42 -class channel_on(channel):
43 pass
44
45 -class channel_off(channel):
46 pass
47 48 49 # use these for multi-compartment cell model:
50 -class soma(compartment):
51 pass
52
53 -class dendr_compartment(compartment):
54 pass
55
56 -class neurite_compartment(compartment):
57 pass
58
59 -class neuron(compatODEComponent):
60 pass
61 62 63 # coupling
64 -class synapse(compatODELeafComponent):
65 pass
66
67 -class exc_synapse(synapse):
68 pass
69
70 -class inh_synapse(synapse):
71 pass
72 73 # point-neuron network
74 -class pnnetwork(compatODEComponent):
75 pass
76 77 # regular network
78 -class network(compatODEComponent):
79 pass
80 81 soma.compatibleContainers=(neuron, pnnetwork) 82 neuron.compatibleSubcomponents=(soma, dendr_compartment, neurite_compartment, 83 synapse) 84 neuron.compatibleContainers=(network,) 85 synapse.compatibleContainers=(neuron, pnnetwork) 86 pnnetwork.compatibleSubcomponents=(soma, synapse) 87 network.compatibleSubcomponents=(neuron,) 88 89 # helpful defaults exported for end-users 90 global voltage, V 91 voltage = 'V' 92 V = Var(voltage) 93 94 95 # ---------------------------------------------------------------------------- 96 ## Object factory functions, and related utilities 97
98 -def disconnectSynapse(syn, mspec):
99 """Disconnect synapse object syn from ModelSpec object mspec in 100 which it has been declared, and remove its target cell's synaptic 101 channels.""" 102 targsfound = [] 103 for targname in syn.connxnTargets: 104 targsfound.extend(searchModelSpec(mspec, targname)) 105 if len(targsfound) == len(syn.connxnTargets): 106 if remain(targsfound, syn.connxnTargets) == []: 107 # all targets were found 108 mspec.remove(targsfound) 109 else: 110 raise ValueError("synapse targets were not found in model spec" 111 "argument") 112 else: 113 raise ValueError("synapse targets were not found in model spec" 114 "argument") 115 syn.delConnxnTarget(targsfound) 116 mspec.remove(syn)
117 118
119 -def connectWithSynapse(synname, syntypestr, source_cell, dest_cell, 120 dest_compartment_name="", 121 threshfun=None, alpha=None, beta=None, 122 threshfun_d=None, alpha_d=None, beta_d=None, 123 adapt_typestr=None, vrev=None, g=None, 124 noauxs=True, subclass=channel):
125 """Make a chemical or electrical synapse between two neurons. For 126 gap junctions this function returns None, otherwise will return the 127 synapse gating variable object. 128 129 Valid source_cell and dest_cell is a neuron or, for point neurons only, 130 a soma is allowed. 131 The optional argument dest_compartment_name is a declared name in 132 dest_cell, or by default will be the soma. 133 134 The standard assumed rate equation for a chemical synapse conductance is 135 136 s' = (1-alpha) * s - beta * s 137 138 Use the _d versions of threshfun, alpha, and beta, for an adapting chemical 139 synapse type's adapting variable. 140 """ 141 if isinstance(source_cell, compartment): 142 addedNeuron1Name = "" 143 soma1name = source_cell.name 144 elif isinstance(source_cell, neuron): 145 addedNeuron1Name = source_cell.name + '.' 146 soma1name = source_cell._componentTypeMap['soma'][0].name 147 else: 148 raise TypeError("Invalid cell type to connect with synapse") 149 if isinstance(dest_cell, compartment): 150 addedNeuron2Name = "" 151 if dest_compartment_name == "": 152 comp2name = dest_cell.name 153 else: 154 raise ValueError("Cannot specify destination compartment name " 155 "when dest_cell has type soma") 156 elif isinstance(dest_cell, neuron): 157 addedNeuron2Name = dest_cell.name + '.' 158 if dest_compartment_name == "": 159 comp2name = dest_cell._componentTypeMap['soma'][0].name 160 else: 161 comp_objlist = dest_cell._componentTypeMap['compartment'] 162 comp_names = [o.name for o in comp_objlist] 163 # check legitimacy of dest_compartment_name in dest_cell 164 if dest_cell._componentTypeMap['soma'][0].name == dest_compartment_name: 165 comp2name = dest_compartment_name 166 elif dest_compartment_name in comp_names: 167 comp2name = dest_compartment_name 168 else: 169 raise ValueError("Destination compartment name " \ 170 + dest_compartment_name \ 171 + "not found in destination cell") 172 else: 173 raise TypeError("Invalid cell type to connect with synapse") 174 175 if addedNeuron1Name != "": 176 # we were given a neuron for source_cell, so should use its name 177 # to help distinguish multiple occurrences of svar name at time 178 # of registration in a network component. 179 svarname1 = source_cell.name+"_"+soma1name 180 else: 181 svarname1 = soma1name 182 if addedNeuron2Name != "": 183 # we were given a neuron for dest_cell, so should use its name 184 # to help distinguish multiple occurrences of svar name at time 185 # of registration in a network component. 186 svarname2 = dest_cell.name+"_"+comp2name 187 else: 188 svarname2 = comp2name 189 190 if syntypestr == 'gap': 191 # make the first 1/2 of the channel going from source -> dest 192 channel_s12 = makeSynapseChannel(synname+'_1', None, 193 comp2name+'.'+voltage, 194 syntypestr, vrev=soma1name+'.'+voltage, 195 g=g, noauxs=noauxs, 196 subclass=subclass, 197 gamma1={voltage: (synname+'_1',soma1name+'.'+voltage), 198 synname+'_1': (voltage,)}) 199 channel_name = nameResolver(channel_s12, dest_cell) 200 # update with globalized name 201 targetChannel = addedNeuron2Name+comp2name+'.'+channel_name 202 channel_s12.name = channel_name 203 # make the other 1/2 channel going in the other direction 204 channel_s21 = makeSynapseChannel(synname+'_2', None, 205 soma1name+'.'+voltage, 206 syntypestr, vrev=comp2name+'.'+voltage, 207 g=g, noauxs=noauxs, 208 subclass=subclass, 209 gamma1={voltage: (synname+'_2',comp2name+'.'+voltage), 210 synname+'_2': (voltage,)}) 211 channel_name = nameResolver(channel_s21, source_cell) 212 # update with globalized name 213 targetChannel = addedNeuron1Name+soma1name+'.'+channel_name 214 channel_s21.name = channel_name 215 else: 216 # chemical synapses 217 svar = Var('s_'+svarname1+'_'+svarname2) 218 svarname = nameResolver(svar, dest_cell) 219 channel_s12 = makeSynapseChannel(synname, 220 addedNeuron1Name+synname+'.'+svarname, 221 comp2name+'.'+voltage, 222 syntypestr, vrev=vrev, g=g, noauxs=noauxs, 223 subclass=subclass, 224 gamma1=args(voltage=(synname,), 225 synname=(voltage,))) 226 channel_name = nameResolver(channel_s12, dest_cell) 227 # update with globalized name 228 targetChannel = addedNeuron2Name+comp2name+'.'+channel_name 229 channel_s12.name = channel_name 230 # check that this source synaptic gating variable hasn't been made before, 231 # otherwise re-use that variable? 232 if threshfun_d == alpha_d == beta_d == None: 233 syn12 = makeSynapse(synname, svarname, 234 addedNeuron1Name+soma1name+'.'+voltage, 235 syntypestr, threshfun=threshfun, alpha=alpha, 236 beta=beta, targetchannel=targetChannel, noauxs=noauxs) 237 else: 238 syn12 = makeAdaptingSynapse(synname, svarname, 'd', 239 addedNeuron1Name+soma1name+'.'+voltage, 240 syntypestr, adapt_typestr, 241 threshfun=threshfun, alpha=alpha, beta=beta, 242 threshfun_d=threshfun_d, alpha_d=alpha_d, beta_d=beta_d, 243 targetchannel=targetChannel, noauxs=noauxs) 244 if addedNeuron2Name == "": 245 # we were given a soma for dest_cell 246 dest_cell.add(channel_s12) 247 else: 248 comp = dest_cell.components[comp2name] 249 dest_cell.remove(comp2name) 250 comp.add(channel_s12) 251 dest_cell.add(comp) 252 if syntypestr == 'gap': 253 source_cell.addConnxnTarget(synname+'_1') 254 source_cell.add(channel_s21) 255 dest_cell.addConnxnTarget(synname+'_2') 256 return None 257 else: 258 source_cell.addConnxnTarget(synname) 259 if addedNeuron1Name != "": 260 # we were given a neuron class for source_cell, so we have to 261 # add the synaptic variable component to that cell 262 source_cell.add(syn12) 263 return syn12
264 265 266 # ---------------------------------------------------------------------------- 267 268
269 -def makeSynapseChannel(name, gatevarname=None, voltage=voltage, typestr=None, 270 vrev=None, g=None, parlist=None, noauxs=True, subclass=channel, 271 gamma1=None, gamma2=None):
272 """Make a chemical or electrical (gap junction) synapse channel in a soma. 273 To select these, the typestr argument must be one of 'exc', 'inh', or 'gap'. 274 For a user-defined chemical synapse use a different string name for typestr. 275 For gap junctions, use the pre-synaptic membrane potential's soma ModelSpec 276 object or its full hierarchical name string for the argument vrev. 277 278 Returns a channel-type object by default, or some user-defined subclass if desired. 279 """ 280 channel_s = subclass(name) 281 if vrev is None: 282 if typestr == 'gap': 283 raise ValueError("Must provide vrev for gap junctions") 284 elif typestr == 'exc': 285 vrevpar = Par('0', 'vrev') 286 elif typestr == 'inh': 287 vrevpar = Par('-75', 'vrev') 288 elif isinstance(typestr, str): 289 vrevpar = Par('vrev') 290 else: 291 raise TypeError("Invalid type for synapse") 292 elif typestr == 'gap': 293 if isinstance(vrev, str): 294 vrevpar = Var(vrev) 295 else: 296 vrevpar = vrev 297 elif isinstance(vrev, str): 298 vrevpar = Par(vrev, 'vrev') 299 elif isinstance(vrev, float) or isinstance(vrev, int): 300 vrevpar = Par(repr(vrev), 'vrev') 301 else: 302 # allow for ModelSpec function of other variables 303 vrevpar = vrev 304 v_pot = Var(voltage) # this is to provide the local name only 305 if g is None: 306 gpar = Par('g') 307 elif isinstance(g, str): 308 gpar = Par(g, 'g') 309 else: 310 gpar = Par(repr(g), 'g') 311 if typestr == 'gap': 312 condQ = gpar 313 if gatevarname is not None: 314 raise ValueError("gatevarname must be None for gap junctions") 315 else: 316 gatevar = Var(gatevarname) 317 condQ = gpar*gatevar 318 319 if noauxs: 320 I = Var(condQ*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 321 if typestr == 'gap': 322 # don't include vrevpar, which is really an already existing 323 # V variable in a different cell 324 channel_s.add([gpar,I]) 325 else: 326 channel_s.add([gpar,vrevpar,I]) 327 else: 328 cond = Var(condQ, 'cond', [0,Inf], specType='ExpFuncSpec') 329 shunt = Var(cond*vrevpar, 'shunt', specType='ExpFuncSpec') 330 I = Var(cond*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 331 if typestr == 'gap': 332 arg_list = [] 333 else: 334 arg_list = [gatevarname] 335 if not isinstance(vrevpar, Par): 336 if parlist is not None: 337 avoid_pars = [str(q) for q in parlist if isinstance(q, Par)] 338 else: 339 avoid_pars = [] 340 arg_list.extend(remain(vrevpar.freeSymbols, avoid_pars)) 341 cond_fn = Fun(condQ(), arg_list, 'dssrt_fn_cond') 342 shunt_fn = Fun(condQ()*vrevpar, arg_list, 'dssrt_fn_shunt') 343 I_fn = Fun(condQ()*(v_pot-vrevpar), arg_list+[voltage], 344 'dssrt_fn_I') 345 channel_s.add([gpar,I,cond,shunt,cond_fn,shunt_fn,I_fn]) 346 if typestr != 'gap': 347 # don't include vrevpar, which is really an already existing 348 # V variable in a different cell 349 channel_s.add([vrevpar]) 350 if gamma1 is None: 351 gamma1 = {} 352 if gamma2 is None: 353 gamma2 = {} 354 channel_s.gamma1 = gamma1 355 channel_s.gamma2 = gamma2 356 return channel_s
357 358
359 -def makeExtInputCurrentChannel(name, noauxs=True, subclass=channel, 360 gamma1=None, gamma2=None):
361 """External input signal used directly as a current. Supply the Input having 362 coordinate name 'ext_input'. 363 364 Returns a channel-type object by default, or some user-defined subclass if desired. 365 """ 366 channel_Iext = subclass(name) 367 # Input doesn't need a name for the internal signal (confused ModelSpec as an 368 # unresolved name!) 369 inp = Input('', 'ext_input', [0,Inf], specType='ExpFuncSpec') 370 I_ext = Var(-inp, 'I', specType='ExpFuncSpec') 371 if noauxs: 372 channel_Iext.add([inp, I_ext]) 373 else: 374 arg_list = ['ext_input'] 375 I_cond = Var('0', 'cond', [0,Inf], specType='ExpFuncSpec') # no conductance for this plain current 376 I_shunt = Var(inp, 'shunt', specType='ExpFuncSpec') # technically not a shunt but acts like a V-independent shunt! 377 cond_fn = Fun(inp(), arg_list, 'dssrt_fn_cond') 378 shunt_fn = Fun(inp()*vrevpar, arg_list, 'dssrt_fn_shunt') 379 I_fn = Fun(inp()*(v_pot-vrevpar), arg_list+[voltage], 380 'dssrt_fn_I') 381 channel_Iext.add([inp, I_ext, I_cond, I_shunt, 382 cond_fn, shunt_fn, I_fn]) 383 if gamma1 is not None: 384 raise ValueError("gamma1 must be None") 385 channel_Iext.gamma1 = {} 386 if gamma2 is None: 387 gamma2 = {} 388 channel_Iext.gamma2 = gamma2 389 return channel_Iext
390 391
392 -def makeExtInputConductanceChannel(name, voltage=voltage, 393 g=None, vrev=None, parlist=None, 394 noauxs=True, subclass=channel, 395 gamma1=None, gamma2=None):
396 """External input signal used as a conductance. Supply the Input having 397 coordinate name 'ext_input'. 398 399 Returns a channel-type object by default, or some user-defined subclass if desired. 400 """ 401 if g is None: 402 gpar = Par('g') 403 elif isinstance(g, str): 404 gpar = Par(g, 'g') 405 else: 406 gpar = Par(repr(g), 'g') 407 if vrev is None: 408 vrevpar = Par('vrev') 409 elif isinstance(vrev, str): 410 try: 411 val = float(vrev) 412 except ValueError: 413 # allow for ModelSpec reference to or function of other variables 414 vrevpar = QuantSpec('vrev', vrev) 415 else: 416 vrevpar = Par(vrev, 'vrev') 417 elif isinstance(vrev, float) or isinstance(vrev, int): 418 vrevpar = Par(repr(vrev), 'vrev') 419 else: 420 # allow for ModelSpec function of other variables 421 vrevpar = vrev 422 if hasattr(vrevpar, 'name'): 423 declare_list = [gpar, vrevpar] 424 else: 425 declare_list = [gpar] 426 channel_Iext = subclass(name) 427 v_pot = Var(voltage) # this is to provide the local name only 428 # Input doesn't need a name for the internal signal (confused ModelSpec as an 429 # unresolved name!) 430 inp = Input('', 'ext_input', [0,Inf], specType='ExpFuncSpec') 431 I_ext = Var(gpar*inp*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 432 channel_Iext.add(declare_list) 433 channel_Iext.add([inp, I_ext]) 434 if not noauxs: 435 I_cond = Var(gpar*inp, 'cond', [0,Inf], specType='ExpFuncSpec') # no conductance for this plain current 436 I_shunt = Var(I_cond*vrevpar, 'shunt', specType='ExpFuncSpec') # technically not a shunt but acts like a V-independent shunt! 437 arg_list = ['ext_input'] 438 if not isinstance(vrevpar, Par): 439 if parlist is not None: 440 avoid_pars = [str(q) for q in parlist if isinstance(q, Par)] 441 else: 442 avoid_pars = [] 443 arg_list.extend(remain(vrevpar.freeSymbols, avoid_pars)) 444 cond_fn = Fun(gpar*inp(), arg_list, 'dssrt_fn_cond') 445 shunt_fn = Fun(gpar*inp()*vrevpar, arg_list, 'dssrt_fn_shunt') 446 I_fn = Fun(gpar*inp()*(v_pot-vrevpar), arg_list+[voltage], 447 'dssrt_fn_I') 448 channel_Iext.add([I_cond, I_shunt, cond_fn, shunt_fn, I_fn]) 449 if gamma1 is not None: 450 raise ValueError("gamma1 must be None") 451 channel_Iext.gamma1 = {} 452 if gamma2 is None: 453 gamma2 = {} 454 channel_Iext.gamma2 = gamma2 455 return channel_Iext
456 457
458 -def makeFunctionConductanceChannel(name, parameter_name, 459 func_def_str, 460 voltage=voltage, 461 g=None, vrev=None, parlist=None, 462 noauxs=True, subclass=channel, 463 gamma1=None, gamma2=None):
464 """Explicit function waveform used as a conductance, e.g. for an alpha-function 465 post-synaptic response. Creates a time-activated event to trigger the waveform 466 based on a named parameter. The function will be a function of (local/relative) 467 time t only. 468 469 Returns a channel-type object by default, or some user-defined subclass if desired. 470 """ 471 if g is None: 472 gpar = Par('g') 473 elif isinstance(g, str): 474 gpar = Par(g, 'g') 475 else: 476 gpar = Par(repr(g), 'g') 477 if vrev is None: 478 vrevpar = Par('vrev') 479 elif isinstance(vrev, str): 480 try: 481 val = float(vrev) 482 except ValueError: 483 # allow for ModelSpec reference to or function of other variables 484 vrevpar = QuantSpec('vrev', vrev) 485 else: 486 vrevpar = Par(vrev, 'vrev') 487 elif isinstance(vrev, float) or isinstance(vrev, int): 488 vrevpar = Par(repr(vrev), 'vrev') 489 else: 490 # allow for ModelSpec function of other variables 491 vrevpar = vrev 492 if hasattr(vrevpar, 'name'): 493 declare_list = [gpar, vrevpar] 494 else: 495 declare_list = [gpar] 496 func_def = QuantSpec('f', func_def_str) 497 # define parameters based on remaining free names in definition 498 declare_list.extend([Par(0, pname) for pname in remain(func_def.freeSymbols, 499 ['t'])]) 500 fun_par = Par(0, parameter_name) 501 inp = Fun('if(t > %s, %s, 0)' % (parameter_name, func_def_str), 502 ['t'], 'func', [-Inf,Inf], specType='ExpFuncSpec') 503 v_pot = Var(voltage) # this is to provide the local name only 504 I_ext = Var(gpar*'func(t)'*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 505 channel_Ifunc = subclass(name) 506 declare_list.extend([fun_par, inp, I_ext]) 507 channel_Ifunc.add(declare_list) 508 if not noauxs: 509 I_cond = Var(gpar*'func(t)', 'cond', [0,Inf], specType='ExpFuncSpec') # no conductance for this plain current 510 I_shunt = Var(I_cond*vrevpar, 'shunt', specType='ExpFuncSpec') # technically not a shunt but acts like a V-independent shunt! 511 arg_list = ['t'] 512 if not isinstance(vrevpar, Par): 513 if parlist is not None: 514 avoid_pars = [str(q) for q in parlist if isinstance(q, Par)] 515 else: 516 avoid_pars = [] 517 arg_list.extend(remain(vrevpar.freeSymbols, avoid_pars)) 518 cond_fn = Fun(gpar*'func(t)', arg_list, 'dssrt_fn_cond') 519 shunt_fn = Fun(gpar*'func(t)'*vrevpar, arg_list, 'dssrt_fn_shunt') 520 I_fn = Fun(gpar*'func(t)'*(v_pot-vrevpar), arg_list+[voltage], 521 'dssrt_fn_I') 522 channel_Ifunc.add([I_cond, I_shunt, cond_fn, shunt_fn, I_fn]) 523 if gamma1 is not None: 524 raise ValueError("gamma1 must be None") 525 channel_Ifunc.gamma1 = {} 526 if gamma2 is None: 527 gamma2 = {} 528 channel_Ifunc.gamma2 = gamma2 529 return channel_Ifunc
530 531
532 -def makeBiasChannel(name, I=None, noauxs=True, subclass=channel, 533 gamma1=None, gamma2=None):
534 """Constant bias / applied current "channel". 535 536 Returns a channel-type object by default, or some user-defined subclass if desired. 537 """ 538 channel_Ib = subclass(name) 539 if I is None: 540 Ibias = Par('Ibias') 541 elif isinstance(I, str): 542 Ibias = Par(I, 'Ibias') 543 else: 544 Ibias = Par(repr(I), 'Ibias') 545 Ib = Var(-Ibias, 'I', specType='ExpFuncSpec') 546 if noauxs: 547 channel_Ib.add([Ibias,Ib]) 548 else: 549 Ibias_cond = Var('0', 'cond', [0,Inf], specType='ExpFuncSpec') # no conductance for this plain current 550 Ibias_shunt = Var(Ibias, 'shunt', specType='ExpFuncSpec') # technically not a shunt but acts like a V-independent shunt! 551 cond_fn = Fun('0', [], 'dssrt_fn_cond') 552 shunt_fn = Fun(Ibias, [], 'dssrt_fn_shunt') 553 I_fn = Fun(Ibias, [], 'dssrt_fn_I') 554 channel_Ib.add([Ibias,Ibias_cond,Ibias_shunt,Ib, cond_fn, shunt_fn, I_fn]) 555 if gamma1 is not None: 556 raise ValueError("gamma1 must be None") 557 channel_Ib.gamma1 = {} 558 if gamma2 is None: 559 gamma2 = {} 560 channel_Ib.gamma2 = gamma2 561 return channel_Ib
562 563
564 -def makeChannel_halfact(name,voltage=voltage,s=None,isinstant=False,sinf=None, 565 taus=None,spow=1,s2=None,isinstant2=False,sinf2=None,taus2=None, 566 spow2=1,vrev=None,g=None, 567 parlist=None, noauxs=True, subclass=channel, 568 gamma1=None, gamma2=None, nonlocal_variables=None):
569 """Make an ionic membrane channel using the steady state and rate function formalism. 570 571 i.e., that the gating variable s has a differential equation in the form: 572 s' = (sinf(v) - v)/taus(v) 573 The channel may have up to two gating variables, each of which is given by an ODE. 574 575 If either gating variable has its corresponding 'isinstant' argument set to True, then 576 that variable is set to be instananeous (algebraic, not an ODE), i.e. s = sinf(v). The 577 taus function will then be ignored. 578 579 The resulting channel current will be of the form 580 name.I = g * s^spow * s2^spow2 * (voltage - vrev) 581 582 Provide any additional Par or Fun objects necessary for the complete definition of the 583 channel kinetics in the optional parlist argument. 584 585 nonlocal_variables list (optional) provides string names of any dynamic variables 586 referenced in the declared specifications that are not local to this channel. 587 588 Returns a channel-type object by default, or some user-defined subclass if desired. 589 """ 590 if g is None: 591 gpar = Par('g') 592 elif isinstance(g, str): 593 gpar = Par(g, 'g') 594 else: 595 gpar = Par(repr(g), 'g') 596 if vrev is None: 597 vrevpar = Par('vrev') 598 elif isinstance(vrev, str): 599 try: 600 val = float(vrev) 601 except ValueError: 602 # allow for ModelSpec reference to or function of other variables 603 vrevpar = QuantSpec('vrev', vrev) 604 else: 605 vrevpar = Par(vrev, 'vrev') 606 elif isinstance(vrev, float) or isinstance(vrev, int): 607 vrevpar = Par(repr(vrev), 'vrev') 608 else: 609 # allow for ModelSpec function of other variables 610 vrevpar = vrev 611 if hasattr(vrevpar, 'name'): 612 declare_list = [gpar, vrevpar] 613 else: 614 declare_list = [gpar] 615 thischannel = subclass(name) 616 v_pot = Var(voltage) # this is to provide the local name only 617 if nonlocal_variables is None: 618 nonlocal_variables = [] 619 # deal with gating variable, if present 620 if s is None: 621 condQ = gpar 622 assert taus==sinf==None 623 assert spow==1, "power should be unused when gating variable is None" 624 # override 625 spow=0 626 else: 627 if taus==sinf==None: 628 assert isinstant == False 629 gatevar = Par(s) 630 else: 631 if isinstance(sinf, (Quantity, QuantSpec)): 632 sinf_term = Var(sinf(), name=s+'inf') 633 elif isinstance(sinf, str): 634 sinf_term = Var(sinf, name=s+'inf') 635 else: 636 raise TypeError("sinf argument must be a string, Quantity or a QuantSpec") 637 if isinstant: 638 if noauxs: 639 gatevar = Var(sinf_term(), name=s, domain=[0,1], 640 specType='ExpFuncSpec') 641 else: 642 taus_term = Var('0', name='tau'+s) 643 taus_fn = Fun( '0', [v_pot]+nonlocal_variables, 644 'dssrt_fn_tau'+s ) 645 sinf_fn = Fun( sinf_term(), [v_pot]+nonlocal_variables, 646 'dssrt_fn_'+s+'inf' ) 647 gatevar = Var(sinf_term(), name=s, domain=[0,1], 648 specType='ExpFuncSpec') 649 else: 650 if isinstance(taus, (Quantity, QuantSpec)): 651 taus_term = Var(taus(), name='tau'+s) 652 elif isinstance(taus, str): 653 taus_term = Var(taus, name='tau'+s) 654 else: 655 raise TypeError("taus argument must be a string, Quantity or a QuantSpec") 656 # temporary declaration of symbol s (the argument string) as a Var 657 s_ = Var(s) 658 if noauxs: 659 gatevar = Var( (sinf_term() - s_) / taus_term(), name=s, 660 domain=[0,1], specType='RHSfuncSpec') 661 else: 662 taus_fn = Fun( taus_term(), [v_pot]+nonlocal_variables, 663 'dssrt_fn_tau'+s ) 664 sinf_fn = Fun( sinf_term(), [v_pot]+nonlocal_variables, 665 'dssrt_fn_'+s+'inf' ) 666 gatevar = Var( (sinf_term() - s_) / taus_term(), name=s, 667 domain=[0,1], specType='RHSfuncSpec') 668 if not noauxs: 669 declare_list.extend([taus_term,sinf_term,taus_fn,sinf_fn]) 670 declare_list.append(gatevar) 671 condQ = gpar*makePowerSpec(gatevar,spow) 672 673 # deal with second gating variable, if present 674 if s2 is not None: 675 assert s is not None, "Must use first gating variable name first!" 676 assert spow2 != 0, "Do not use second gating variable with spow2 = 0!" 677 s2_ = Var(s2) 678 if isinstance(sinf2, QuantSpec): 679 sinf2_term = Var(sinf2(), name=s2+'inf', specType='ExpFuncSpec') 680 elif isinstance(sinf2, str): 681 sinf2_term = Var(sinf2, name=s2+'inf', specType='ExpFuncSpec') 682 else: 683 raise TypeError("sinf2 argument must be a string or a QuantSpec") 684 if isinstant2: 685 if noauxs: 686 gatevar2 = Var(sinf2_term(), name=s2, domain=[0,1], 687 specType='ExpFuncSpec') 688 else: 689 taus2_term = Var('0', name='tau'+s2) 690 taus2_fn = Fun( '0', [v_pot]+nonlocal_variables, 691 'dssrt_fn_tau'+s2 ) 692 sinf2_fn = Fun( sinf2_term(), [v_pot]+nonlocal_variables, 693 'dssrt_fn_'+s2+'inf' ) 694 gatevar2 = Var(sinf2_term(), name=s2, domain=[0,1], 695 specType='ExpFuncSpec') 696 else: 697 if isinstance(taus2, QuantSpec): 698 taus2_term = Var(taus2(), name='tau'+s2, specType='ExpFuncSpec') 699 elif isinstance(taus2, str): 700 taus2_term = Var(taus2, name='tau'+s2, specType='ExpFuncSpec') 701 else: 702 raise TypeError("taus2 argument must be a string or a QuantSpec") 703 if noauxs: 704 gatevar2 = Var( (sinf2_term() - s2_) / taus2_term(), name=s2, 705 domain=[0,1], specType='RHSfuncSpec') 706 else: 707 taus2_fn = Fun( taus2_term(), [v_pot]+nonlocal_variables, 708 'dssrt_fn_tau'+s2 ) 709 sinf2_fn = Fun( sinf2_term(), [v_pot]+nonlocal_variables, 710 'dssrt_fn_'+s2+'inf' ) 711 gatevar2 = Var( (sinf2_term() - s2_) / taus2_term(), name=s2, 712 domain=[0,1], specType='RHSfuncSpec') 713 declare_list.append(gatevar2) 714 if not noauxs: 715 declare_list.extend([taus2_term, sinf2_term, taus2_fn, sinf2_fn]) 716 condQ = condQ * makePowerSpec(gatevar2, spow2) 717 718 if noauxs: 719 I = Var(condQ*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 720 declare_list.append(I) 721 else: 722 cond = Var(condQ, 'cond', [0,Inf], specType='ExpFuncSpec') 723 shunt = Var(cond*vrevpar, 'shunt', specType='ExpFuncSpec') 724 I = Var(cond*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 725 if s is None: 726 arg_list = [] 727 else: 728 arg_list = [s] 729 if s2 is not None: 730 arg_list.append(s2) 731 # vrev may be a Par, a singleton reference to another parameter already defined, 732 # which counts the same; or else a compound reference (function of parameters) 733 if not isinstance(vrevpar, Par) and vrevpar.freeSymbols != [str(vrevpar)]: 734 if parlist is not None: 735 avoid_pars = [str(q) for q in parlist if isinstance(q, Par)] 736 else: 737 avoid_pars = [] 738 arg_list.extend(remain(vrevpar.freeSymbols, avoid_pars)) 739 condQ_var_free = intersect(condQ.freeSymbols, nonlocal_variables) 740 cond_var_free = intersect(cond.freeSymbols, nonlocal_variables) 741 cond_fn = Fun(condQ(), arg_list+condQ_var_free, 'dssrt_fn_cond') 742 shunt_fn = Fun(cond()*vrevpar, arg_list+cond_var_free, 'dssrt_fn_shunt') 743 I_fn = Fun(cond()*(v_pot-vrevpar), arg_list+[voltage]+cond_var_free, 744 'dssrt_fn_I') 745 declare_list.extend([I,cond,shunt, cond_fn, shunt_fn, I_fn]) 746 if gamma1 is None: 747 gamma1 = {} 748 if gamma2 is None: 749 gamma2 = {} 750 thischannel.gamma1 = gamma1 751 thischannel.gamma2 = gamma2 752 if parlist is not None: 753 declare_list.extend(parlist) 754 thischannel.add(declare_list) 755 return thischannel
756 757
758 -def makeChannel_rates(name,voltage=voltage, 759 s=None,isinstant=False,arate=None,brate=None,spow=1, 760 s2=None,isinstant2=False,arate2=None,brate2=None,spow2=1, 761 vrev=None,g=None,parlist=None,noauxs=True,subclass=channel, 762 gamma1=None, gamma2=None, nonlocal_variables=None):
763 """Make an ionic membrane channel using the forward and backward rate formalism. 764 765 i.e., that the gating variable s has a differential equation in the form: 766 s' = sa(v) * (1-s) - sb(v) * s 767 The channel may have up to two gating variables, each of which is given by an ODE. 768 769 If either gating variable has its corresponding 'isinstant' argument set to True, then 770 that variable is set to be instananeous (algebraic, not an ODE), 771 i.e. s = sinf(v) = a(v) / (a(v) + b(v)) 772 773 The resulting channel current will be of the form 774 name.I = g * s^spow * s2^spow2 * (voltage - vrev) 775 776 Provide any additional Par or Fun objects necessary for the complete definition of the 777 channel kinetics in the optional parlist argument. 778 779 nonlocal_variables list (optional) provides string names of any dynamic variables 780 referenced in the declared specifications that are not local to this channel. 781 782 Returns a channel-type object by default, or some user-defined subclass if desired. 783 """ 784 if g is None: 785 gpar = Par('g') 786 elif isinstance(g, str): 787 gpar = Par(g, 'g') 788 else: 789 gpar = Par(repr(g), 'g') 790 if vrev is None: 791 vrevpar = Par('vrev') 792 elif isinstance(vrev, str): 793 try: 794 val = float(vrev) 795 except ValueError: 796 # allow for ModelSpec reference to or function of other variables 797 vrevpar = QuantSpec('vrev', vrev) 798 else: 799 vrevpar = Par(vrev, 'vrev') 800 elif isinstance(vrev, float) or isinstance(vrev, int): 801 vrevpar = Par(repr(vrev), 'vrev') 802 else: 803 # allow for ModelSpec function of other variables 804 vrevpar = vrev 805 if hasattr(vrevpar, 'name'): 806 declare_list = [gpar, vrevpar] 807 else: 808 declare_list = [gpar] 809 thischannel = subclass(name) 810 v_pot = Var(voltage) # this is to provide the local name only 811 if nonlocal_variables is None: 812 nonlocal_variables = [] 813 # deal with gating variable, if present 814 if s is None: 815 condQ = gpar 816 assert arate==brate==None 817 assert spow==1, "power should be unused when gating variable is None" 818 # override 819 spow=0 820 else: 821 if arate==brate==None: 822 assert isinstant == False 823 gatevar = Par(s) 824 else: 825 if isinstance(arate, (Quantity, QuantSpec)): 826 aterm = Var(arate(), name='a'+s) 827 elif isinstance(arate, str): 828 aterm = Var(arate, name='a'+s) 829 else: 830 raise TypeError("arate argument must be a string or QuantSpec") 831 if isinstance(brate, (Quantity, QuantSpec)): 832 bterm = Var(brate(), name='b'+s) 833 elif isinstance(brate, str): 834 bterm = Var(brate, name='b'+s) 835 else: 836 raise TypeError("brate argument must be a string or QuantSpec") 837 # temporary declaration of symbol s (the argument string) as a Var 838 s_ = Var(s) 839 if isinstant: 840 if noauxs: 841 gatevar = Var( aterm()/(aterm()+bterm()), name=s, 842 domain=[0,1], specType='ExpFuncSpec') 843 else: 844 taus_term = Var('0', name='tau'+s) 845 sinf_term = Var( aterm/(aterm+bterm), name=s+'inf') 846 taus_fn = Fun( '0', [voltage]+nonlocal_variables, 847 'dssrt_fn_tau'+s ) 848 sinf_fn = Fun( aterm()/(aterm()+bterm()), 849 [voltage]+nonlocal_variables, 850 'dssrt_fn_'+s+'inf' ) 851 gatevar = Var( aterm()/(aterm()+bterm()), name=s, 852 domain=[0,1], specType='ExpFuncSpec') 853 else: 854 if noauxs: 855 gatevar = Var( aterm()*(1-s_)-bterm()*s_, name=s, 856 domain=[0,1], specType='RHSfuncSpec') 857 else: 858 taus_term = Var( 1/(aterm+bterm), name='tau'+s) 859 sinf_term = Var( aterm*taus_term, name=s+'inf') 860 taus_fn = Fun( 1/(aterm()+bterm()), [voltage]+nonlocal_variables, 861 'dssrt_fn_tau'+s ) 862 sinf_fn = Fun( aterm()/(aterm()+bterm()), 863 [voltage]+nonlocal_variables, 864 'dssrt_fn_'+s+'inf' ) 865 gatevar = Var( aterm()*(1-s_)-bterm()*s_, name=s, 866 domain=[0,1], specType='RHSfuncSpec') 867 if not noauxs: 868 declare_list.extend([taus_term,sinf_term,aterm,bterm, 869 taus_fn, sinf_fn]) 870 declare_list.append(gatevar) 871 condQ = gpar*makePowerSpec(gatevar,spow) 872 873 # deal with second gating variable, if present 874 if s2 is not None: 875 assert s is not None, "Must use first gating variable name first!" 876 assert spow2 != 0, "Do not use second gating variable with spow2 = 0!" 877 if isinstance(arate2, (Quantity, QuantSpec)): 878 aterm2 = Var(arate2(), name='a'+s2) 879 elif isinstance(arate2, str): 880 aterm2 = Var(arate2, name='a'+s2) 881 else: 882 raise TypeError("arate2 argument must be a string, Quantity or QuantSpec") 883 if isinstance(brate2, (Quantity, QuantSpec)): 884 bterm2 = Var(brate2(), name='b'+s2) 885 elif isinstance(brate2, str): 886 bterm2 = Var(brate2, name='b'+s2) 887 else: 888 raise TypeError("brate2 argument must be a string, Quantity or QuantSpec") 889 s2_ = Var(s2) 890 if isinstant2: 891 if noauxs: 892 gatevar2 = Var( aterm2()/(aterm2()+bterm2()), name=s2, 893 domain=[0,1], specType='ExpFuncSpec') 894 else: 895 taus2_term = Var( '0', name='tau'+s2) 896 sinf2_term = Var( aterm2()/(aterm2()+bterm2()), name=s2+'inf') 897 taus2_fn = Fun( '0', [voltage]+nonlocal_variables, 898 'dssrt_fn_tau'+s2 ) 899 sinf2_fn = Fun( aterm2()/(aterm2()+bterm2()), 900 [voltage]+nonlocal_variables, 901 'dssrt_fn_'+s2+'inf' ) 902 gatevar2 = Var( aterm2()/(aterm2()+bterm2()), name=s2, 903 domain=[0,1], specType='ExpFuncSpec') 904 else: 905 if noauxs: 906 gatevar2 = Var(aterm2()*(1-s2_)-bterm2()*s2_, name=s2, 907 domain=[0,1], specType='RHSfuncSpec') 908 else: 909 taus2_term = Var( 1/(aterm2()+bterm2()), name='tau'+s2) 910 sinf2_term = Var( aterm2()*taus2_term(), name=s2+'inf') 911 taus2_fn = Fun( 1/(aterm2()+bterm2()), [voltage]+nonlocal_variables, 912 'dssrt_fn_tau'+s2 ) 913 sinf2_fn = Fun( aterm2()/(aterm2()+bterm2()), 914 [voltage]+nonlocal_variables, 915 'dssrt_fn_'+s2+'inf' ) 916 gatevar2 = Var(aterm2()*(1-s2_)-bterm2()*s2_, name=s2, 917 domain=[0,1], specType='RHSfuncSpec') 918 condQ = condQ * makePowerSpec(gatevar2, spow2) 919 declare_list.append(gatevar2) 920 if not noauxs: 921 declare_list.extend([taus2_term, sinf2_term, aterm2, bterm2, 922 taus2_fn, sinf2_fn]) 923 924 if noauxs: 925 I = Var(condQ*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 926 declare_list.append(I) 927 else: 928 cond = Var(condQ, 'cond', [0,Inf], specType='ExpFuncSpec') 929 shunt = Var(cond*vrevpar, 'shunt', specType='ExpFuncSpec') 930 I = Var(cond*(v_pot-vrevpar), 'I', specType='ExpFuncSpec') 931 if s is None: 932 arg_list = [] 933 else: 934 arg_list = [s] 935 if s2 is not None: 936 arg_list.append(s2) 937 # vrev may be a Par, a singleton reference to another parameter already defined, 938 # which counts the same; or else a compound reference (function of parameters) 939 if not isinstance(vrevpar, Par) and vrevpar.freeSymbols != [str(vrevpar)]: 940 if parlist is not None: 941 avoid_pars = [str(q) for q in parlist if isinstance(q, Par)] 942 else: 943 avoid_pars = [] 944 arg_list.extend(remain(vrevpar.freeSymbols, avoid_pars)) 945 condQ_var_free = intersect(condQ.freeSymbols, nonlocal_variables) 946 cond_var_free = intersect(cond.freeSymbols, nonlocal_variables) 947 cond_fn = Fun(condQ(), arg_list+condQ_var_free, 'dssrt_fn_cond') 948 shunt_fn = Fun(cond()*vrevpar, arg_list+cond_var_free, 'dssrt_fn_shunt') 949 I_fn = Fun(cond()*(v_pot-vrevpar), arg_list+[voltage]+cond_var_free, 950 'dssrt_fn_I') 951 declare_list.extend([I,cond,shunt, cond_fn, shunt_fn, I_fn]) 952 if gamma1 is None: 953 gamma1 = {} 954 if gamma2 is None: 955 gamma2 = {} 956 thischannel.gamma1 = gamma1 957 thischannel.gamma2 = gamma2 958 if parlist is not None: 959 declare_list.extend(parlist) 960 thischannel.add(declare_list) 961 return thischannel
962 963
964 -def makeSoma(name, voltage=voltage, channelList=None, C=None, noauxs=True, 965 subclass=soma, channelclass=channel):
966 """Build a soma type of "compartment" from a list of channels and a membrane 967 capacitance C, using the local voltage name. 968 969 Specify a sub-class for the channel to target specific types to summate, 970 in case the user wants to have runtime control over which are "active" in the 971 dynamics (mainly related to the DSSRT toolbox). 972 973 Returns a soma-type object by default, or some user-defined subclass if desired. 974 """ 975 channel_cls_name = className(channelclass) 976 v = Var(QuantSpec(str(voltage), '-for(%s,I,+)/C'%channel_cls_name, 'RHSfuncSpec'), 977 domain=[-150,100]) 978 if C is None: 979 C = Par('C') 980 elif isinstance(C, str): 981 gpar = Par(C, 'C') 982 else: 983 C = Par(repr(C), 'C') 984 985 asoma = soma.__new__(subclass) 986 asoma.__init__(name) 987 declare_list = copy(channelList) 988 if noauxs: 989 declare_list.extend([v,C]) 990 else: 991 tauv = Var(QuantSpec('tauv', '1/(C*for(%s,cond,+))'%channel_cls_name)) 992 vinf = Var(QuantSpec('vinf', 'for(%s,shunt,+)'%channel_cls_name + \ 993 '/for(%s,cond,+)'%channel_cls_name)) 994 # for each cond in channel, append signature of dssrt_fn_cond 995 arg_list_c = [] 996 arg_list_s = [] 997 c_sum = [] 998 s_sum = [] 999 dssrt_inputs = {} 1000 vstr = str(voltage) 1001 for c in channelList: 1002 def process_targ_name(inp): 1003 if inp == vstr: 1004 return inp 1005 else: 1006 return c.name+'.'+inp
1007 for targ, inputs in c.gamma1.items(): 1008 proc_targ = process_targ_name(targ) 1009 if proc_targ not in dssrt_inputs: 1010 dssrt_inputs[proc_targ] = args(gamma1=[], gamma2=[]) 1011 dssrt_inputs[proc_targ].gamma1.extend([process_targ_name(i) \ 1012 for i in inputs]) 1013 for targ, inputs in c.gamma2.items(): 1014 proc_targ = process_targ_name(targ) 1015 if proc_targ not in dssrt_inputs: 1016 dssrt_inputs[proc_targ] = args(gamma1=[], gamma2=[]) 1017 dssrt_inputs[proc_targ].gamma2.extend([process_targ_name(i) \ 1018 for i in inputs]) 1019 cond = c['dssrt_fn_cond'] 1020 shunt = c['dssrt_fn_shunt'] 1021 cond_sig = cond.signature 1022 shunt_sig = shunt.signature 1023 arg_list_c.extend(cond_sig) 1024 arg_list_s.extend(shunt_sig) 1025 c_sum.append( c.name + '.' + cond.name+'('+','.join(cond_sig)+')' ) 1026 s_sum.append( c.name + '.' + shunt.name+'('+','.join(shunt_sig)+')' ) 1027 arg_list1 = makeSeqUnique(arg_list_c) 1028 arg_list1.sort() 1029 arg_list2 = makeSeqUnique(arg_list_c+arg_list_s) 1030 arg_list2.sort() 1031 tauv_fn = Fun('1/(C*(%s))'%'+'.join(c_sum), arg_list1, 1032 'dssrt_fn_tau'+vstr) 1033 vinf_fn = Fun('(%s)/(%s)'%('+'.join(s_sum), '+'.join(c_sum)), 1034 arg_list2, 'dssrt_fn_'+vstr+'inf') 1035 declare_list.extend([v,C,vinf,tauv,tauv_fn,vinf_fn]) 1036 asoma.dssrt_inputs = dssrt_inputs 1037 asoma.add(declare_list) 1038 return asoma 1039 1040 1041
1042 -def makeDendrite(name, voltage=voltage, channelList=None, C=None, noauxs=True, 1043 subclass=dendr_compartment, channelclass=channel):
1044 """Build a dendrite type of "compartment" from a list of channels and a membrane 1045 capacitance C, using the local voltage name. 1046 1047 Specify a sub-class for the channel to target specific types to summate, 1048 in case the user wants to have runtime control over which are "active" in the 1049 dynamics (mainly related to the DSSRT toolbox). 1050 1051 Returns a dendrite-type object by default, or some user-defined subclass if desired. 1052 """ 1053 channel_cls_name = className(channelclass) 1054 v = Var(QuantSpec(str(voltage), '-for(%s,I,+)/C'%channel_cls_name, 'RHSfuncSpec'), 1055 domain=[-150,100]) 1056 if C is None: 1057 C = Par('C') 1058 elif isinstance(C, str): 1059 gpar = Par(C, 'C') 1060 else: 1061 C = Par(repr(C), 'C') 1062 1063 acomp = dendr_compartment.__new__(subclass) 1064 acomp.__init__(name) 1065 declare_list = copy(channelList) 1066 if noauxs: 1067 declare_list.extend([v,C]) 1068 else: 1069 tauv = Var(QuantSpec('tauv', '1/(C*for(%s,cond,+))'%channel_cls_name)) 1070 vinf = Var(QuantSpec('vinf', 'for(%s,shunt,+)'%channel_cls_name + \ 1071 '/for(%s,cond,+)'%channel_cls_name)) 1072 declare_list.extend([v,C,vinf,tauv]) 1073 acomp.add(declare_list) 1074 return acomp
1075 1076
1077 -def makeNeurite(name, voltage=voltage, channelList=None, C=None, noauxs=True, 1078 subclass=dendr_compartment, channelclass=channel):
1079 """Build a neurite type of "compartment" from a list of channels and a membrane 1080 capacitance C, using the local voltage name. 1081 1082 Specify a sub-class for the channel to target specific types to summate, 1083 in case the user wants to have runtime control over which are "active" in the 1084 dynamics (mainly related to the DSSRT toolbox). 1085 1086 Returns a neurite-type object by default, or some user-defined subclass if desired. 1087 """ 1088 channel_cls_name = className(channelclass) 1089 v = Var(QuantSpec(str(voltage), '-for(%s,I,+)/C'%channel_cls_name, 'RHSfuncSpec'), 1090 domain=[-150,100]) 1091 if C is None: 1092 C = Par('C') 1093 elif isinstance(C, str): 1094 gpar = Par(C, 'C') 1095 else: 1096 C = Par(repr(C), 'C') 1097 1098 acomp = neurite_compartment.__new__(subclass) 1099 acomp.__init__(name) 1100 declare_list = copy(channelList) 1101 if noauxs: 1102 declare_list.extend([v,C]) 1103 else: 1104 tauv = Var(QuantSpec('tauv', '1/(C*for(%s,cond,+))'%channel_cls_name)) 1105 vinf = Var(QuantSpec('vinf', 'for(%s,shunt,+)'%channel_cls_name + \ 1106 '/for(%s,cond,+)'%channel_cls_name)) 1107 declare_list.extend([v,C,vinf,tauv]) 1108 acomp.add(declare_list) 1109 return acomp
1110 1111
1112 -def makePointNeuron(name, voltage=voltage, channelList=None, synapseList=None, 1113 C=None, noauxs=True):
1114 """Factory function for single compartment neurons ("point neurons"). 1115 1116 Returns a neuron type object from a list of ion channels 1117 and synapses (if any are defined). The ion channels will internally be built into 1118 a soma compartment type on the fly. 1119 """ 1120 asoma = makeSoma('soma', voltage, channelList, C=C, noauxs=noauxs) 1121 n = neuron(name) 1122 n.add(asoma) 1123 if synapseList is not None: 1124 n.add(synapseList) 1125 return n
1126 1127
1128 -def makePointNeuronNetwork(name, componentList):
1129 """Factory function returning a pnnetwork type object from a list of compatible 1130 components (somas and synapses). 1131 """ 1132 net = pnnetwork(name) 1133 net.add(componentList) 1134 return net
1135 1136
1137 -def makeNeuronNetwork(name, neuronList):
1138 """Factory function returning a network type object from a list of compatible 1139 components (neurons). The neurons can be single or multi-compartmental. 1140 1141 Currently untested! 1142 """ 1143 net = network(name) 1144 net.add(neuronList) 1145 return net
1146 1147
1148 -def makeSynapse(name, gatevar, precompartment, typestr, threshfun=None, 1149 alpha=None, beta=None, targetchannel=None, 1150 evalopt=True, noauxs=True):
1151 """Make a chemical synapse channel object. 1152 """ 1153 if targetchannel is None: 1154 raise TypeError("Must provide name of synaptic channel object in " 1155 "target cell's compartment") 1156 if typestr == 'exc': 1157 if alpha is None: 1158 alpha = 10. 1159 if beta is None: 1160 beta = 0.5 1161 syn = exc_synapse(name) 1162 elif typestr == 'inh': 1163 if alpha is None: 1164 alpha = 1. 1165 if beta is None: 1166 beta = 0.1 1167 syn = inh_synapse(name) 1168 elif typestr == "": 1169 syn = synapse(name) 1170 else: 1171 raise ValueError("Invalid type of synapse specified") 1172 s_ = Var(gatevar) 1173 if alpha is None: 1174 a = Par('alpha') 1175 elif isinstance(alpha, str): 1176 gpar = Par(alpha, 'alpha') 1177 else: 1178 a = Par(repr(alpha), 'alpha') 1179 if beta is None: 1180 b = Par('beta') 1181 elif isinstance(beta, str): 1182 gpar = Par(beta, 'beta') 1183 else: 1184 b = Par(repr(beta), 'beta') 1185 if threshfun is None: 1186 f = Fun('0.5+0.5*tanh(v/4.)', ['v'], 'thresh') 1187 else: 1188 assert isinstance(threshfun, tuple), \ 1189 "threshfun must be pair (vname, funbody)" 1190 if isinstance(threshfun[1], QuantSpec): 1191 funbody = threshfun[1].specStr 1192 elif isinstance(threshfun[1], str): 1193 funbody = threshfun[1] 1194 else: 1195 raise TypeError("threshold function must be a string or a " 1196 "QuantSpec") 1197 if threshfun[0] not in funbody: 1198 print "Warning: voltage name %s does not appear in function body!"%threshfun[0] 1199 f = Fun(funbody, [threshfun[0]], 'thresh') 1200 assert len(f.signature) == 1, \ 1201 'threshold function must be a function of a single argument (voltage)' 1202 if evalopt: 1203 alpha_ = a*f.eval(**{f.signature[0]: precompartment}) 1204 else: 1205 alpha_ = a*f(precompartment) 1206 syn.add(f) 1207 s = Var(alpha_*(1-s_)-b*s_, name=gatevar, 1208 specType='RHSfuncSpec', domain=[0,1]) 1209 syn.add([s,a,b]) 1210 syn.addConnxnTarget(targetchannel) 1211 if not noauxs: 1212 sinf = Var( alpha_/(alpha_+b), name=gatevar+'inf') 1213 taus = Var( 1/(alpha_+b), name='tau'+gatevar) 1214 sinf_fn = Fun( alpha_/(alpha_+b), [precompartment], 'dssrt_fn_'+gatevar+'inf') 1215 taus_fn = Fun( 1/(alpha_+b), [precompartment], 'dssrt_fn_tau'+gatevar) 1216 syn.add([sinf,taus,sinf_fn,taus_fn]) 1217 return syn
1218 1219
1220 -def makeAdaptingSynapse(name, gatevar, adaptvar, precompartment, typestr, 1221 adapt_typestr, 1222 threshfun=None, alpha=None, beta=None, 1223 threshfun_d=None, alpha_d=None, beta_d=None, 1224 targetchannel=None, evalopt=True, noauxs=True):
1225 """Make an adapting chemical synapse channel object. 1226 """ 1227 if targetchannel is None: 1228 raise TypeError("Must provide name of synaptic channel object in " 1229 "target cell's compartment") 1230 if typestr == 'exc': 1231 if alpha is None: 1232 alpha = 10. 1233 if beta is None: 1234 beta = 0.5 1235 syn = exc_synapse(name) 1236 elif typestr == 'inh': 1237 if alpha is None: 1238 alpha = 1. 1239 if beta is None: 1240 beta = 0.1 1241 syn = inh_synapse(name) 1242 elif typestr == "": 1243 syn = synapse(name) 1244 else: 1245 raise ValueError("Invalid type of synapse specified") 1246 # synaptic variable 1247 s_ = Var(gatevar) 1248 if alpha is None: 1249 a = Par('alpha') 1250 elif isinstance(alpha, str): 1251 a = Par(alpha, 'alpha') 1252 else: 1253 a = Par(repr(alpha), 'alpha') 1254 if beta is None: 1255 b = Par('beta') 1256 elif isinstance(beta, str): 1257 b = Par(beta, 'beta') 1258 else: 1259 b = Par(repr(beta), 'beta') 1260 if threshfun is None: 1261 f = Fun('0.5+0.5*tanh(v/4.)', ['v'], 'thresh') 1262 else: 1263 assert isinstance(threshfun, tuple), \ 1264 "threshfun must be pair (vname, funbody)" 1265 if isinstance(threshfun[1], QuantSpec): 1266 funbody = threshfun[1].specStr 1267 elif isinstance(threshfun[1], str): 1268 funbody = threshfun[1] 1269 else: 1270 raise TypeError("threshold function must be a string or a " 1271 "QuantSpec") 1272 assert threshfun[0] in funbody, \ 1273 "voltage name %s does not appear in function body!"%threshfun[0] 1274 f = Fun(funbody, [threshfun[0]], 'thresh') 1275 assert len(f.signature) == 1, \ 1276 'threshold function must be a function of a single argument (voltage)' 1277 if evalopt: 1278 alpha_ = a*f.eval(**{f.signature[0]: precompartment}) 1279 else: 1280 alpha_ = a*f(precompartment) 1281 syn.add(f) 1282 # adaptive variable 1283 if adapt_typestr == 'f': 1284 at = '' 1285 elif adapt_typestr == 'd': 1286 at = '-' 1287 else: 1288 raise ValueError("Invalid type for adapting synapse: use 'f' or 'd'") 1289 d_ = Var(adaptvar) 1290 if alpha_d is None: 1291 a_d = Par('alpha_d') 1292 elif isinstance(alpha_d, str): 1293 a_d = Par(alpha_d, 'alpha_d') 1294 else: 1295 a_d = Par(repr(alpha_d), 'alpha_d') 1296 if beta_d is None: 1297 b_d = Par('beta_d') 1298 elif isinstance(beta_d, str): 1299 b_d = Par(beta_d, 'beta_d') 1300 else: 1301 b_d = Par(repr(beta_d), 'beta_d') 1302 if threshfun_d is None: 1303 f_d = Fun('0.5+0.5*tanh(%sv/4.)'%atype, ['v'], 'thresh') 1304 else: 1305 assert isinstance(threshfun_d, tuple), \ 1306 "threshfun must be pair (vname, funbody)" 1307 if isinstance(threshfun_d[1], QuantSpec): 1308 funbody_d = threshfun_d[1].specStr 1309 elif isinstance(threshfun_d[1], str): 1310 funbody_d = threshfun_d[1] 1311 else: 1312 raise TypeError("threshold function must be a string or a " 1313 "QuantSpec") 1314 assert threshfun_d[0] in funbody_d, \ 1315 "voltage name %s does not appear in function body!"%threshfun[0] 1316 f_d = Fun(funbody_d, [threshfun_d[0]], 'thresh') 1317 assert len(f_d.signature) == 1, \ 1318 'threshold function must be a function of a single argument (voltage)' 1319 if evalopt: 1320 alpha_d_ = a_d*f_d.eval(**{f_d.signature[0]: precompartment}) 1321 else: 1322 alpha_d_ = a_d*f_d(precompartment) 1323 syn.add(f_d) 1324 d = Var(alpha_d_*(1-d_)-b_d*d_, name=adaptvar, 1325 specType='RHSfuncSpec', domain=[0,1]) 1326 s = Var(alpha_*(d_-s_)-b*s_, name=gatevar, 1327 specType='RHSfuncSpec', domain=[0,1]) 1328 syn.add([s,a,b,d,a_d,b_d]) 1329 syn.addConnxnTarget(targetchannel) 1330 if not noauxs: 1331 sinf = Var( alpha_*d_/(alpha_+b), name=gatevar+'inf') 1332 taus = Var( 1/(alpha_+b), name='tau'+gatevar) 1333 dinf = Var( alpha_d_/(alpha_d_+b_d), name=adaptvar+'inf') 1334 taud = Var( 1/(alpha_d_+b_d), name='tau'+adaptvar) 1335 sinf_fn = Fun( ) 1336 taus_fn = Fun( ) 1337 syn.add([sinf,taus,dinf,taud,sinf_fn,taus_fn]) 1338 return syn
1339 1340 # ----------------------------------------------------------------------------- 1341 1342 ## Helper functions 1343
1344 -def makePowerSpec(v, p):
1345 if isinstance(p, int): 1346 if p < 0: 1347 raise ValueError("power should be > 0") 1348 elif p==0: 1349 return "1" 1350 elif p > 4: 1351 return Pow(v, p) 1352 else: 1353 res = v 1354 for i in range(p-1): 1355 res = res * v 1356 return res 1357 else: 1358 return Pow(v, p)
1359 1360 1361 # the following functions will be removed if the calling order for Quantities 1362 # is changed (as it ought to be some day). 1363
1364 -def makePar(parname, val=None):
1365 if val is None: 1366 return Par(parname) 1367 elif isinstance(val, str): 1368 gpar = Par(val, parname) 1369 else: 1370 return Par(repr(val), parname)
1371 1372
1373 -def makeFun(funname, sig, defn):
1374 return Fun(defn, sig, funname)
1375