reproducing offcial code

[tesla-cat]

Hi, I am trying to reproduce the simulations with my own simple python code

for the sims/ch2/neuron/ example, I have:

This looks quite different from the official example, I wonder where am I doing it wrong?

it is produced by the following code:

from typing import Dict, List
import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

class Recorder:
    def __init__(s):
        s.records: Dict[str, List] = {}

    def record(s, k, v):
        if k not in s.records: s.records[k] = []
        s.records[k].append(v)
    
    def show(s):
        scale = 1.8
        plt.figure(figsize=(2*scale, len(s.records)*scale))
        for i, (k, v) in enumerate(s.records.items()):
            plt.subplot(len(s.records), 1, i+1)
            plt.plot(v)
            plt.title(k)   
        plt.tight_layout() 
        plt.show()

recorder = Recorder()

#============================================

def init_act_func(gamma, sigma):
    x = np.linspace(-2, 2, 1000)
    x2 = gamma * np.where(x < 0, 0, x)
    y = x2 / (x2+1)
    gaussian = np.exp(-x**2/(2*sigma**2))
    gaussian /= np.sum(gaussian)
    y2 = np.convolve(y, gaussian, mode='same')
    mask = np.gradient(y2) >= 0
    x = x[mask]; y = y[mask]; y2 = y2[mask]
    act_func = interp1d(x, y2)
    if 0:
        plt.plot(x, y)
        plt.plot(x, act_func(x))
        plt.show()
    return act_func

#======================================

gbe = 0.3; gbi = 1.0; gl = 0.3
Ee = 1.0; Ei = 0.25; El = 0.3
Vmr = 0.3
dt_vm = 0.3
Theta = 0.5 
act_func = init_act_func(gamma=30, sigma=0.01)
mode = 'spike'

class Neuron:
    def __init__(s):
        s.Vm = Vmr
        s.y = 0

    def update(s, ge, gi=0):
        ge *= gbe; gi *= gbi
        Inet = ge*(Ee-s.Vm) + gi*(Ei-s.Vm) + gl*(El-s.Vm)
        s.Vm += dt_vm * Inet
        
        if mode=='spike':
            if s.Vm > Theta: s.y = 1; s.Vm = Vmr
            else: s.y = 0
        elif mode=='rate':
            # Vm_eq = (ge*Ee + gi*Ei + gl*El) / (ge + gi + gl)
            ge_thr = (gi*(Ei-Theta) + gl*(El-Theta)) / (Theta-Ee)
            ys = act_func(ge - ge_thr)
            s.y += dt_vm * (ys - s.y) 

        recorder.record('ge', ge)
        recorder.record('Inet', Inet)
        recorder.record('Vm', s.Vm)
        recorder.record('act', s.y)

if __name__=='__main__':
    neuron = Neuron()
    for t in range(10):
        neuron.update(0)
    for t in range(10, 160):
        neuron.update(1)
    for t in range(160, 200):
        neuron.update(0)
    recorder.show()