Debugging

Even if you write clear and readable code, even if you fully understand your codes, env if you are very familiar with your model, weird bugs will inevitably appear and you will need to debug them in some way.

Fortunately, BrainPy supports debugging with pdb module or breakpoint (The latest version of BrainPy removes the support of debugging in IDEs). That is to say, you do not need to resort to using bunch of print statements to see what’s happening in their code. On the contrary, you can work with Python’s interactive source code debugger to see the state of any variable in your model.

For the variables you are interested in, you just need to add the pdb.set_trace() or breakpoint() after the code line.

In this section, let’s take the HH neuron model as an example to illustrate how to debug your model within BrainPy.

[1]:
import brainpy as bp
import numpy as np
import pdb

If you want to debug your model, we would like to recommond you to open the show_code=True.

[2]:
bp.profile.set(show_code=True, jit=False)

Here, the HH neuron model is defined as:

[3]:
E_Na = 50.
E_K = -77.
E_leak = -54.387
C = 1.0
g_Na = 120.
g_K = 36.
g_leak = 0.03
V_th = 20.
noise = 1.

ST = bp.types.NeuState(
    {'V': -65., 'm': 0.05, 'h': 0.60,
     'n': 0.32, 'spike': 0., 'input': 0.}
)

@bp.integrate
def int_m(m, _t, V):
    alpha = 0.1 * (V + 40) / (1 - np.exp(-(V + 40) / 10))
    beta = 4.0 * np.exp(-(V + 65) / 18)
    return alpha * (1 - m) - beta * m

@bp.integrate
def int_h(h, _t, V):
    alpha = 0.07 * np.exp(-(V + 65) / 20.)
    beta = 1 / (1 + np.exp(-(V + 35) / 10))
    return alpha * (1 - h) - beta * h

@bp.integrate
def int_n(n, _t, V):
    alpha = 0.01 * (V + 55) / (1 - np.exp(-(V + 55) / 10))
    beta = 0.125 * np.exp(-(V + 65) / 80)
    return alpha * (1 - n) - beta * n

@bp.integrate
def int_V(V, _t, m, h, n, I_ext):
    I_Na = (g_Na * np.power(m, 3.0) * h) * (V - E_Na)
    I_K = (g_K * np.power(n, 4.0))* (V - E_K)
    I_leak = g_leak * (V - E_leak)
    dVdt = (- I_Na - I_K - I_leak + I_ext)/C
    return dVdt, noise / C

def update(ST, _t):
    m = np.clip(int_m(ST['m'], _t, ST['V']), 0., 1.)
    h = np.clip(int_h(ST['h'], _t, ST['V']), 0., 1.)
    n = np.clip(int_n(ST['n'], _t, ST['V']), 0., 1.)
    V = int_V(ST['V'], _t, m, h, n, ST['input'])

    pdb.set_trace()

    spike = np.logical_and(ST['V'] < V_th, V >= V_th)
    ST['spike'] = spike
    ST['V'] = V
    ST['m'] = m
    ST['h'] = h
    ST['n'] = n
    ST['input'] = 0.

HH = bp.NeuType(ST=ST,
                name='HH_neuron',
                steps=update,
                mode='vector')

In this example, we add pdb.set_trace() after the variables \(m\), \(h\), \(n\) and \(V\).

Then we can create a neuron group, and try to run this neuron model:

[4]:
group = bp.NeuGroup(HH, geometry=1, monitors=['spike'])
group.run(1000., inputs=('input', 10.))
def NeuGroup0_input_step(ST, input_inp,):
  # "input" step function of NeuGroup0
  ST[5] += input_inp



def NeuGroup0_monitor_step(ST, _i, mon_ST_spike,):
  # "monitor" step function of NeuGroup0
  mon_ST_spike[_i] = ST[4]



