Standards in compartmental modeling / Libraries of mechanisms

Dear all,

Just to give some background. I recently moved my compartmental modeling work to brian2, most of this material was originally implemented in Neuron and shared publicly (e.g. on ModelDB). It works great (thanks !!), but I feel that we would benefit from discussing the way to write mechanisms, trying to establish standards. The equation string is very convenient to design custom models, but it gets very long (and poorly readable) when adding many biophysical mechanisms. Having a Neuron-like way to include mechanism would make sense in this research topic (I think).

I share below a draft implementation, it’s focused on the mechanisms used to study dendritic integration in cortical cells (e.g. see Branco et al., Science 2010). I am far from sure that this is a good way to do it. I’d be interested in having some feedback on this.

Ideally, we could also have a way to share mechanisms, we could set up a library of biophysical mechanisms (part of the brian2 repo or somewhere else).
What do you think ?

Hi @yzerlaut. This looks like you’ve done some really cool work! I don’t think we’ll include it as part of Brian itself because our philosophy is to only include core stuff as much as possible now, and put everything else in separate packages. It would be great if you made a package out of this though!

@dan this looks like a very clever and handy approach. As soon as one departs from a single-compartment model with a few currents, the equations stuff become horribly messy. Maybe we can spin up a library for Brian extensions?

1 Like

Hi Yann, I just realized that I did not switch on notifications for this category so I’m only now seeing your post :upside_down_face: This is really some impressive work!
I agree with Dan that we are probably not going to integrate this in Brian itself directly. One reason is that a list of models like this is of course always arbitrary and we want to avoid discussions about, say, which variant of the Hodgkin-Huxley model we should include. The second reason is that there is a danger with built-in mechanisms (or mechanisms distributed databases): since they are easy to use, most users will use them regardless of whether they are appropriate or not, and without understanding the details of their implementation. An example is the number of NEURON models that simply use insert hh for a “biologically-realistic neuron model”, forgetting that this is a model of a squid axon, or the use of insert pas + insert hh, even though the hh mechanism already contains a passive current. Of course you could argue that this is just a mistake by the user, but software design choices can make certain types of mistakes more or less likely. If you are interested in a longer form of this argument, have a look at this blog post by @romainbrette.

I hope all that doesn’t sound too dismissive of @yzerlaut’s work! I do agree with @yzerlaut and @rth that the equation string approach quickly reaches its limits with complex multi-compartment models. And I think recently Brian is used by more and more researchers doing this kind of work, so I’d love to improve its capabilities in this domain.

From my point of view, our best way forward consists of two parallel lines of work:

  1. Improve Brian to make writing extensions like the one presented here more easy, and in general allow for more structured equations that make them easier to read. E.g. we discussed having support to add suffixes or a templating syntax to the Equations object directly, which would be a bit more elegant (additional checks, etc.) than using string formatting. See this github issue (but note that the proposal changed quite a lot in the discussion). Another thing that would be nice to have is macros in equations so that you can make it clearer when e.g. different gating variables use the same rate equation just with different parameters (see this github issue).
  2. Have some kind of standardized way to make models like the one in @yzerlaut’s code accessible. For me, this does not necessarily have to be in the form of a library from which you can import code, but it could also be simply presenting them on a website in a standardized and copy&paste-friendly way. So similar to what @yzerlaut already has in his code, we’d have a short description + paper link for each model, together with the model code and maybe a plot showing the model in action. We discussed something like this 7 years ago, actually… We could also autogenerate some of this from NeuroML descriptions and/or provide NeuroML descriptions in addition to the Brian code. @Vigneswaran has been working on automatic model descriptions over this summer, I could imagine he would be also interested in contributing to efforts like this.

What do you think about all that? As you can see from the linked github issues that have been around forever, we are lacking the resources do implement all this by ourselves, but I think there are some great opportunities for contributions here. Everything around Equations is relatively accessible compared to the code generation parts of the code, I think.

I have a few random thoughts on the topic. Hopefully, some of them can be useful.

Re 1: What makes the equations messy is that there is no wave to separate them into logical blocks. What @yzerlaut 's code does, is just that separation. Even if equation parser can accept variables with name separation like d nap.m/dt = messiness can be significantly reduced. But the true modularity will be better.

