Update NeuroML2 version of muscle model to match Neuron version

[pgleeson]

Will try to match behaviour at https://github.com/openworm/muscle_model/tree/master/neuron_implementation

[vellamike]

@pgleeson I suggest we do this next time we're in the same room - It would be very educational for me.

[vellamike]

@pgleeson I suggest we do this next time we're in the same room - It would be very educational for me.

[pgleeson]

Ok, it depends when you need it done by, as we'll probably not meet for a few weeks. I'll do some initial restructuring of the files in https://github.com/openworm/muscle_model/tree/master/NeuroML2 in the time being.

What was the decision to go with the Na, K and Ca N type channels based on? It's a fast and a slow K channel and a Ca chan in the Boyle and Cohen paper.

[vellamike]

Those mechanisms are mis-named really. The kinetics were optimized automatically so exactly how you interpret the simulation in terms of ionic mechanisms requires some thinking. Looking at this, I need to clean up the code a bit and make some stuff more clear.

On 3 January 2014 15:51, Padraig Gleeson notifications@github.com wrote:

Ok, it depends when you need it done by, as we'll probably not meet for a few weeks. I'll do some initial restructuring of the files in https://github.com/openworm/muscle_model/tree/master/NeuroML2 in the time being.

What was the decision to go with the Na, K and Ca N type channels based on? It's a fast and a slow K channel and a Ca chan in the Boyle and Cohen paper.

— Reply to this email directly or view it on GitHubhttps://github.com/openworm/OpenWorm/issues/169#issuecomment-31530800 .

[slarson]

Seems to @rayner and I that it makes sense to test the current instantiation of this model against the paper. He's gotten it to run. We've had a look at the activity of the model and compared it against the paper.

This figure shows the overall behavior of the real muscle and the model:

We tested the current implementation in the repo and we find that for the 20ms stimulus time used in the figure that the results are not the same. The system seems to charge up to produce a "spike" no matter what (even with no input) but in the figure some amount of input is needed to produce activity.

What we need is a system that at rest does not "spike". This means that the channels are not the same as in the paper.

Another figure in the paper looks at the channels:

Here I-V curves are generated for the main channels and compared against real data. In order for us to determine if our versions of those channels are correct, we need to be able to create our own I-V curves.

@pgleeson or @vellamike -- do you have a magical bit of code somewhere that nicely shows the I-V plots within the NEURON / NeuroML ecosystem that we could apply here, or should we try to make it up ourselves?

[pgleeson]

You could look at some of the python scripts under development here: https://github.com/OpenSourceBrain/BlueBrainProjectShowcase/tree/master/Channelpedia.

There is one script NML2ChannelAnalyse.py which can generate a LEMS file for plotting some of the key graphs associated with a channel, e.g. rate activation curves, conductance responses to voltage clamps, etc. See below for an example

LEMS plots in grey around the outside, some generated data from the equivalent channel on Channelpedia on the inside.

You could also browse Chanelpedia for similar/interesting channels and use test_all.py to generate & test NML2 equivalents (work in progress, so test visually against plots on Chanelpedia...).

Regarding the I-V curves above, these are plots are harder to generate via the scripts above, as they involve the steady state of the current vs V, which is not a readily accessible quantity from the LEMS (while minf, hinf & the instantaneous currents are available, minf^3*hinf is never calculated). However, you could do some post processing on the m_inf.lems.dat files that are produced.

[slarson]

@pgleeson very helpful, thank you! This would be something good for @rayner to have a look at and play with. We've also been looking again at the wiki page for the muscle model and realizing that the graphs I asked for at the top were already created via jnml over here:

http://opensourcebrain.org/projects/muscle_model/wiki

@pgleeson what is the recommended status of jnml today? Should it be subsumed into jlems now? Should it be accessed mainly via python scripts? Is the python neuroml library able to do what jnml did before? I'm getting confused between all the cross products between Python/Java and LEMS/NeuroML.

[slarson]

@VahidGh is going to have a look at this too and see if he can start modeling channel EGL-19, an L-Type Calcium channel, which pertains to the muscle cell. @VahidGh feel free to throw in a screen shot in here or ask additional questions if you run into trouble.

[pgleeson]

Just has a look at that wiki page and saw it needed some updates (as well as a minor change to the example scripts). A new pull of the muscle model (and a rebuild or update of jNeuroML) should give you working examples which match the traces in the latest wiki.