def NeuGroup0_update(ST, _t,):
 # "update" step function of NeuGroup0
 _int_m_m = ST[1]
 _int_m__t = _t
 _int_m_V = ST[0]
 _int_m_alpha = 0.1 * (_int_m_V + 40) / (1 - np.exp(-(_int_m_V + 40) / 10))
 _int_m_beta = 4.0 * np.exp(-(_int_m_V + 65) / 18)
 _dfm_dt = _int_m_alpha * (1 - _int_m_m) - _int_m_beta * _int_m_m
 _int_m_m = 0.1*_dfm_dt + _int_m_m
 _int_m_res = _int_m_m
 m = np.clip(_int_m_res, 0.0, 1.0)

 _int_h_h = ST[2]
 _int_h__t = _t
 _int_h_V = ST[0]
 _int_h_alpha = 0.07 * np.exp(-(_int_h_V + 65) / 20.0)
 _int_h_beta = 1 / (1 + np.exp(-(_int_h_V + 35) / 10))
 _dfh_dt = _int_h_alpha * (1 - _int_h_h) - _int_h_beta * _int_h_h
 _int_h_h = 0.1*_dfh_dt + _int_h_h
 _int_h_res = _int_h_h
 h = np.clip(_int_h_res, 0.0, 1.0)

 _int_n_n = ST[3]
 _int_n__t = _t
 _int_n_V = ST[0]
 _int_n_alpha = 0.01 * (_int_n_V + 55) / (1 - np.exp(-(_int_n_V + 55) / 10))
 _int_n_beta = 0.125 * np.exp(-(_int_n_V + 65) / 80)
 _dfn_dt = _int_n_alpha * (1 - _int_n_n) - _int_n_beta * _int_n_n
 _int_n_n = 0.1*_dfn_dt + _int_n_n
 _int_n_res = _int_n_n
 n = np.clip(_int_n_res, 0.0, 1.0)

 _int_V_V = ST[0]
 _int_V__t = _t
 _int_V_m = m
 _int_V_h = h
 _int_V_n = n
 _int_V_I_ext = ST[5]
 _int_V_I_Na = g_Na * np.power(_int_V_m, 3.0) * _int_V_h * (_int_V_V - E_Na)
 _int_V_I_K = g_K * np.power(_int_V_n, 4.0) * (_int_V_V - E_K)
 _int_V_I_leak = g_leak * (_int_V_V - E_leak)
 _int_V_dVdt = (-_int_V_I_Na - _int_V_I_K - _int_V_I_leak + _int_V_I_ext) / C
 _dfV_dt = _int_V_dVdt
 _V_dW = _normal_like(_int_V_V)
 _dgV_dt = noise / C
 _int_V_V = _int_V_V + 0.316227766016838*_V_dW*_dgV_dt + 0.1*_dfV_dt
 _int_V_res = _int_V_V
 V = _int_V_res

 pdb.set_trace()

 spike = np.logical_and(ST[0] < V_th, V >= V_th)
 ST[4] = spike
 ST[0] = V
 ST[1] = m
 ST[2] = h
 ST[3] = n
 ST[5] = 0.0



def step_func(_t, _i, _dt):
  NeuGroup0_runner.input_step(NeuGroup0.ST["_data"], NeuGroup0_runner.input_inp,)
  NeuGroup0_runner.update(NeuGroup0.ST["_data"], _t,)
  NeuGroup0_runner.monitor_step(NeuGroup0.ST["_data"], _i, NeuGroup0.mon["spike"],)

> c:\users\oujag\codes\projects\brainpy\docs\advanced(52)update()

ipdb> p m
array([0.05123855])
ipdb> p n
array([0.31995744])
ipdb> p h
array([0.59995445])
ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(53)update()

ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(54)update()

ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(55)update()

ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(56)update()

ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(57)update()

ipdb> n
> c:\users\oujag\codes\projects\brainpy\docs\advanced(58)update()

ipdb> p ST
array([[-6.41827214e+01],
       [ 5.12385538e-02],
       [ 5.99954448e-01],
       [ 3.19957442e-01],
       [ 0.00000000e+00],
       [ 1.00000000e+01]])
ipdb> q
---------------------------------------------------------------------------
BdbQuit                                   Traceback (most recent call last)
<ipython-input-4-703915fa7809> in <module>
      1 group = bp.NeuGroup(HH, geometry=1, monitors=['spike'])
----> 2 group.run(1000., inputs=('input', 10.))

~\codes\projects\BrainPy\brainpy\core\base.py in run(self, duration, inputs, report, report_percent)
    584         else:
    585             for run_idx in range(run_length):
--> 586                 step_func(_t=times[run_idx], _i=run_idx, _dt=dt)
    587
    588         if profile.run_on_gpu():

? in step_func(_t, _i, _dt)

? in update(ST, _t)

? in update(ST, _t)

~\Miniconda3\envs\py38\lib\bdb.py in trace_dispatch(self, frame, event, arg)
     86             return # None
     87         if event == 'line':
---> 88             return self.dispatch_line(frame)
     89         if event == 'call':
     90             return self.dispatch_call(frame, arg)

~\Miniconda3\envs\py38\lib\bdb.py in dispatch_line(self, frame)
    111         if self.stop_here(frame) or self.break_here(frame):
    112             self.user_line(frame)
--> 113             if self.quitting: raise BdbQuit
    114         return self.trace_dispatch
    115

BdbQuit: