Adding conductance dynamics to ConductanceSystem

[AndreaRamirezH]

Hi @wsphillips,

I'm trying to model a neuron with conductances regulated by an integral control rule that uses the calcium concentration, [Ca]. See https://www.cell.com/neuron/pdf/S0896-6273(14)00292-X.pdf The new equations that I want to add are eqs (3) in the paper: D(g_i) = 1/tau_g (mRNA_i - g_i) D(mRNA_i) = 1/tau_i (Ca_{target} - [Ca])

I think that ConductanceSystem could contain the ODESystem with the equations for D(g) and D(mRNA) as well as D(gating variables) but I don't know how to implement these changes. What would be the best way to incorporate them into the code ?

Thanks!

[vi-hart]

Here is an implementation which should create the system for Figure 1B in the paper. The parameter values come from the supplemental material.

Currently there is no specialized constructor for adding a custom ODESystem into a ConductanceSystem but directly using the structure definition works.

Currently there isn't direct support for MTK unit validation #5 , so I am converting units "by hand" with ustrip.

import Unitful: ustrip, μM, μS, mS, ms, s, cm
import Conductor: GatingVariable, Linear

# Define variables & parameters
V = MembranePotential()  # mV
Ca2 = 109.2 * exp(V/12.5)  # μM
Ca_tgt = 1  # μM
@variables mRNA(t) g(t)
@parameters τg = ustrip(ms, 3600s) 
@parameters τm = ustrip(μM*ms/mS, 9.6e5μM*s/μS)

# Equation 3
eqs = [
    D(g) ~ (mRNA - g)/τg,
    D(mRNA) ~ (Ca2 - Ca_tgt)/τm
]

# Construct the ionchannel with conductance dynamics
name = :CaCh
CaCh = ConductanceSystem(
    g,
    Calcium,
    GatingVariable[],
    Set(V),
    ODESystem(eqs; name=name),  # control rule
    Linear,
    nothing,
    name
)

@named LCh = IonChannel(Leak; max_g=0.1μS/cm^2)

reversals = Equilibria([Calcium => 50mV, Leak => -85mV])
channels = [CaCh, LCh]

@named neuron = CompartmentSystem(V, channels, reversals)

[AndreaRamirezH]

OK, this works but only for those two leak conductances. What if I want to have the gating variables for all the ion channels though? Maybe modifying the equations within the function ConductanceSystem, like pushing the two new regulation equations to eqs, and adding g_i(t) and mRNA_i(t) as states of the ODE system, no?

[wsphillips]

I will check this later today

[wsphillips]

ok I see what's going on. We can fix this with a patch. We just need to allow max_g to be of type Num instead of expecting a real number value.

[wsphillips]

@AndreaRamirezH see #21 for a possible fix.

That being said: You can sidestep things like this by using the default struct constructors for ConductanceSystem (as @pjh6654 mentioned). What's missing from her solution was that you'd have to replicate all the other work normally done by the convenient front-end via IonChannel etc. So in addition to making your ODESystem with your extra equations, you would have to write ALL the equations since you're manually defining the object.

The solution in #21 right now extends the IonChannel constructor to accept ODESystem extensions and allow max_g to be of type Num. This means you can do something like:

using Conductor, IfElse, Unitful, ModelingToolkit
import Unitful: mV, mS, cm
import Conductor: t

Vₘ = MembranePotential(-65mV) # set default V₀ = -65mV
@variables gbar(t) b(t)
@parameters z

@named extras = ODESystem([gbar ~ 2b * 5z])

nav_kinetics = [
    Gate(AlphaBeta,
         IfElse.ifelse(Vₘ == -40.0, 1.0, (0.1*(Vₘ + 40.0))/(1.0 - exp(-(Vₘ + 40.0)/10.0))),
         4.0*exp(-(Vₘ + 65.0)/18.0), 3, name = :m)
    Gate(AlphaBeta,
         0.07*exp(-(Vₘ+65.0)/20.0),
         1.0/(1.0 + exp(-(Vₘ + 35.0)/10.0)), name = :h)]

@named NaV = IonChannel(Sodium, nav_kinetics, max_g = gbar, extensions = [extras]) 

Something we might do in the future is break up "typical" parts of the system construction process into an internal API that can be used to roll your own system constructor for custom applications. Note if you want even MORE flexibility, you can back up to the AbstractConductanceSystem type and even define your own data layout.

[AndreaRamirezH]

Hi @wsphillips, sorry it took me so long to reply.

Thank you for taking a look and directing me to the fix, it works well :) Thanks @pjh6654 also for posting your answer!

[wsphillips]

No problem. I'm glad it was a useful patch. I've been busy the past couple months but I plan to resume active work on this package either this week or next week. Feel free to follow up with more issues as you find them.