Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Make StochKv3.mod NEURON 9.0a compatible #35

Merged
merged 7 commits into from
Jun 18, 2024
Merged

Make StochKv3.mod NEURON 9.0a compatible #35

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
7a55267
Make StochKv3.mod NEURON9 compatible
darshanmandge Feb 28, 2024
ff15792
black fix
darshanmandge Feb 28, 2024
6a40e72
Update README
darshanmandge Feb 28, 2024
e364cb7
limit NEURON to 8.0.2
anilbey Jun 18, 2024
33f714f
relax NEURON upper version until 8.2.4
anilbey Jun 18, 2024
2db583a
update precision of opt test
anilbey Jun 18, 2024
512981e
undo threshold_current ground truth update
anilbey Jun 18, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 9 additions & 8 deletions README.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,9 +12,9 @@ A universal workflow for creation, validation, and generalization of detailed ne
Introduction
---------

The biophysically detailed electrical neuron model (e-model) is one of the central tools in neuroscience. Here we present a demo of the single neuron e-model creation, validation, and generalization described in `A universal workflow for creation, validation, and generalization of detailed neuronal models (Reva, Rossert et al) <https://doi.org/10.1016/j.patter.2023.100855>`_.
The biophysically detailed electrical neuron model (e-model) is one of the central tools in neuroscience. Here we present a demo of the single neuron e-model creation, validation, and generalization described in `A universal workflow for creation, validation, and generalization of detailed neuronal models (Reva, Rossert et al) <https://doi.org/10.1016/j.patter.2023.100855>`_.

This demo is built on the example of L5PC of the SSCx in juvenile rat.
This demo is built on the example of L5PC of the SSCx in juvenile rat.

Structure of the codebase
~~~~~~~~~~~~~~~~~~~~~~~~~
Expand All @@ -34,7 +34,7 @@ When you using methods or code from this repository for your research, we ask yo

A universal workflow for creation, validation, and generalization of detailed neuronal models" (Reva, Rossert et al).

.. code-block::
.. code-block::