Obviously, there may be different opinions here. I rather prefer to have cleaner code, that is easy to read/development than prevent possible mistakes from a user. What can be different from old NEURON’s paradigm, is that Brain can perform a better analysis of modules and warns a user. For example, if Brain will organize equations into a hierarchical tree, there may be a possibility to check equation similarity, therefore if a user insert pas and insert hh, Brain could simply indicate: “pas.g and hh.pas.g seem the same variables”. NEURON cannot do that just because of internal architecture, I guess. I should admit, I did such mistakes a lot in the past, I wish, NEURON could catch them.

Re 2

we’d have a short description + paper link for each model

That will be something like ModeDB - a useful, but a very old-school tool. Everyone will have to make its own way and reads through lots of line of code and (maybe) understands what is inside. So the quality of the code will be on the author of the code. In contrast, a carefully designed standard library can help to avoid user mistakes, simplify code, and speed up new module development a lot, IMHO.

1 Like

Hi. Thanks for your feedback @rth. Before I give my 2 cents, just a general remark because I think it might not have been clear enough from my earlier post: I am giving here my opinion (and to some extent of the Brian core development team) about what our ideal solution would look like. We certainly do not want to discourage users from using @yzerlaut’s work (or even worse, discourage @yzerlaut from packaging/distributing it). If there is something that we can do in the short term (my outlined approach is rather a long-term one), please let us know. Discussing the library here is already a first step, but we could also mention it on for example.

Ok, back to the points that @rth raised.

I think we are already not far from a solution for this. We can have one Equations object per channel, for example, and in the end add them all together. The main problem is the “scoping” issue, e.g. how to avoid clashes between two gating variables “m”, for example. In the discussion on github I linked to earlier, our approach was to deal with this via suffixes/prefixes, similar to @yzerlaut’s code. I see the appeal of the point syntax as well, but that would be a much bigger change (not only for the equation parsing, but also things like attribute access). I have to think about this a bit more.

There’s certainly some nice analysis we could do, but I don’t think this is straightforward. In NEURON’s example, the variables are called i, g, and e for the pas mechanism, but il, gl, and el in hh. They do have the same functional form in both, but that will be the case for a lot of gating variables across mechanisms as well.

What I have in mind is different from ModelDB. ModelDB is about complete models that were used in publications. What I had in mind was rather something that presents a certain neuron/synapse model, or a single channel type. So rather closer to something like Channelpedia. But then, it might not be worth the effort (or possibly) to duplicate the work that went into these initiatives. Maybe focussing on some automatic import/translation from NeuroML (which we already have to some extent) makes more sense.

The thing is, do we have any standard models? And should it come with parameters as well, even though we know that these strongly differ across neurons of the same type? And will the standard model change when a new paper comes out that updates the model? As you can see, I am sceptical that this will be manageable :smile: And I guess there’s a reason that NEURON does not come with a standard library after all these years…

Again, all this is my opinion and I am very interested in hearing other view points. I will resurrect the old github issue I linked to above and use @yzerlaut’s code as a reference for the kind of problem we’d have to solve.

Dear all.
Thank you very much for the feedback ! Very much appreciated ! (ans sorry for the late reply).

This is a very interesting discussion, I can totally see that it points to old discussions and I agree with you @dan and @mstimberg, it doesn’t need to be included within the software. My motivation was really just about finding a good/standard/robust way to implement the mechanisms, and find the way to make this implementation as “Brian2-ish” as possible.

So @mstimberg your proposal sounds very nice to me, working on those two parallel lines would be great. 1) some (minor) tweaks in the Equation object, e.g. the EquationTemplate of your Issue’s discussion with @dan seems very close to what I would like to have and 2) just releasing the standards of way to implement those mechanisms. For me it can be of any form: a dedicated website, even just the “from_papers/” directory of the repo would be good already. And yes, benefiting from the NeuroML material would indeed be very valuable !
For point 1), there might also be indeed be some conceptual work to do as well. E.g. you pointed out:

The main problem is the “scoping” issue, e.g. how to avoid clashes between two gating variables “m”,
Even more tricky things should be caught, like: “[Ca2+]-dependent currents should be inserted only with a model for Calcium dynamics (pump, …)”, or the redundancy problem of @rth, etc … It might be good to have a place to list the constraints that we would like to be caught in an ideal implementation, maybe we could start an “Issue” for that ?

On my side, I’ll set up a github repository repository with this material, I will add a few additional things that might be useful for cellular physiology in the cortex (e.g. pointing to nicely formatted morphology datasets that load straight into brian2, …).

I will post the link of the repo here if useful for anyone (I’ll try to do this soon).