Regarding @slarson's questions: What is the recommended status of jnml today? Stable, but still in active development. The version from the development branches (or the precompiled jar that can be installed with svn checkout svn://svn.code.sf.net/p/neuroml/code/jNeuroMLJar) is always the most well tested and the one I use for all my work. Should it be subsumed into jlems now? I often use these (jneuroml/jnml/jlems) interchangeably, but jNeuroML effectively wraps a standalone package jLEMS, and adds the NeuroML ComponentType definitions and some import/export options (e.g. to NEURON), and other bells and whistles like validating NML files (try jnml -h). Both are in active development and jnml/jlems can be used synonymously, except when jlems is used to run a pure, non neuronal LEMS file. When Geppetto talks of using jLEMS, it usually means jLEMS & some of the other NeuroML packages that go to make up jNeuroML. Should it be accessed mainly via python scripts? Not necessarily, the command line jnml utility is still the main way of running a LEMS file. Python scripts using libNeuroML are very useful (e.g. in c302) for generating NeuroML2 files, but these are still normally run with jnml. There is a pure Python library for parsing & executing LEMS, PyLEMS, but this is slower than jnml & doesn't have all the options like code generation etc. Is the python neuroml library able to do what jnml did before? Mostly, but these don't really replace jnml. libNeuroML is mainly for reading & writing NeuroML, and can be used with PyLEMS for a pure Python environment for all your NeuroML2/LEMS needs (see libNeuroML and PyLEMS: using Python to combine procedural and declarative modeling approaches in computational neuroscience). However jNeuroML will usually be required when you want to run your generated model or convert it to other simulator formats.

[rayner]

I forked the muscle_model and BlueBrainProjectShowcase repositories and added a few changes to help visualise the current/voltage relationship for the ion channels in the muscle cell.

To use it:

1. Install jNeuroML
2. Clone my forks of the muscle_model and BlueBrainProjectShowcase repositories.
3. Using the NML2ChannelAnalyse.py script in BlueBrainProjectShowcase/Channelpedia/, create a LEMS file from the muscle_model NeuroML files: e.g. python NML2ChannelAnalyse.py [muscle_model install dir]/NeuroML2/k_fast.channel.nml k_fast
4. Run the simulation using jNeuroML: jnml LEMS_Test_k_fast.xml. This will now create data files of channel currents over time for each of the clamp voltages in the simulation (e.g. i_100.lems.dat for 100mV, i_min80.lems.dat for -80mV).
5. Now you can run:

python i_plot.py filename.lems.dat (to plot the current against time, replicating one of the graphs produced by running the LEMS file through jNeuroML)

python iv_analyse.py (to find the peak current values from all the current data files, and create a plot of peak current against voltage)

To do:

• Find out why the ca_boyle channel doesn't produce a working LEMS file (the current variable appears to contain an invalid value). The NeuroML file contains a warning that it doesn't yet work as intended, which may be related.
• Calculate steady-state currents as well as peak ones (fit curve to the current values for each clamp voltage, then find where the gradient becomes zero?)
• Find out if there's any way to produce the I/V plots when running the sim through jNeuroML in the first place, so we can view it alongside all the other graphs.

[rayner]

N.B. dependencies for the BlueBrainProjectShowcase code are listed in https://github.com/rayner/BlueBrainProjectShowcase/blob/master/.travis.yml

[pgleeson]

Haven't forgotten about this... It should be possible to get some kind of I-V curve out of NML2ChannelAnalyse.py, but this will require some refactoring/tidying up (& maybe moving to separate repo of this analysis script). Slowly working on it...

By the way, the ca_boyle channel does seem to produce a working LEMS file from NML2ChannelAnalyse.py when I use the original muscle model repo & BlueBrainProjectShowcase version

[slarson]

@theobrown Have a look at this issue. Give a shot at building @rayner 's versions of this and see where you get

[slarson]

Update:

Our IV plot for ca_boyle doesn't look like the one from the paper. So we've been looking at the currents that come out for different levels of voltage:

two things about this. First, we aren't plotting for voltage below zero -- we need to do that. But even for the values above zero, we aren't seeing the features we would expect from the boyle paper:

One thing we need to check in on is the scaling of the Y-axis -- are we off by orders of magnitude or is this the right scale?

Second, currents seem to be getting increasingly negative as voltage is increased, but instead they should be getting more positive.

As a next step we should compare the parameters we are using for this channel with what are in the paper.

[pgleeson]

I've added scripts analyse_k_fast.sh and analyse_k_slow.sh to make it easier to see the properties of the k channels. This requires the latest code from https://github.com/OpenSourceBrain/BlueBrainProjectShowcase/tree/master/Channelpedia.

