Gap junction network not seeming to conduct

[swienbar]

Hi,
I am trying to make a network of passive tubes connected all on one side by gap junctions. I am hoping to have one central tube and the rest on the periphery.
When I attempt to connect them with gap junctions, I can record the gap junction voltage 'vgap' and that is as expected, but it doesn't seem to flow through to the other side of the gap junction. Am I missing a step here?

Code: Select all

from neuron import h
from neuron.units import ms, mV
from neuron import gui
import matplotlib.pyplot as plt
from math import sin, cos, pi
h.nrn_load_dll("C:\\Users\\Sophia Wienbar\\Documents\\ConeModel\\mod\\nrnmech.dll")

class ConeObj:
    def __init__(self, gid, cellLength):
        self.gid = gid
        self.cellLength = cellLength
        self.__morphology()
        self.__conductances()
    
    def __morphology(self):
        self.tube = h.Section(name='tube', cell=self)
        self.tube.L = self.cellLength
        self.tube.diam = 1
        self.all = [self.tube]
    
    def __conductances(self):
        for sec in self.all:
            sec.Ra = 100    # Axial resistance in Ohm * cm
            sec.cm = 1      # Membrane capacitance in micro Farads / cm^2
        self.tube.insert('pas')
        for seg in self.tube:
            seg.pas.g = 0.001  # Passive conductance in S/cm2
            seg.pas.e = -65    # Leak reversal potential mV 
    
    def __repr__(self):
        return 'ConeObj[{}]'.format(self.gid)

def ConeNetwork(N, cellLength, Gj):
    # construct a circle of N - 1 cones around a central cone
    rad = 2 # um, radius of circle
    cones = []
    gapJ_toCenter = [0]*(N-1)
    gapJ_fromCenter = [0]*(N-1)
    # append central cone
    cones.append(ConeObj(0, cellLength))
    # make peripheral cones
    for i in range(N-1):
        cones.append(ConeObj(i+1, cellLength))
    
    # connect the central cone to the outer cones via gap junctions
    if N > 1:
        for j in range(N-1):
            print(j)
            # connect central cone to peripheral
            cones[0].tube.push()
            gapJ_fromCenter[j] = h.gap(cones[0].tube(0.999))
            gapJ_fromCenter[j].g = Gj # uS
            h.pop_section()
            cones[j+1].tube.push()
            gapJ_toCenter[j] = h.gap(cones[j+1].tube(0.999))
            gapJ_toCenter[j].g = Gj # uS
            h.pop_section()
            
            gapJ_fromCenter[j]._ref_vgap = cones[0].tube(0.999)._ref_v
            gapJ_toCenter[j]._ref_vgap = cones[j+1].tube(0.999)._ref_v
    return cones, [gapJ_fromCenter, gapJ_toCenter]

[Cones, GapJ] = ConeNetwork(6, 500, 10) # 6 cones that are 500 um long with Gj 10 nS

# define stimulus at one end
stim = h.IClamp(Cones[0].tube(0))
stim.delay = 2 # ms
stim.dur = 1 # ms
stim.amp = 0.1 # nA

CenterBeginV = h.Vector().record(Cones[0].tube(0)._ref_v)
CenterEndV = h.Vector().record(Cones[0].tube(1)._ref_v)
CenterGJ = h.Vector().record(GapJ[0][1]._ref_vgap)
NeighEndV = h.Vector().record(Cones[1].tube(1)._ref_v)
NeighBeginV = h.Vector().record(Cones[1].tube(0)._ref_v)
t = h.Vector().record(h._ref_t)

# begin stim
h.finitialize(-65 * mV)
h.continuerun(10 * ms)

# plot results
fig = plt.figure()
plt.plot(t, CenterBeginV, 'b')
plt.plot(t, CenterEndV, 'r')
plt.plot(t, CenterGJ, 'k')
plt.plot(t, NeighEndV, 'g')
plt.plot(t, NeighBeginV, 'm')
plt.ylabel('Voltage (mV)')
plt.xlabel('Time (ms)')
plt.show()

I am using the following mod file for the gap junctions

Code: Select all

NEURON {
  SUFFIX gap
  POINT_PROCESS gap
  RANGE g, i
  POINTER vgap
  NONSPECIFIC_CURRENT i
}
PARAMETER { g = 1 (nS)}
ASSIGNED {
  v (millivolt)
  vgap (millivolt)
  i (nanoamp)
}
BREAKPOINT { i = (vgap - v)*g*(1e-3) }

Thank you so much in advance!

[ramcdougal]

Note that vgap in your mod file is supposed to be the voltage of the OTHER half of the gap junction... but here, all the fromCenters are on cones[0] and the other side is also set to be cones[0], so your vgap and v are the same, hence no current flow. Likewise for toCenter.

That is, change it to:

Code: Select all

            gapJ_fromCenter[j]._ref_vgap = cones[j+1].tube(0.999)._ref_v
            gapJ_toCenter[j]._ref_vgap = cones[0].tube(0.999)._ref_v

As an aside, you don't need either of:

Code: Select all

cones[0].tube.push()
h.pop_section()

The second undoes the first, but other than that, they serve no role here as all section information is covered in the segment passed to h.gap. (Likewise for the corresponding cones[j+1] push and pop.) Basically, except when dealing with legacy code, there's no reason to use a section's push method or h.pop_section.

[swienbar]

That worked! Thank you so much!

[ted]

Very nice, but there still may be a reason to push and pop the section stack: specifically, spatial discretization. A 1 um section with length 500 um, cm 1 uf/cm2, and Ra 100 ohm cm must be represented by many segments for the sake of accuracy in space. And this will also affect accuracy in time.

The d_lambda rule, described in

Hines ML, Carnevale NT. (2001). NEURON: a tool for neuroscientists. The Neuroscientist, 7(2), 123-135.

generally gives good results, and is easy to implement. Comment to myself and whoever else cares: a preprint of this paper is currently available here http://www.neuron.yale.edu/neuron/stati ... e_rev2.pdf but it should also be made available at a more permanent site, e.g. should be downloadable from nrn.readthedocs.io

Edit added a bit later: And so it is (albeit via a somewhat more obscure URL)! See https://nrn.readthedocs.io/en/latest/_d ... ntists.pdf

WRT implementation: execute the following
AFTER all sections have been created and their anatomical and biophysical properties have been specified
but
BEFORE any point processes (including gap junctions) or NetCons have been attached

Or one might prefer to include similar* code in every cell class definition so that it is executed every time a new cell instance is created. (*--similar but affecting only the new instance's sections. This is left to the reader as an exercise.)

Code: Select all

# for each section assign the smallest odd nseg such that 
# nseg > L / (0.1 * AC length constant at 100 Hz)
# this uses lambda_f(), a function defined in stdlib.hoc
# which is loaded when one imports gui from neuron
# or executes h.load_file("stdrun.hoc") or h.load_file("stdlib.hoc")
for tmpsec in h.allsec():
  tmpsec_id = h.this_section(sec=tmpsec)
  h.push_section(tmpsec_id) # so lambda_f uses current section's Ra, cm, diam
  tmpsec.nseg = 1 + 2*int((tmpsec.L/(0.1*h.lambda_f(100))+0.999)/2.0)
  h.pop_section()