Compute neuron transfer function

low_firing_rate_dynamics_sweep.py

def compute_transfer_function(us, xs, dt=dt, smoothing=100.0):
    Us = np.fft.fftshift(np.fft.fft(us))
    Xs = np.fft.fftshift(np.fft.fft(xs))
    fs = np.fft.fftshift(np.fft.fftfreq(len(us), dt))

    wnd = np.exp(-np.square(fs) / np.square(smoothing))
    wnd /= np.sum(wnd) * dt

    UXs = np.fft.fftshift(np.fft.fft(Us * np.conj(Xs)))
    XXs = np.fft.fftshift(np.fft.fft(Xs * np.conj(Xs)))

    As = np.fft.ifft(np.fft.ifftshift(UXs * wnd))
    Bs = np.fft.ifft(np.fft.ifftshift(XXs * wnd))

    return fs, As / Bs


def run_single_experiment(args):
    # Unpack the arguments
    i, j, k, rate = args
    seed = 58219 * i + 213
    neuron_type = neuron_types[j]

    n_neurons = int(100000 / rate)  # 10 Hz ==> 10,000; 100 Hz ==> 1000

    with nengo.Network(seed=seed) as model:
        ens_us = nengo.Node(
            nengo.processes.WhiteSignal(period=T,
                                        high=bandwidth,
                                        rms=0.5,
                                        seed=seed))

        ens_xs = nengo.Ensemble(n_neurons=n_neurons,
                                neuron_type=neuron_type(),
                                dimensions=1,
                                max_rates=nengo.dists.Uniform(
                                    0.5 * rate, rate),
                                seed=seed)

        nengo.Connection(ens_us, ens_xs, synapse=None)

        p_us = nengo.Probe(ens_us, synapse=None)
        p_xs = nengo.Probe(ens_xs, synapse=None)

    with nengo.Simulator(model, dt=dt, progress_bar=False) as sim:
        sim.run(T)

    ts = sim.trange()
    us = sim.data[p_us][:, 0]
    xs = sim.data[p_xs][:, 0]

    fs, Hs = compute_transfer_function(us, xs, dt)
    i0 = np.min(np.where(fs >= -bandwidth)[0])
    i1 = np.max(np.where(fs <= bandwidth)[0])

    return i, j, k, Hs[i0:i1]