@article{reva2023universal,
title={A universal workflow for creation, validation, and generalization of detailed neuronal models},
Expand Down Expand Up @@ -86,7 +86,7 @@ The visualization of the bAP/EPSP validations can be found in `validation.ipynb

The morphologies for these validations are located in the `input/morphologies <https://github.com/BlueBrain/SSCxEModelExamples/tree/main/validation/input/morphologies>`_ folder.

To run bAP/EPSP validations use::
To run bAP/EPSP validations use::

python main.py att_conf.json

Expand All @@ -98,10 +98,10 @@ Somatic validations are located in the `somatic_validation <https://github.com/B

Note that this is the only step that does not use the ``requirements.txt`` in the main directory.

`somatic-val-requirements.txt <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/somatic_validation/somatic-val-requirements.txt>`_ needs to be installed and the mechanisms need to be compiled with the following command before running the notebooks::
`somatic-val-requirements.txt <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/somatic_validation/somatic-val-requirements.txt>`_ needs to be installed and the mechanisms need to be compiled with the following command before running the notebooks::

nrnivmodl mechanisms

nrnivmodl mechanisms

First, e-features for the validations have to be extracted from the chosen patch clamp protocol. To extract e-features use `feature-extraction.ipynb <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/somatic_validation/feature-extraction.ipynb>`_, the results of this extraction can be found in the `somatic_validation/L5TPC <https://github.com/BlueBrain/SSCxEModelExamples/tree/main/somatic_validation/L5TPC>`_ folder. To run and visualize results of the somatic validation run `somatic-validation.ipynb <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/somatic_validation/somatic-validation.ipynb>`_.

4. Generalization
Expand Down Expand Up @@ -136,6 +136,7 @@ Instead of manually compiling the mechanisms for each step of the pipeline, the
Requirements
------------

For the exact reproducibility of paper results, we recommend `NEURON>=7.8.0,<8.0.0` and repository code release version `1.0.1 <https://github.com/BlueBrain/SSCxEModelExamples/releases/tag/1.0.1>`. Users with other versions may encounter some small discrepancies in results.
The `requirements.txt <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/requirements.txt>`_ at the main directory should be used for all steps except for the somatic validations.
Install `somatic-val-requirements.txt <https://github.com/BlueBrain/SSCxEModelExamples/blob/main/somatic_validation/somatic-val-requirements.txt>`_ before running the somatic validation notebooks or tests.

Expand All @@ -147,7 +148,7 @@ This project/research was supported by funding to the Blue Brain Project, a rese
License
-------

This work is licensed under `Creative Commons (CC BY) 4.0 <https://creativecommons.org/licenses/by/4.0/>`_
This work is licensed under `Creative Commons (CC BY) 4.0 <https://creativecommons.org/licenses/by/4.0/>`_

For MOD files for which the original source is available on ModelDB, any specific licenses on mentioned on ModelDB, or the generic License of ModelDB apply.

Expand Down
87 changes: 59 additions & 28 deletions model_management/mm_run_minimal/emodels_dir/singlecell-optimization/mechanisms/StochKv3.mod
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
TITLE StochKv2.mod
TITLE StochKv3.mod

COMMENT
----------------------------------------------------------------
Expand Down Expand Up @@ -28,11 +28,11 @@ NEURON {
USEION k READ ek WRITE ik
RANGE N, eta, gk, gamma, deterministic, gkbar, ik
RANGE N0, N1, n0_n1, n1_n0
GLOBAL ninf, linf, ltau, ntau, an, bn, al, bl
GLOBAL P_an, P_bn, P_al, P_bl
RANGE ninf, linf, ltau, ntau, an, bn, al, bl
RANGE P_an, P_bn, P_al, P_bl
GLOBAL vmin, vmax
:BBCOREPOINTER rng
POINTER rng
BBCOREPOINTER rng
:POINTER rng
}

UNITS {
Expand Down Expand Up @@ -99,15 +99,22 @@ for comparison with Pr to decide whether to activate the synapse or not
ENDCOMMENT

VERBATIM
#ifndef NRN_VERSION_GTEQ_8_2_0
#include "nrnran123.h"
extern int cvode_active_;

#include <stdlib.h>
#include <stdio.h>
#include <math.h>

#ifndef CORENEURON_BUILD
double nrn_random_pick(void* r);
void* nrn_random_arg(int argpos);
#endif
#define RANDCAST
#else
#define RANDCAST (Rand*)
#endif

ENDVERBATIM
: ----------------------------------------------------------------
Expand All @@ -117,6 +124,10 @@ INITIAL {
if (cvode_active_ && !deterministic) {
hoc_execerror("StochKv2 with deterministic=0", "cannot be used with cvode");
}

if( usingR123 ) {
nrnran123_setseq((nrnran123_State*)_p_rng, 0, 0);
}
ENDVERBATIM

eta = (gkbar / gamma) * (10000)
Expand Down Expand Up @@ -248,9 +259,12 @@ PROCEDURE trates(v (mV)) {
FUNCTION strap(x) {
if (x < 0) {
strap = 0
VERBATIM
fprintf (stderr,"skv.mod:strap: negative state");
ENDVERBATIM
: This function gets executed in the DERIVATIVE and BREAKPOINT blocks
: which should get vectorized and fprintf doesn't allow this
: Also Intel classic compilers generate invalid code when #pragma omp simd is used
:VERBATIM
: fprintf (stderr,"skv.mod:strap: negative state");
:ENDVERBATIM
} else {
strap = x
}
Expand All @@ -273,11 +287,12 @@ PROCEDURE setRNG() {
VERBATIM
// For compatibility, allow for either MCellRan4 or Random123. Distinguish by the arg types
// Object => MCellRan4, seeds (double) => Random123
#if !NRNBBCORE
#ifndef CORENEURON_BUILD
usingR123 = 0;
if( ifarg(1) && hoc_is_double_arg(1) ) {
nrnran123_State** pv = (nrnran123_State**)(&_p_rng);
uint32_t a2 = 0;
uint32_t a3 = 0;

if (*pv) {
nrnran123_deletestream(*pv);
Expand All @@ -286,7 +301,10 @@ VERBATIM
if (ifarg(2)) {
a2 = (uint32_t)*getarg(2);
}
*pv = nrnran123_newstream((uint32_t)*getarg(1), a2);
if (ifarg(3)) {
a3 = (uint32_t)*getarg(3);
}
*pv = nrnran123_newstream3((uint32_t)*getarg(1), a2, a3);
usingR123 = 1;
} else if( ifarg(1) ) {
void** pv = (void**)(&_p_rng);
Expand All @@ -302,46 +320,56 @@ ENDVERBATIM
FUNCTION urand() {

VERBATIM
double value;
double value = 0.0;
if( usingR123 ) {
value = nrnran123_dblpick((nrnran123_State*)_p_rng);
} else if (_p_rng) {
value = nrn_random_pick(_p_rng);
#ifndef CORENEURON_BUILD
value = nrn_random_pick(RANDCAST _p_rng);
#endif
} else {
value = 0.5;
// see BBPBGLIB-972
value = 0.0;
}
_lurand = value;
Comment on lines -311 to 333
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"The changes to mod files like StochKv3 are all technical. There is no change in intrinsic behaviour / results "

Just to remain "true" to the above quote: this change to value variable is the one something that we were a bit unsure about. In the JIRA linked issue we explained why this was necessary and how this was the issue in load-balancing / ExperimentalMechComplex usage with NEURON.

We have ran simulations and results are same and now this change is in production. But if you really want to keep exact old code, you can set the value back to 0.5.

Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just to add some more information if someone will revisit this in the future:

FUNCTION urand() {VERBATIM
    double value;
    if( usingR123 ) {
        value = nrnran123_dblpick((nrnran123_State*)_p_rng);
    } else if (_p_rng) {
#ifndef CORENEURON_BUILD
        value = nrn_random_pick(RANDCAST _p_rng);
#endif
    } else {
        
        value = 0.5;
    }
    _lurand = value;
ENDVERBATIM
}

In the above urand() function we are picking a random number. In all practical usage, we would have an RNG stream setup either using Random123 or MCellRan4. The last else block returns a value 0.5 (or 0) when a user hasn't set up the RNG variable i.e. StochKV is not properly instantiated. This won't be the case if someone is really using StochKV.

So I would say let's ignore my initial comment about this change. Using the new MOD file is perfectly OK.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the detailed explanation @pramodk!

ENDVERBATIM
}

VERBATIM
/*
static void bbcore_write(double* x, int* d, int* xx, int* offset, _threadargsproto_) {
if (d) {
uint32_t* di = ((uint32_t*)d) + *offset;
// temporary just enough to see how much space is being used
if (!_p_rng) {
di[0] = 0; di[1] = 0;
di[0] = 0; di[1] = 0, di[2] = 0;
}else{
nrnran123_State** pv = (nrnran123_State**)(&_p_rng);
nrnran123_getids(*pv, di, di+1);
nrnran123_getids3(*pv, di, di+1, di+2);
// write stream sequence
char which;
nrnran123_getseq(*pv, di+3, &which);
di[4] = (int)which;
}
//printf("StochKv2.mod %p: bbcore_write offset=%d %d %d\n", _p, *offset, d?di[0]:-1, d?di[1]:-1);
//printf("StochKv3.mod %p: bbcore_write offset=%d %d %d\n", _p, *offset, d?di[0]:-1, d?di[1]:-1);
}
*offset += 2;
*offset += 5;
}
static void bbcore_read(double* x, int* d, int* xx, int* offset, _threadargsproto_) {
assert(!_p_rng);
uint32_t* di = ((uint32_t*)d) + *offset;
if (di[0] != 0 || di[1] != 0)
if (di[0] != 0 || di[1] != 0|| di[2] != 0)
{
nrnran123_State** pv = (nrnran123_State**)(&_p_rng);
*pv = nrnran123_newstream(di[0], di[1]);
#if !NRNBBCORE
if(*pv) {
nrnran123_deletestream(*pv);
}
#endif
*pv = nrnran123_newstream3(di[0], di[1], di[2]);
nrnran123_setseq(*pv, di[3], (char)di[4]);
}
//printf("StochKv2.mod %p: bbcore_read offset=%d %d %d\n", _p, *offset, di[0], di[1]);
*offset += 2;
//printf("StochKv3.mod %p: bbcore_read offset=%d %d %d\n", _p, *offset, di[0], di[1]);
*offset += 5;
}
*/
ENDVERBATIM

: Returns random numbers drawn from a binomial distribution
Expand Down Expand Up @@ -471,13 +499,16 @@ VERBATIM
FUNCTION bbsavestate() {
bbsavestate = 0
VERBATIM
#ifdef ENABLE_SAVE_STATE
#ifndef CORENEURON_BUILD
// TODO: since N0,N1 are no longer state variables, they will need to be written using this callback
// provided that it is the version that supports multivalue writing
/* first arg is direction (-1 get info, 0 save, 1 restore), second is value*/
double *xdir, *xval, *hoc_pgetarg();
double *xdir, *xval;
#ifndef NRN_VERSION_GTEQ_8_2_0
double *hoc_pgetarg();
long nrn_get_random_sequence(void* r);
void nrn_set_random_sequence(void* r, int val);
#endif
xdir = hoc_pgetarg(1);
xval = hoc_pgetarg(2);
if (_p_rng) {
Expand All @@ -498,13 +529,13 @@ VERBATIM
xval[0] = (double) seq;
xval[1] = (double) which;
} else {
xval[0] = (double)nrn_get_random_sequence(_p_rng);
xval[0] = (double)nrn_get_random_sequence(RANDCAST _p_rng);
}
} else{
if( usingR123 ) {
nrnran123_setseq( (nrnran123_State*)_p_rng, (uint32_t)xval[0], (char)xval[1] );
} else {
nrn_set_random_sequence(_p_rng, (long)(xval[0]));
nrn_set_random_sequence(RANDCAST _p_rng, (long)(xval[0]));
}
}
}
Expand Down
1 change: 1 addition & 0 deletions model_management/mm_run_minimal/emodels_dir/singlecell-optimization/setup/__init__.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
To view a copy of this license, visit https://creativecommons.org/licenses/by/4.0/legalcode or send a letter to Creative Commons, 171
Second Street, Suite 300, San Francisco, California, 94105, USA.
"""

from . import evaluator
from . import template

Expand Down
21 changes: 10 additions & 11 deletions model_management/mm_run_minimal/emodels_dir/singlecell-optimization/setup/evaluator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -405,14 +405,16 @@ def define_protocols(
"type" in protocol_definition
and protocol_definition["type"] == "RatSSCxThresholdDetectionProtocol"
):
protocols_dict[
"ThresholdDetection"
] = protocols.RatSSCxThresholdDetectionProtocol(
"IDRest",
step_protocol_template=read_step_protocol(
"Threshold", protocol_definition["step_template"], recordings
),
prefix=prefix,
protocols_dict["ThresholdDetection"] = (
protocols.RatSSCxThresholdDetectionProtocol(
"IDRest",
step_protocol_template=read_step_protocol(
"Threshold",
protocol_definition["step_template"],
recordings,
),
prefix=prefix,
)
)
else:
stimuli = []
Expand Down Expand Up @@ -488,7 +490,6 @@ def calculate_score(self, responses, trace_check=False):


class eFELFeatureExtra(eFELFeature):

"""eFEL feature extra"""

SERIALIZED_FIELDS = (
Expand Down Expand Up @@ -647,7 +648,6 @@ def calculate_score(self, responses, trace_check=False):


class SingletonWeightObjective(EFeatureObjective):

"""Single EPhys feature"""

def __init__(self, name, feature, weight):
Expand Down Expand Up @@ -775,7 +775,6 @@ def define_fitness_calculator(


class MultiEvaluator(bpopt.evaluators.Evaluator):

"""Multiple cell evaluator"""

def __init__(
Expand Down
13 changes: 5 additions & 8 deletions model_management/mm_run_minimal/emodels_dir/singlecell-optimization/setup/protocols.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,6 @@


class RatSSCxMainProtocol(ephys.protocols.Protocol):

"""Main protocol to fit RatSSCx neuron ephys parameters.

Pseudo code:
Expand Down Expand Up @@ -193,7 +192,6 @@ def run(self, cell_model, param_values, sim=None, isolate=None):


class RatSSCxRinHoldcurrentProtocol(ephys.protocols.Protocol):

"""IDRest protocol to fit RatSSCx neuron ephys parameters"""

def __init__(
Expand Down Expand Up @@ -385,7 +383,6 @@ def voltage_base(self, current, cell_model, param_values, sim=None, short=False)


class RatSSCxThresholdDetectionProtocol(ephys.protocols.Protocol):

"""IDRest protocol to fit RatSSCx neuron ephys parameters"""

def __init__(
Expand Down Expand Up @@ -641,7 +638,6 @@ def search_spike_threshold(


class StepProtocol(ephys.protocols.SweepProtocol):

"""Protocol consisting of step and holding current"""

def __init__(
Expand All @@ -665,9 +661,11 @@ def __init__(

super(StepProtocol, self).__init__(
name,
stimuli=step_stimuli + [holding_stimulus]
if holding_stimulus is not None
else step_stimuli,
stimuli=(
step_stimuli + [holding_stimulus]
if holding_stimulus is not None
else step_stimuli
),
recordings=recordings,
cvode_active=cvode_active,
)
Expand Down Expand Up @@ -752,7 +750,6 @@ def step_amplitude(self):


class StepThresholdProtocol(StepProtocol):

"""Step protocol based on threshold"""

def __init__(
Expand Down
12 changes: 6 additions & 6 deletions model_management/mm_run_minimal/emodels_dir/singlecell-optimization/setup/template.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -122,12 +122,12 @@ def define_parameters(params_filename, stage=None, past_params=None):
dist_param_names = definition["parameters"]
else:
dist_param_names = None
distributions[
distribution
] = ephys.parameterscalers.NrnSegmentSomaDistanceScaler(
name=distribution,
distribution=definition["fun"],
dist_param_names=dist_param_names,
distributions[distribution] = (
ephys.parameterscalers.NrnSegmentSomaDistanceScaler(
name=distribution,
distribution=definition["fun"],
dist_param_names=dist_param_names,
)
)

params_definitions = definitions["parameters"]
Expand Down
Loading
Loading