Also, the latest commit adds a custom gate for Ca channel inactivation, plus an internal decaying Ca pool. Not all the parameters are correct, but it's producing better plots for Figure 2 in the paper:

Discuss this more during the hangout.

[slarson]

@pgleeson this is really really terrific. @rayner let me know if you've had a chance to play with this yet. Would be great to see the next generation I/V curve after these changes.

[rayner]

When I run the analyse_ca_boyle.sh script, I get the error:

SEVERE: (ERROR) checking DerivedVariable tau val=0 in ComponentType, name=customHGate (Custom gate for h) extends gateHHtauInf alpha (none) k (concentration) ca_half (concentration) instances (none) Expression dimension does not match delcared dimension expression gives ExprDimensional[dimensionless] but delcared is Dimension[time t=1]

I can work around this by changing the DerivedVariable definition for "tau" so that it has dimensions of "time", but then I get this instead:

SEVERE: (ERROR) Error when initialising StateType, id=h, Type=customHGate org.lemsml.jlems.core.run.RuntimeError: Problem getting exposed var caConc in: net1[networkWithTemperature]

I'm not sure where the missing variable exposure is coming from.

I suspect there might be something I need to update. I've already tried upgrading jNeuroML to version 0.5.2, though. I'll keep looking.

In other news, the Boyle & Cohen code compiles and runs. It's necessary to create a "results" directory under muscle_model/BoyleCohen2008/, though, otherwise it seems to silently discard the results files.

[pgleeson]

The analyse_ca_boyle.sh script predates the updated ca channel, which contains the [Ca2+] dependence. There is no support yet for channels which depend on Ca in the NML2ChannelAnalyse.py script. A solution would be to add a settable caConc parameter in vClampedCell here: https://github.com/OpenSourceBrain/BlueBrainProjectShowcase/blob/master/Channelpedia/LEMS_Test_template.xml (see here for an example of caConc usage https://github.com/NeuroML/NeuroML2/blob/master/NeuroML2CoreTypes/Cells.xml#L369)

You can give that a try if you like, but if you have analyse_k_fast.sh working & can generate the plots above, you have everything that's currently implemented.

[rayner]

Yes, analyse_k_fast.sh and analyse_k_slow.sh are working. Thanks for the explanation - I'll have a go at adding a caConc parameter to the template as you suggest.

[VahidGh]

@slarson playing with these codes, I could generate currently available ca_boyle, etc plots, but as pgleeson explained, there is no complete support for Ca channels yet. For egl-19 which codes l-type voltage-dependent Ca channel in muscle cells, I've almost extracted experimental data from this paper which also has been referenced in the Boyle & Cohen paper. But for the final tests, I need to plot the I-V curve to compare with figures in the paper, which is not possible right now. Should I start working with another channel or be of help on Ca and I-V curves plotting issues?

[slarson]

@VahidGh great! @rayner can you help @VahidGh generate the I-V curves for tha Ca channel? I've been slow on updating the readme but I have seen that we can do this.

In general @VahidGh, if you are able to get stuff working we could use your help in documenting your process for how you got it working and then offering that as a pull request to the README and/or the docs so others can follow your path.

[VahidGh]

@slarson Actually I could generate the I-V curves for Ca channel with data from both papers with a little python script. but the results are not the same (trying to find the problem):

@pgleeson converting k_f from 5 to -5 did not solve the problem. resulting figures were completely deformed, instead playing with V_0.5_f (from 25.2 to 5.2) generated a better result:

[pgleeson]

Nice work @VahidGh. As discussed with @rayner, we have all the original matlab scripts for the model here. I'd check through them to make sure they reproduce the figure from the paper. I've just updated the iv.m, and it seems to match what you need:

The Python version of this model https://github.com/openworm/muscle_model/tree/master/BoyleCohen2008/PythonSupport/Main_Version can be updated to match the matlab version (i.e. make a iv.py, based on your script, but reusing the input_vars.py).

The next step would be to make sure the NeuroML parameters and the input.csv values are the same. I've started this in compareToNeuroML2.py and it would be great if someone could help with this.

[VahidGh]

@pgleeson Thank you, for your suggestions. I added the iv_all.py with the help of current codes and using the input.csv values as you suggested. Finally I could generate these plots (for all voltage values and skipping optimization, etc for now):

But as you can see, the peak current values for Ca, does not match the matlab version (although I couldn't find how does the matlab version generating those plots (if the plots generated from the iv.m script with just an e^2 for I_Ca)).

[slarson]

This effort continues over at https://github.com/openworm/muscle_model