IClamp and ImpedanceRatio in python

When Python is the interpreter, what is a good
design for the interface to the basic NEURON
concepts.

Moderator: hines

Post Reply
rmaina
Posts: 2
Joined: Wed Feb 09, 2022 1:53 pm

IClamp and ImpedanceRatio in python

Post by rmaina »

Hello,

I'm relatively new to neuron and attempting to inject a current into a custom cell and record the transfer impedance along the axon as one would with the ImpedanceRatio tool in the Neuron GUI but in python. Attached is my python code, I built the same cell in the neuron GUI and recorded the transfer impedance and the outputs do not match. I am using this as a sanity check my next step is to Import an array from a file, convert it to a vector and play that vector on the Iclamp for the next step. Any help on this would be greatly appreciated.

Code: Select all

import neuron as ne
import nrn as n
import numpy as np

class cell:
    def __init__(self):
        self.ne = ne
        self.n = n


    def make_cell(self):
        soma = self.soma = self.n.Section()
        axon = self.axon = self.n.Section()
        dend = self.dend = self.n.Section()

        soma.L = 100
        soma.diam = 50
        soma.nseg = 5
        soma.Ra = 100
        soma.cm = 1

        axon.L = 1500
        axon.diam = 10
        axon.nseg = 50
        axon.Ra = 100
        axon.cm = 1

        dend.L = 3400
        dend.diam = 50
        dend.nseg =70
        dend.Ra = 100
        dend.cm = 1 

        axon.connect(soma,1)
        dend.connect(soma,0)

        soma.insert('hh')
        axon.insert('hh')
        dend.insert('hh')
        
        self.vect = ne.h.Vector()
        soma.push()
        self.vect.record(self.ne.h.ref(soma(0.5).v))
        self.ne.h.pop_section()

    def run(self,tstop,v_init=-65):
        print ("Before sim " + str(self.soma(0.5).v))
        self.ne.h.finitialize(v_init)
        self.ne.h.fcurrent()
        self.ne.h.dt = 1          
        z_transfer = []
        
        #load file 
        in_wav = np.loadtxt('ModulatedWave.txt')
            
        #IClamp 
        stim = self.ne.h.IClamp(self.soma(0.5))
        stim.delay = 0
        stim.dur = 1000
        #stim.dur = 1e9             #for vector input 
        stim.amp = 0.5
        
        #read transfer impedance
        imp = self.ne.h.Impedance()
        imp.loc(0.5, sec=self.soma)     #location of input electrode
        
        #input current
        #i_vect = ne.h.Vector(np.array(in_wav))
        t_list = list(range(0,1000))
        
        while self.ne.h.t < tstop:
            #ne.h.VecStim.play(stim._ref_amp, i_vect , 1)
            imp.compute(int(self.ne.h.t), 1)
            z_transfer.append(imp.transfer(0.41, self.axon))
            self.ne.h.fadvance()
            
        print ("After sim " + str(self.soma(0.5).v))
  
        z_transfer = np.array(z_transfer)
        output_f = open('channel_result.txt', 'w')
        np.savetxt(output_f, z_transfer)
        output_f.close()

c = cell()
c.make_cell()
c.run(1000)
ted
Site Admin
Posts: 6287
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: IClamp and ImpedanceRatio in python

Post by ted »

It looks like you're off to a good start, but let me suggest a slightly different strategy about what to import and how to structure code that sets up and exercises NEURON models. This borrows heavily from the "ball and stick" tutorial that starts at https://neuron.yale.edu/neuron/docs/bal ... del-part-1, and happens to contain many other useful tips. The tutorial shows a way to define cell classes which will produce cells that are ready for simulation on serial or parallel hardware.

For some time now it has not been necessary to import either neuron or nrn. It is generally sufficient to

Code: Select all

from neuron import h # gets NEURON's computational engine
h.load_file('stdrun.hoc') # gets high-level simulation control functions
or do this

Code: Select all

from neuron import h, gui
which gets everything those two lines do, PLUS NEURON's built-in GUI tools AND the NEURON Main Menu toolbar, which are handy for "almost code-free" model specification and debugging.

