Page 1 of 1

Simultaneous equations in python

Posted: Wed Aug 30, 2017 5:53 am
by catharina
Hi,

I want to construct a network of pyramdial cells and interneuron. Therefore I have chosen to use a cell template similar to this one https://senselab.med.yale.edu/modeldb/S ... del=222321. Right now I am trying to construct cell templates in python that are similar to the original ones but use the same mechanisms. However I am having trouble including the mechanisms for the Na-K pump nakpump.mod. This mechanism needs a stability function that recaclulates sodium and potassium conductances and the maximal flux of the pump. This is done by defining simultaneous equations. The function is called once before the simulation is run. Here is the code of the original implementation in hoc:

Code: Select all

proc pyramstab() {
	access psoma
	rm=$1
	init()
	vrest = -70

	dr_na = vrest-ena
	dr_k = vrest-ek

	eqinit()
	depvar gna, gk
	eqn gk:: gk = 1/rm - gna
	eqn gna:: gna = (-dr_k*(3/2)*gk)/dr_na
	solve()
	forall gk_leak=gk
	forall gna_leak=gna

	p_k=(1+km_k_nakpump/ko)^(-2)
	p_na=(1+km_na_nakpump/nai)^(-3)
	p_tot=p_k*p_na
	forall imax_nakpump = (-gna_leak*dr_na)/(3*p_tot)
}
I have tried to write my own version of this function in python. However because the eqn, eqinit() and solve() is not available in python, I do solve the equations and call my function at every time step before I call h.steprun(). This is the code:

Code: Select all

def stab(celllist):
    rm=5e3
    vrest = -70

    for nn in range(len(celllist)): # for each cell calculate and assign gk_leak, gna_leak and imax_nakpump
        dr_na = vrest-celllist[nn].soma.ena
        dr_k = vrest-celllist[nn].soma.ek

        gk=(2*dr_na)/((3*dr_k-2*dr_na)*rm) # equations solved for gk and gna, see EqnSolve.nb
        gna=(3*dr_k)/((3*dr_k-2*dr_na)*rm)
        for ss in celllist[nn].sections:
            ss.gk_leak=gk
            ss.gna_leak=gna

        p_k=(1+celllist[nn].soma.km_k_nakpump/celllist[nn].soma.ko)**(-2)
        p_na=(1+celllist[nn].soma.km_na_nakpump/celllist[nn].soma.nai)**(-3)
        p_tot=p_k*p_na
        for ss in celllist[nn].sections:
            ss.imax_nakpump = (-ss.gna_leak*dr_na)/(3*p_tot)
    return
The problem is that this appaerently does not work the same way as the implementation in hoc does. Now, I am aware that it is because the equations are not defined as simultanous ones in NEURON and hence are not solved parallel to the cable equation. Does anyone know how I can define simultaneous equations in python or how I could work around it? If anyone knows an example where this was done, it would also help me?

I would appreciate any help!

Best, Cathi

Re: Simultaneous equations in python

Posted: Wed Aug 30, 2017 11:31 am
by catharina
It took me a while but I finally found a way to work around it.

That's the code I add in python:

Code: Select all

h('objref celllist,N')
h.celllist=h.List()
for elem in celllist:
    h.celllist.append(elem)

h('proc hstab() {\n\
    for nn=0,1 {\n\
    	rm=5e3\n\
	vrest = -70\n\
    //print nn\n\
	dr_na = vrest-celllist.object(nn).soma.ena\n\
	dr_k = vrest-celllist.object(nn).soma.ek\n\
	eqinit()\n\
	depvar gna, gk\n\
	eqn gk:: gk = 1/rm - gna\n\
	eqn gna:: gna = (-dr_k*(3/2)*gk)/dr_na\n\
	solve()\n\
    for ii=0,2 {\n\
    celllist.o(nn).sections[ii].gk_leak=gk\n\
    celllist.o(nn).sections[ii].gna_leak=gna}\n\
	p_k=(1+celllist.object(nn).soma.km_k_nakpump/celllist.object(nn).soma.ko)^(-2)\n\
	p_na=(1+celllist.object(nn).soma.km_na_nakpump/celllist.object(nn).soma.nai)^(-3)\n\
	p_tot=p_k*p_na\n\
    for ii=0,2 {\n\
    celllist.o(nn).sections[ii].imax_nakpump= (-celllist.o(nn).sections[ii].gna_leak*dr_na)/(3*p_tot)\n\
}}}')
It's not nice but it works. =D Still I think it would be kind of nice to know how to do that directly from python.

Re: Simultaneous equations in python

Posted: Wed Aug 30, 2017 10:20 pm
by ted
hoc's simultaneous equation feature is not required to deal with a system of two linear equations in two unknowns--a bit of algebra is sufficient to solve first for gk, then for gna (or vice-versa).

Re: Simultaneous equations in python

Posted: Thu Aug 31, 2017 5:45 am
by catharina
Yes, that's what I tried first, see the second block of code I posted. However this doesn't give the same results as the hoc functions. So in this case the stability functions in hoc do ensure that the resting potential stays constant, but the python function doesn't.

Re: Simultaneous equations in python

Posted: Thu Aug 31, 2017 9:22 am
by ted
catharina wrote:However this doesn't give the same results as the hoc functions. So in this case the stability functions in hoc do ensure that the resting potential stays constant, but the python function doesn't.
You're saying that your model implementation in Python produces the same results as the original hoc implementation does, as long as your Python implementation uses the "stability functions" implemented in hoc?