# Convert code from Brian1 to Brian2

**URL:** <https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704>\
**Category:** Support\
**Created:** [30 May 2022 14:58 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704 "2022-05-30T14:58:28Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![serge25](https://avatars.discourse-cdn.com/v4/letter/s/4af34b/32.png) [@serge25](https://brian.discourse.group/u/serge25)\
**Post date:** [30 May 2022 14:58 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/1 "2022-05-30T14:58:28Z")

</div>

# Description of problem

I’m trying to convert Brian1 program ([https://senselab.med.yale.edu/modeldb/ShowModel?model=189346#tabs-1](https://senselab.med.yale.edu/modeldb/ShowModel?model=189346#tabs-1)) to Brian2.  
I have error "‘The identifier “stim” could not be resolved.’ and I don’t know how to solve it.

# Minimal code to reproduce problem

```auto
from brian2 import *
import sys, subprocess, os
import numpy as np
import matplotlib
from matplotlib import lines
import pylab as pl

defaultclock.dt = 0.05 *ms

# Number of neurons
M = 3000
Ne = int(M)
Ni = int(0.25*M)
p = 0. # neurons are not connected.

# ext_max: maximal external input
# EK: potassium reversal potential
ext_max = 8 *kHz
EK = [-90*mV, -75*mV, -60*mV]

for vk_n in EK:

    # Synaptic reversal potential
    ve = 25 *mV
    vi = -70 *mV
    
    # Excitatory ML parameters
    Ce = 5 *nF 
    gke = 22 *usiemens 
    gnae= 20 *usiemens 
    gle = 2 *usiemens

    v1e = -7 *mV
    v2e = 15 *mV
    v3e = -8 *mV 
    v4e = 15 *mV
    v3te = 0 *mV
    v4te = 15 *mV
    
    vke = vk_n 
    vnae= 60 *mV
    vle = -60 *mV
    phie = 0.55
    
    # Inhibitory ML parameters
    Ci = 1.0 *nF
    gki = 7.0 *usiemens
    gnai= 16.0 *usiemens
    gli = 7.0 *usiemens 

    v1i = -7.35 *mV 
    v2i = 23.0 *mV 
    v3i = -12.0 *mV
    v4i = 12.0 *mV 
    v3ti = -15.0 *mV
    v4ti = 10.0 *mV 

    vki = vk_n 
    vnai= 60 *mV
    vli = -60 *mV
    phii = 1.0

    # simulation time
    simt = 2000 *ms    

    # synaptic decay time constants
    tauext = 3.0 *ms       
    taue = 3.0 *ms 
    taui = 7.0 *ms 

    # synaptic coupling strength (neurons are not connected.)
    Jee = 0.0 *nsiemens  
    Jie = 0.0 *nsiemens  
    Jei = 0.0 *nsiemens  
    Jii = 0.0 *nsiemens  
    # synaptic strength of external inputs
    win = 200 *nsiemens
    
    # spike threshold, refractory period
    v_the = 15 *mV
    v_thi = 10 *mV
    rfe = 2 *ms
    rfi = 1 *ms

    # External inputs to exc and inh neurons.
    # Given a time vector tvec = [t0, t1, ..., tn] and a rate vector rvec = [r0, r1, ..., rn],
    # the function stim(t) generates a piecewise linear function that connects (t0,r0), (t1,r1), ..., (tn,rn).
    tvec = [0.*ms, simt] 
    rvec = [0.*Hz, ext_max]

    def heavi(x):
        if x > 0:
            return 1
        elif x <=0:
            return 0
        
    def stim(t):
        r = 0.*Hz
        for j in range(len(tvec)-1):
            r = r + (rvec[j] + (rvec[j+1]-rvec[j])/(tvec[j+1]-tvec[j])*(t-tvec[j]))*(heavi(t-tvec[j]) - heavi(t-tvec[j+1]))
        return r

    ratee = 'stim(t)'
    ratei = 'stim(t)'

    # neuron models
    
    # excitatory ML neurons
    eqse = Equations('''

    dv/dt = ( - gnae*minf*(v - vnae) - gke*w*(v - vke) - gle*(v - vle) \
             - gext*(v - ve) - ge*(v - ve) - gi*(v - vi) )/Ce : volt

    dw/dt = phie*(winf - w)/tauw : 1

    dgext/dt = -gext/tauext : siemens

    dge/dt = -ge/taue : siemens

    dgi/dt = -gi/taui : siemens

    minf = .5*(1 + tanh((v - v1e)/v2e)) : 1

    winf = .5*(1 + tanh((v - v3e)/v4e)) : 1

    tauw = 1/cosh((v - v3te)/(2*v4te))*ms : second
        
    ''')

    # inhibitory ML neurons
    eqsi = Equations('''

    dv/dt = ( - gnai*minf*(v - vnai) - gki*w*(v - vki) - gli*(v - vli) \
             - gext*(v - ve) - ge*(v - ve) - gi*(v - vi))/Ci : volt

    dw/dt = phii*(winf - w)/tauw : 1

    dgext/dt = -gext/tauext : siemens 

    dge/dt = -ge/taue : siemens

    dgi/dt = -gi/taui : siemens

    minf = .5*(1 + tanh((v - v1i)/v2i)) : 1

    winf = .5*(1 + tanh((v - v3i)/v4i)) : 1

    tauw = 1/cosh((v - v3ti)/(2*v4ti))*ms : second

    ''')

    Pe = NeuronGroup(Ne, model=eqse,
         threshold='v>v_the',
         refractory = rfe)

    Pi = NeuronGroup(Ni, model=eqsi,
         threshold='v>v_thi',
         refractory = rfi)

    # Recurrent connections
    Cee = Synapses(Pe, Pe, on_pre='ge += Jee')
    Cee.connect(p=p)
    Cie = Synapses(Pe, Pi, on_pre='ge += Jie')
    Cie.connect(p=p)
    Cei = Synapses(Pi, Pe, on_pre='ge += Jei')
    Cei.connect(p=p)
    Cii = Synapses(Pi, Pi, on_pre='ge += Jii')
    Cii.connect(p=p)
    

    # External Poisson spikes
    inpute = PoissonGroup(Ne, rates = ratee)
    inputi = PoissonGroup(Ni, rates = ratei)

    # Connect the input to neurons.
    input_co1 = Synapses(inpute, Pe, on_pre='gext=win')    
    input_co2 = Synapses(inputi, Pi, on_pre='gext=win')
    
 # input_co1 = IdentityConnection(inpute, Pe, 'gext', weight=win)    
 # input_co2 = IdentityConnection(inputi, Pi, 'gext', weight=win)
        
    # Initialization
    Pe.v = -40*mV 
    Pe.w = 0.1 
    Pi.v = -40*mV 
    Pi.w = 0.1

    # Record population activity
    Me = StateMonitor(Pe, ['v'], record=True)
    Mi = StateMonitor(Pi, ['v'], record=True)

    Re = PopulationRateMonitor(Pe)
    Ri = PopulationRateMonitor(Pi)

    # Run simulations
    net = Network(Pe,Pi,inpute,inputi,input_co1,input_co2,Cee,Cei,Cie,Cii,Me,Mi,Re,Ri)

    net.run(simt)

    # Display results.
    
    #===== inhibitory rate =====#
    pl.figure(1,figsize=(10,8))
    axfig2 = pl.subplot(111)
    plot(np.arange(0,8,8./len(Ri.times)),Ri.rate,label=r'$\mathregular{E_{K}: }$'+str(vk_n))

    # remove boundary
    axfig2.spines['top'].set_visible(False)
    axfig2.spines['right'].set_visible(False)
    axfig2.get_xaxis().tick_bottom()
    axfig2.get_yaxis().tick_left()

    # axis
    pl.xlim([0.,8.])
    pl.ylim([0,550])
    pl.xticks(np.arange(0,9,2))
    pl.yticks([100,200,300,400,500])
    pl.xlabel('input rate (kHz)',fontsize=30)
    pl.ylabel('inh rate (Hz)',fontsize=30)

    # legend
    pl.legend(loc='best',fontsize=24,frameon=False)
    leg = pl.gca().get_legend()
    llines=leg.get_lines()
    pl.setp(llines,linewidth=3)

    #===== excitatory rate =====#
    pl.figure(2,figsize=(10,8))
    axfig2 = pl.subplot(111)
    plot(np.arange(0,8,8./len(Re.times)),Re.rate,label=r'$\mathregular{E_{K}: }$'+str(vk_n))

    # remove boundary
    axfig2.spines['top'].set_visible(False)
    axfig2.spines['right'].set_visible(False)
    axfig2.get_xaxis().tick_bottom()
    axfig2.get_yaxis().tick_left()

    # axis
    pl.xlim([0.,8.])
    pl.ylim([0,550])
    pl.xticks(np.arange(0,9,2))
    pl.yticks([100,200,300,400,500])
    pl.xlabel('input rate (kHz)',fontsize=30)
    pl.ylabel('exc rate (Hz)',fontsize=30)

    # legend
    pl.legend(loc='best',fontsize=24,frameon=False)
    leg = pl.gca().get_legend()
    llines=leg.get_lines()
    pl.setp(llines,linewidth=3)

    matplotlib.rcParams.update({'font.size':24})

    
    reinit_default_clock(t=0*ms)
    clear(True)

pl.show()

```

# What you have aready tried

# Expected output (if relevant)

# Actual output (if relevant)

# Full traceback of error (if relevant)

---

<div class="post-metadata">

**Author:** ![adam](https://yyz2.discourse-cdn.com/free1/user_avatar/brian.discourse.group/adam/32/197_2.png) [@adam](https://brian.discourse.group/u/adam)\
**Post date:** [30 May 2022 17:48 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/2 "2022-05-30T17:48:50Z")

</div>

Hi serge,  
I’m not sure, but the documentation on user-defined functions talks about preparing functions for using in brian  
[https://brian2.readthedocs.io/en/stable/advanced/functions.html#user-provided-functions](https://brian2.readthedocs.io/en/stable/advanced/functions.html#user-provided-functions)  
so you could try the `@checkunits(...)` or `stim = Function(stim, ...)` approach.

---

<div class="post-metadata">

**Author:** ![serge25](https://avatars.discourse-cdn.com/v4/letter/s/4af34b/32.png) [@serge25](https://brian.discourse.group/u/serge25)\
**Post date:** [31 May 2022 18:12 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/3 "2022-05-31T18:12:08Z")

</div>

@adam Thank you. But at this code user function use time of simulation `def stim(t):`. Is it right using of time in simulator for user function?

---

<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:** [1 June 2022 10:03 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/4 "2022-06-01T10:03:02Z")

</div>

Hi @serge25. @adam’s answer is correct, defining your own functions for use with Brian 2 is slightly more complex compared to Brian 1 due to the code generation approach. Referring to the time `t` as an argument is ok.  
That said, using Python functions was often only necessary in the first place due to Brian 1’s limited expressiveness – with Brian 2 you can express more things directly in the equations and therefore might not have to resort to user-defined functions. From your code:

```python
    # External inputs to exc and inh neurons.
    # Given a time vector tvec = [t0, t1, ..., tn] and a rate vector rvec = [r0, r1, ..., rn],
    # the function stim(t) generates a piecewise linear function that connects (t0,r0), (t1,r1), ..., (tn,rn).
    tvec = [0.*ms, simt] 
    rvec = [0.*Hz, ext_max]

```

Here, the function `stim(t)` is quite general, but with the given `tvec` and `rvec` all it does is to linearly go from 0Hz to `ext_max` over the time course of the simulation. If that is all you need, you can directly write:

```python
ratee = '(t/simt)*ext_max'
ratei = '(t/simt)*ext_max'

```

which will do the same thing without referring to any user-defined function (and therefore also automatically work with C++ code generation, etc.).

---

<div class="post-metadata">

**Author:** ![serge25](https://avatars.discourse-cdn.com/v4/letter/s/4af34b/32.png) [@serge25](https://brian.discourse.group/u/serge25)\
**Post date:** [1 June 2022 16:10 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/5 "2022-06-01T16:10:01Z")

</div>

@mstimberg thank you for your answer. Is it right to reproduce this code `input_co1 = IdentityConnection(inpute, Pe, 'gext', weight=win)` from Brian1 like ` input_co1 = Synapses(inpute, Pe, on_pre='gext=win') input_co1.connect() input_co2 = Synapses(inputi, Pi, on_pre='gext=win') input_co2.connect()` for brian2?

---

<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 June 2022 08:35 UTC](https://brian.discourse.group/t/convert-code-from-brian1-to-brian2/704/6 "2022-06-02T08:35:15Z")

</div>

```python
input_co1 = IdentityConnection(inpute, Pe, 'gext', weight=win)

```

should translate into

```python
input_co1 = Synapses(inpute, Pe, on_pre='gext += win')
input_co1.connect(j='i')

```

(`IdentityConnection` means a one-to-one connection pattern, and each spike _increases_ the postsynaptic `gext`)  
Also note that we have [detailed documentation](https://brian2.readthedocs.io/en/stable/introduction/brian1_to_2/index.html) for the Brian 1 → 2 transition, which might be helpful.