And just sharing some thoughts on the discussion about the Frankenstein models. Personally, I totally agree with the opinion of @rbrette. In my scientific work in general, I would personally prefer to stick to study-specific low-dimensional models. But other investigators put more emphasis on “you can not look at this phenomenon without the full high-dimensional picture, so please use the Frankenstein model ! And yes, if only this T-channel model is available, plug in the thalamic relay cell model in your Layer 2/3 V1 pyramidal dendrite model, even if the parameters are wrong, at least you’ll get a degree of the complexity”. My experience is that the field of theoretical neuroscience is divided into two groups of people: those that fear high-dimensionality (usually people with maths/physics background) and those that fear the lack of complexity (usually more biology backgrounds). To some extent, I find it constructive to have those antagonistic approaches in the field, so it could make sense to have a simulator that would allow those people to interact easily, but that’s just my opinion. And anyway, the only way to discuss/refine/discard the Frankenstein models is to benefit from good implementations :wink:

Hi everyone (and thanks to @yzerlaut for the additional thoughts on this subject). I just opened a pull request to discuss adding a more structured syntax for equations:
I think having this would make @rth and @yzerlaut at least a little bit happy, right :slight_smile:?
I’d be very grateful for feedback, please let me know what you think.

1 Like

This looks fantastic !
To my opinion, reducing this (relatively complex) model to a set of gating_var in the implementation thanks to the templates is a great idea. It really makes sense that this gating_var object should be emphasized as the (or one of the) building block of classical biophysical modeling. The more I think about it, the more I think that it is a much better unit to work with than the “mechanism” ! Very elegant solution… :ok_hand:
And yes, this implementation also addresses all the limitations that I saw in the large Equation strings: it drastically improves readability, it minimizes the potential mistakes in setting parameters in the questions, etc…
Looking forward to play with this ! (I will follow your PR and work on my repo with this implementation).
So yeah, I am really happy ! :wink: Thanks a lot @mstimberg and everyone involved !


@mstimberg, WoW! I completely agree with @yzerlaut, it looks like an elegant and well thought solution for the problem. Nesting templates are perfect here, really well done :clap: What isn’t clear yet, how programming access to a gating_var will be organized in the model.

1 Like

Thanks for all the praise, but please don’t forget that the current code in the pull request is more of a proof-of-concept and it still needs quite a bit of work (error checking, handling corner cases, etc.).

Could you give some more detail on what you mean by this? Do you mean the reference by name to m & h in the line

ina = Equations('ina = gnabar*m**3*h*(ENa-v) : amp') + m + h


Oh, and if you have an opinion on the Equations vs. EquationTemplate question (see PR), please let me know!

@mstimberg, sorry for the lack of clarity.
Yes, you are right. I wonder, will this automatically separate gating variables, by channel and add some prefix/suffix to variable name (prefix is better, Neuron’s hoc suffix notation is very confusing for beginners).

For example can I use the same m and h names in ina and ical channels?

ina = Equations('ina = gnabar*m**3*h*(ENa-v) : amp') +\
               rate_expression=pos_sigmoid(midpoint=-38., scale=7.),
               forward_rate=exp_voltage_dep(magnitude=5., midpoint=-60, scale=18.),
               reverse_rate=neg_exp_voltage_dep(magnitude=36., midpoint=-60, scale=25.),
               tau_base=0.04*ms, tau_scale=10*ms) + \
               rate_expression=neg_sigmoid(midpoint=-65., scale=6.),
               forward_rate=exp_voltage_dep(magnitude=7., midpoint=-60., scale=11.),
               reverse_rate=neg_exp_voltage_dep(magnitude=10., midpoint=-60., scale=25.),
               tau_base=0.6*ms, tau_scale=100*ms)
ical = Equations('ical = gcalbar*m*h*(ECa-v) : amp') +\
             ) + \

Re Equations vs. EquationTemplate : I think, if you plane only name substitution, without comprehensive analysis, straightforward {...} notation is good. However, a separate template class will be useful for modularity.

Hi again. I thought about the prefix/suffix issue, but I think for now I’d rather not implement anything automatic for this, because it would also mean to change the way things are accessed later, would make the link between the equations and the names of the variables more indirect, etc. Your example would indeed not work at the moment, but with my latest changes in the pull request you’d only have to change e.g. name='m' to name='m_na' in a single place.

I actually think it is better to just go with Equations for now, I’ve laid out my reasoning in the pull request. Would be happy to hear your opinion on it! I’d suggest we move the discussion about this proposal/pull request to the pull request on github, otherwise it gets a bit confusing.