The question "how much stuff should be included in a cell class definition" often comes up. My opinion is that it's generally best to limit a cell class definition to specifying the anatomical and biophysical properties that are being represented in the model. The resulting cell classes are more modular and easier for others to reuse without first having to remove undesired instrumentation (current or voltage clamps, or code that sets up Vector play or record), graphical displays, or simulation control ("virtual experimental protocol") stuff.

End of sermon.

Here's the first part of something that I wrote based on the example you provided:

Code: Select all

from neuron import h
import numpy as np

class Cell():
  name = 'Cell'
  def __init__(self, gid):
    self._gid = gid
    self._morphology()
    self.all = self.soma.wholetree()
    self._biophysics()
  def __repr__(self):
    return '{}[{}]'.format(self.name, self._gid)
  def _morphology(self):
  # create sections, configure topology, specify geometry
    self.soma = h.Section(name='soma', cell=self)
    self.dend = h.Section(name='dend', cell=self)
    self.axon = h.Section(name='axon', cell=self)
    self.dend.connect(self.soma, 1)
    self.axon.connect(self.soma, 0)
    self.soma.L = 100
    self.soma.diam = 50
    self.dend.L = 3400
    self.dend.diam = 50
    self.axon.L = 1500
    self.axon.diam = 10
  def _biophysics(self):
    for sec in self.all:
      sec.Ra = 100
      sec.cm = 1
      sec.insert('hh') # hh doesn't affect Ra or cm
    self.soma.nseg = 5 # way too many--1, 19, and 17 should do for this model
    self.dend.nseg = 70
    self.axon.nseg = 50

cell = Cell(0)
Notice that
1. the class definition addresses only the anatomical and biophysical properties represented by the model cell
2. _morphology handles topology first, then the geometry of the sections (L and diam)
3. __init__ creates the SectionList (set of sections) all, which contains all of the sections in the cell
4. _biophysics takes advantage of the SectionList all. Also notice that it suggests nseg values that are based on the d_lambda rule (so that no segment is longer than 0.1 length constant at 100 Hz). This is generally sufficient for good accuracy.

Next come statements that set up instrumentation (stimulators, Vector recording, application of the Impedance class):

Code: Select all

# stimulators, Vector record and play, Z class
stim = h.IClamp(cell.soma(0.5))
stim.delay = 0
stim.dur = 1e3
stim.amp = 0.5
# in_wav = np.loadtxt('ModulatedWave.txt')
# needs additional code to use in_wav to drive cell
# set up recording of soma(0.5).v and t
vvec = h.Vector()
vvec.record(cell.soma(0.5)._ref_v)
tvec = h.Vector()
tvec.record(h._ref_t)
# set up calc and recording of transfer impedance
imp = h.Impedance()
imp.loc(0.5, sec = cell.soma)
zxvec = h.Vector()
Finally comes the code that specifies and executes simulation flow control (AKA "the virtual experimental protocol"):

Code: Select all

def dorun(cell, tstop, v_init):
  h.dt = 1 # too long for hh, ok for test
  zvec = h.Vector() # discard pre-existing values
  h.finitialize(v_init)
  print ("Before sim " + str(cell.soma(0.5).v))
  while h.t < tstop:
    # really want to calculate z at 0, 1, 2 . . . Hz?
    imp.compute(int(h.t), 1)
    # if 0.41 is a location of particular interest,
    # should instead break axon into two pieces
    # so that the location remaines fixed regardless of nseg values 
    zxvec.append(imp.transfer(0.41, cell.axon))
    h.fadvance()
  print ("After sim " + str(cell.soma(0.5).v))
  # zxvec is already accessible as a numpy array
  output_f = open('channel_result.txt', 'w')
  np.savetxt(output_f, zxvec) # also save corresponding t?
  output_f.close()

dorun(cell, 1e3, -65.0)
Let me know if you have any questions about this code or the comments.
rmaina
Posts: 2
Joined: Wed Feb 09, 2022 1:53 pm

Re: IClamp and ImpedanceRatio in python

Post by rmaina »

Thanks for the advice, I'll go ahead and implement those changes and give it a shot.
Post Reply