# Reproducibility issue between brian2.5 and 2.7

**URL:** <https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240>\
**Category:** Support\
**Tags:** synapses, reproducibility, custom-event\
**Created:** [19 June 2024 11:19 UTC](https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240 "2024-06-19T11:19:11Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![mreynes](https://avatars.discourse-cdn.com/v4/letter/m/3be4f8/32.png) [@mreynes](https://brian.discourse.group/u/mreynes)\
**Post date:** [19 June 2024 11:19 UTC](https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240/1 "2024-06-19T11:19:11Z")

</div>

Description of problem  
I have a model running perfectly fine on brian 2.5 (python 3.8) and I would like it to run on brian & python latests versions. However, the model does not reproduce at all the results from version 2.5. The miniEPSPs (mini excitatory post synaptic potentials) do not work. Those are defined as custom events, appearing following a poisson process. The rate of arrival is defined in my axosomatic compartment, as a function dependent on the last spike (here, pyramidal spike but this is not important). Depending on that rate, the custom\_event (defining the time of a poisson process) is triggered and this triggers in return changes in the synapses (releasing small potentials).

Minimal code to reproduce problem

The custom\_event (defined when instantiating the axosomatic compartment of our model):

```auto
PY_soma = NeuronGroup(N_PY,Soma_eqs,method='rk4',threshold='v>40*mV',refractory=3*ms,events={'custom_poisson_PY':'rand()<mean_rate_PY*dt','custom_poisson_IN':'rand()<mean_rate_IN*dt'})

```

Equation for the axosomatic compartment of the model where the mean\_rate is defined:

```auto
Soma_eqs = '''

    v = -68*mV*init_timestep + truev*(1-init_timestep) : volt
    init_timestep = init_timedarray(t) : 1
    
    truev = (vgap + (R * s_soma) * g2) / (1 + (R * s_soma) * g1) : volt    
                            
    g1 = g_na *(m_na **3) * h_na * Phi_m + g_k * n_k * Tad + g_nap * m_nap : siemens * meter**-2
    g2 = g_na *(m_na **3) * h_na * Phi_m * E_na + g_k * n_k * Tad * E_k + g_nap * m_nap * E_na + 6.74172*uA*cm**-2 + Iext : amp*meter**-2
   
    I_na = Phi_m * g_na * (m_na **3) * h_na * (v - E_na) : amp * meter**-2
        dm_na/dt = -(m_na - m_nainf) / tau_m_na : 1
        dh_na/dt = -(h_na - h_nainf) / tau_h_na : 1
        alpham_na = (0.182/ms * (v + 25*mV)/mV / (1 - exp(-(v + 25*mV)/9/mV))) * int(abs((v-10*mV)/(-35*mV)) > 1e-6) + (0.182/ms * 9) * int(abs((v-10*mV)/(-35*mV)) < 1e-6) : Hz
        betam_na = (-0.124/ms * (v + 25*mV)/mV / (1 - exp((v + 25*mV)/9/mV))) * int(abs((-v+10*mV)/(35*mV)) > 1e-6) + (0.124/ms * 9) * int(abs((-v+10*mV)/(35*mV)) < 1e-6) : Hz
        alphah_na = (0.024/ms * (v + 40*mV)/mV / (1 - exp(-(v + 40*mV)/5/mV))) * int(abs((v-10*mV)/(-50*mV)) > 1e-6) + (0.024/ms * 5) * int(abs((v-10*mV)/(-50*mV)) < 1e-6) : Hz
        betah_na = (-0.0091/ms * (v + 65*mV)/mV / (1 - exp((v + 65*mV)/5/mV))) * int(abs((-v+10*mV)/(75*mV)) > 1e-6) + (0.0091/ms * 5) * int(abs((-v+10*mV)/(75*mV)) < 1e-6) : Hz
        h_nainf = 1 / (1 + exp((v + 55*mV)/6.2/mV)) : 1
        tau_h_na = (1 / (alphah_na + betah_na)) / Qt : second   
        m_nainf = alpham_na / (alpham_na + betam_na) : 1
        tau_m_na = (1 / (alpham_na + betam_na)) / Qt : second
        
    I_nap = g_nap * m_nap * (v - E_na) : amp * meter**-2 
        dm_nap/dt = -(m_nap - m_napinf)/tau_m_nap : 1
        m_napinf = 0.02 / (1 + exp(-(v + 42*mV)/5/mV)) : 1 
    
    I_k = Tad * g_k * n_k * (v - E_k) : amp * meter**-2
        dn_k/dt = -(n_k - n_kinf) / tau_n_k : 1
        alphan_k = (0.02/mV) * (v - 25*mV) / (1 - exp(-(v - 25*mV)/(9*mV))) : 1
        betan_k = (-0.002/mV) * (v - 25*mV) / (1 - exp((v - 25*mV)/(9*mV))) : 1
        n_kinf = alphan_k / (alphan_k + betan_k) : 1
        tau_n_k = (1*msecond / (alphan_k + betan_k)) / Qt : second
        
    Iext : amp * meter**-2
        
    vgap : volt
     
    g_na : siemens * meter**-2
    g_k : siemens * meter**-2
    g_nap : siemens * meter**-2
    
    x : meter
    
    t_last_spike_PY : second
    
    t_last_spike_IN : second
    
    mean_rate_PY = (log((t-t_last_spike_PY + 50*ms)/(50*ms))/400*kHz)*int((t-t_last_spike_PY)>70*ms) : Hz
    
    mean_rate_IN = (log((t-t_last_spike_PY + 50*ms)/(50*ms))/400*kHz)*int((t-t_last_spike_PY)>70*ms) : Hz
    
    mean_rate_GABAA = (log((t-t_last_spike_IN + 50*ms)/(50*ms))/400*kHz)*int((t-t_last_spike_IN)>70*ms) : Hz

    '''

```

The synapse from where the issue is coming from:

```auto
def syn_ampa(source,target,syntype,surfacetarget,connection_pattern,g_syn,A_mEPSP,poissontype,connection_pattern_mini,synapses_received):
    eq_syn_AMPA = '''_post = (D * g_syn * (1/surfacetarget) * (W/mM) * (v_post - E_syn)) / N_incoming : amp * meter**-2 (summed)
        dW/dt = alpha * (1*mM - W) * T - beta * W : mM (clock-driven)
        T = A_syn * int((t_last_spike + tmax - t) > 0*ms)* int((t - t_last_spike) > 0*ms) : mM
        D : 1
        q1 : second
        g_syn : siemens
        E_syn = 0*mV : volt
        t_last_spike : second
        alpha = 0.94*ms **-1*mM** -1 : second **-1*mM** -1
        beta = 0.18*kHz : Hz
        U = 0.073 : 1
        tau = 700*ms : second
        tmax = 0.3*ms : second
        A_syn = 0.5*mM : mM
        connection_pattern_mini : 1
        surfacetarget : meter**2
        synapses_received : 1
    '''
    eq_syn_AMPA2 = '''_post = (D * A_mEPSP * (1/surfacetarget) * (W_Poisson/mM) * (v_post - E_syn)) / N_incoming : amp * meter**-2 (summed)
        dW_Poisson/dt = alpha * (1*mM - W_Poisson) * T_Poisson - beta * W_Poisson : mM (clock-driven)
        T_Poisson = (A_syn * int((t_last_spike_Poisson_PY + tmax - t) > 0*ms)* int((t - t_last_spike_Poisson_PY) > 0*ms))*int(connection_pattern_mini == 0) + (A_syn * int((t_last_spike_Poisson_IN + tmax - t) > 0*ms)* int((t - t_last_spike_Poisson_IN) > 0*ms))*int(connection_pattern_mini == 1) : mM
        t_last_spike_Poisson_PY : second
        t_last_spike_Poisson_IN : second
        A_mEPSP : siemens
    '''
    pre_code_ampa = '''q1 = ((t - t_last_spike) - tmax) 
        D = 1 - (1 - D*(1 - U)) * exp(-q1/tau)*int((-q1/tau)>-10)*int((-q1/tau)<10)
        t_last_spike = t
        t_last_spike_PY_pre = t
    '''
    poisson_PY = '''
        q1 = ((t - t_last_spike) - tmax) 
        t_last_spike_Poisson_PY = t
        D = (1 - (1 - D*(1 - 0)) * exp(-q1/tau)*int((-q1/tau)>-10)*int((-q1/tau)<10))*int(connection_pattern_mini == 0) + D*int(connection_pattern_mini == 1)
        W0 = W
    '''
    poisson_IN = '''
        q1 = ((t - t_last_spike) - tmax) 
        t_last_spike_Poisson_IN = t
        D = (1 - (1 - D*(1 - 0)) * exp(-q1/tau)*int((-q1/tau)>-10)*int((-q1/tau)<10))*int(connection_pattern_mini == 1) + D*int(connection_pattern_mini == 0)
        W0 = W
    '''
    S=Synapses(source,target,model=(syntype+eq_syn_AMPA)+(poissontype+eq_syn_AMPA2),method='rk4',on_pre={'pre': pre_code_ampa,'mini_PY': poisson_PY,'mini_IN': poisson_IN},on_event={'mini_PY':'custom_poisson_PY','mini_IN':'custom_poisson_IN'})
    if connection_pattern=='':
        S.connect()
    else :
        S.connect(condition=connection_pattern, skip_if_invalid=True)
    S.surfacetarget = surfacetarget
    S.g_syn = g_syn
    S.A_mEPSP = A_mEPSP
    S.connection_pattern_mini = connection_pattern_mini
    S.synapses_received = synapses_received
    # print("AMPA synapses from "+str(source)+" to "+str(target))
    # print(S.N_outgoing_pre)
    # print(S.N_incoming_post)
    return S

```

What you have aready tried  
I have considered it being a ran() error, but the custom\_event actually happen correctly (can be easily tracked using a monitor like `M0=EventMonitor(PY_soma,'custom_poisson_PY',record=True)`). The issue is coming from the ‘integration’ between the custom\_event and the synapse itself: the variables that should be updated when this event happen are simply not (all of them, not just a single one). I have gone through the release notes of both 2.6 and 2.7, and can’t find anything that changed that would impact the synapse.

Expected output (if relevant)  
Update of variables in ‘poisson\_PY’ of the synapse, as the event ‘custom\_poisson\_PY’ happen.

Actual output (if relevant)  
No update of variables

Thank you in advance for your help.

---

<div class="post-metadata">

**Author:** ![mstimberg](https://yyz2.discourse-cdn.com/free1/user_avatar/brian.discourse.group/mstimberg/32/11_2.png) [@mstimberg](https://brian.discourse.group/u/mstimberg)\
**Post date:** [19 June 2024 11:56 UTC](https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240/2 "2024-06-19T11:56:55Z")

</div>

Hi @mreynes. Many thanks for this report, there might have been some inadvertent change in Brian that cause this. When I recently answered another question here on the forum, I stated that " In fact, the `on_event` dictionary should only refer to `pre` and `post`, i.e. you can only set one event type for all `pre` and one event type for all `post` pathways." ([How to activate synaptic pathways on custom events only? - #2 by mstimberg](https://brian.discourse.group/t/how-to-activate-synaptic-pathways-on-custom-events-only/1213/2)). I made this statement after looking at the code, but after having another look at the documentation ([Custom events — Brian 2 2.7.0 documentation](https://brian2.readthedocs.io/en/stable/advanced/custom_events.html#triggering-synaptic-pathways)), it actually seems as if providing custom pathway names should be supported. Since your code uses this feature, and you had it work with Brian 2.5, so it seems to have been broken in an update. I’ll have a closer look soon and keep you posted!

---

<div class="post-metadata">

**Author:** ![mstimberg](https://yyz2.discourse-cdn.com/free1/user_avatar/brian.discourse.group/mstimberg/32/11_2.png) [@mstimberg](https://brian.discourse.group/u/mstimberg)\
**Post date:** [22 June 2024 16:09 UTC](https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240/3 "2024-06-22T16:09:31Z")

</div>

Hi again @mreynes. So it turned out that this was indeed a bug/regression in Brian that slipped in unnoticed! It should be fixed in the latest development version, and I will do a new “proper” release soon that will include the fix. I’d be very grateful if you could test it with the development version, which you can install from the pypi test server like this:

```auto
pip install -i https://test.pypi.org/simple/ Brian2==2.7.0.post5

```

This requires that the dependencies are already installed. If you do the install from scratch, you’d have to add

```auto
--extra-index-url https://pypi.org/simple

```

to the above command. Please let us know whether everything works as before with this version.

---

<div class="post-metadata">

**Author:** ![mstimberg](https://yyz2.discourse-cdn.com/free1/user_avatar/brian.discourse.group/mstimberg/32/11_2.png) [@mstimberg](https://brian.discourse.group/u/mstimberg)\
**Post date:** [2 July 2024 14:53 UTC](https://brian.discourse.group/t/reproducibility-issue-between-brian2-5-and-2-7/1240/4 "2024-07-02T14:53:11Z")

</div>

This fix has been included in Brian 2.7.1: [Release notes — Brian 2 2.7.1 documentation](https://brian2.readthedocs.io/en/stable/introduction/release_notes.html)