Model fitting with refractory period

I’m following the model fitting example of Tutorial: TraceFitter — Brian2modelfitting 0.3 documentation, but I’m using an adaptive exponential model instead of the given HH model. Now I want to introduce also refractory period, but I get an error. Here’s my code:

from brian2modelfitting import *
from brian2 import *
import pandas as pd
inp_trace = pd.read_csv('input_traces_hh.csv', index_col=0).to_numpy()
out_trace = pd.read_csv('output_traces_hh.csv', index_col=0).to_numpy()
area = 20000*umetre**2                                                                                                                                                                                                                    

model = '''dv/dt = (gL*(EL-v)+gL*DeltaT*exp((v-VT)/DeltaT) -w + I)/Cm : volt (unless refractory)                                                                                                                                                                                  
  dw/dt=(a*(v-EL)-w)/tauw : amp                                                                                                                                                                                                                                                   
  gL : siemens (constant)                                                                                                                                                                                                                                                         
  EL : volt (constant)                                                                                                                                                                                                                                                            
  DeltaT :  volt (constant)                                                                                                                                                                                                                                                       
  VT : volt (constant)                                                                                                                                                                                                                                                            
  a : siemens (constant)                                                                                                                                                                                                                                                          
  tauw : second (constant)                                                                                                                                                                                                                                                        
  wadd : ampere (constant)                                                                                                                                                                                                                                                        
  vre : volt (constant)                                                                                                                                                                                                                                                           
  Cm : farad (constant)                                                                                                                                                                                                                                                           
  tref : second (constant)                                                                                                                                                                                                                                                        

opt = NevergradOptimizer()
metric = MSEMetric()
fitter = TraceFitter(model=model,
                     input=inp_trace * amp,
                     param_init={'v': -65*mV},
                     reset='v=vre; w += wadd',
                     threshold='v > 0*mV',

res, error =,
                        gL=[2*psiemens, 300*nsiemens],
traces = fitter.generate_traces()

The code works if I remove “(unless refractory)”, that is, remove the refractory period, but otherwise it gives me a division by zero. Any ideas how to fix this issue? Is refractory period supported at all in brian2modelfitting library?

I also seem to get NaNs in the objective function values when I remove “(unless refractory)”, unless I make the length of the refractory period shorter than the simulation time step.

Hi @Tuoma . Using refractoriness in fitted models should work in principle. The error you are getting is due to a bug, though it seems to be a bug in Brian’s code optimization itself rather than in brian2modelfitting. Brian’s code generation simplifies expression by replacing expressions with fixed values, e.g. 1.0 / 2.0 will be replaced by 0.5. Apparently in this case it tries to replace 1.0 / 0.0 by a value, which leads to the division by zero error. I will have to have a closer look into why this is happening in this specific case, but in the mean time I think you could use the following work-around. Instead of (unless refractory), you use dv/dt = not_refractory*(gL*....), i.e. manually multiply the RHS of the equation by 0 or 1. This is basically what (unless refractory) does, even though in general it does a bit more (e.g. discard changes to the membrane potential by synapses) – but this should not matter here.

Thanks, that seems to work!

1 Like

For future reference, this has now been fixed in Brian’s current development version: Make evaluate_expr more robust by mstimberg · Pull Request #1389 · brian-team/brian2 · GitHub