sg-s / xolotl

A MATLAB neuron simulator. Very fast (written in C++). Flexible (fully object oriented). Immediate (live manipulation in MATLAB). Comes with a powerful parameter optimizer. Get started ➡️
https://go.brandeis.edu/xolotl
GNU General Public License v3.0
45 stars 8 forks source link

clamping currents computed using Ohm's Law, and do not respect custom getCurrent methods #555

Closed JRS92 closed 3 years ago

JRS92 commented 3 years ago

I would think that if you voltage clamp a compartment with a single conductance, the current density reported by integrate() (with output_type = 1) multiplied by the surface area of the compartment should be equal to I_clamp minus the capacitance artifact. However, this does not seem to be the case.

In this figure, I_clamp is plotted on the left, and the current density multiplied by the surface area of the single-conductance compartment is on the right. This is multiple sweeps of a simple voltage ladder.

image

The conductance uses a custom getCurrent() function. Removing that function (inheriting default getCurrent) makes I_clamp match density * surface area.

image

Thus, I believe this is something related to how current densities are calculated? I wonder if this also carries over to calculating the "dominant current" such as in this example: https://xolotl.readthedocs.io/en/master/tutorials/first-neuron/.

sg-s commented 3 years ago

this is interesting -- could you share your custom getCurrent() function?

JRS92 commented 3 years ago

double optiSlowKV5::getCurrent(double V_){ return g * (V_/25)* (exp((V_-E)/25) - 1)/((exp(V_/25) - 1)); }

sg-s commented 3 years ago

i'm still trying to understand the issue, but in the meantime, can you try your model with

x.use_current = 1;

where x is your xolotl model?

JRS92 commented 3 years ago

That does not change the output. Sorry, are you not able to reproduce the problem? I'm able to do this with any conductance...no matter what getCurrent looks like, I_clamp will always show the default g*(V-E).

JRS92 commented 3 years ago

Here is some example code. Swap that getCurrent into liu/ACurrent...or use any 'getCurrent' that you want as long as it's not the same as the default g * (V-E);

x = xolotl;
x.add('compartment','AB','A',.0017);
x.AB.add('liu/ACurrent','gbar',20);
x.dt = 0.1;
x.t_end = 100;
x.V_clamp = 40;
x.output_type = 1;
results = x.integrate;

plot(results.AB.I_clamp,'Color', 'k'); hold on
plot(results.AB.ACurrent.I * x.AB.A,'Color','m'); 
ylim([-0.5 2]); 

image

sg-s commented 3 years ago

aaah, thanks, i understand the problem now:

when voltage clamping, the clamping currents are calculated assuming Ohmic currents and do not respect the getCurrent method in the conductance....right?

this line is where the clamping current is calculated, and as you correctly guessed, it's Ohmic.

thanks for this bug report -- it's a nasty one. i'll fix it today.

sg-s commented 3 years ago

@JRS92 fixed, let me know if this works.

Screen Shot 2021-03-12 at 3 18 45 PM
JRS92 commented 3 years ago

Works great! Thanks