Need help for a Custom Poisson Spiking Neuron with Exponential Tuning Curve

Hi! I have a hardware device that creates Poisson Spikes but its tuning curve is a Rectified Exponential as shown in Fig 1. I want to use Nengo to simulate many of my devices as a neuron population. I have tried to create my own rectified exponential neuron in Nengo by modifying the RectifiedLinear neuron class as shown in the code snippet “Modified RectifiedLinear Neuron type class” and it performs with some limited success. The code in Nengo for a single device is also shown. This implementation doesn’t work well with Nengo’s random initialization of neurons, which is expected since Nengo is still initializing for a RectifiedLinear Neuron. I was wondering if anyone has some tips for implementing a custom neurons in Nengo. If any part of my post is not clear, please reach out to me. Thank you!

image
Fig 1

Nengo Code

import nengo

model = nengo.Network()
with model:
    stim0 = nengo.Node([0])
    max_rate_mat = [5,5]
    encoder_mat = [[1],[-1]]
    intercept_mat = [0.0,0.0]
    a = nengo.Ensemble(n_neurons=2, dimensions=1, max_rates=max_rate_mat, encoders=encoder_mat, intercepts=intercept_mat,
    neuron_type=nengo.PoissonSpiking(nengo.RectifiedLinear()))
    nengo.Connection(stim0, a)

Modified RectifiedLinear Neuron type class

class RectifiedLinear(NeuronType):
    """A rectified linear neuron model.

    Each neuron is modeled as a rectified line. That is, the neuron's activity
    scales linearly with current, unless it passes below zero, at which point
    the neural activity will stay at zero.

    Parameters
    ----------
    amplitude : float
        Scaling factor on the neuron output. Corresponds to the relative
        amplitude of the output of the neuron.
    initial_state : {str: Distribution or array_like}
        Mapping from state variables names to their desired initial value.
        These values will override the defaults set in the class's state attribute.
    """

    negative = False

    amplitude = NumberParam("amplitude", low=0, low_open=True)

    def __init__(self, amplitude=1, initial_state=None):
        super().__init__(initial_state)

        self.amplitude = amplitude

    def gain_bias(self, max_rates, intercepts):
        """Determine gain and bias by shifting and scaling the lines."""
        max_rates = np.array(max_rates, dtype=float, copy=False, ndmin=1)
        intercepts = np.array(intercepts, dtype=float, copy=False, ndmin=1)
        gain = max_rates / (1 - intercepts)
        bias = -intercepts * gain
        return gain, bias

    def max_rates_intercepts(self, gain, bias):
        """Compute the inverse of gain_bias."""
        intercepts = -bias / gain
        max_rates = (1 - intercepts)
        return max_rates, intercepts

    def step(self, dt, J, output):
        """Implement the rectification nonlinearity."""
        output[...] = self.amplitude * (np.exp((np.maximum(0.0, J)))-1)

I’m not quite sure what the problem is that you’re experiencing. For the most part, your code looks correct. There are two main problems that I’ve noticed:

  • Modifying the RectifiedLinear neuron type in place could cause unintended side effects with other models; better to create your own new neuron type.
  • The gain_bias and max_rates_intercepts functions need to be modified for your new neuron type. I think this might be what you were referring to by “This implementation doesn’t work well with Nengo’s random initialization of neurons”. I’ve modified these functions below for your neuron type. Essentially, you want to solve a system of equations for the gain/bias given the max rates/intercepts: f(intercepts) = 0, f(1) = max_rates, where f(x) = exp(max(gain * x + bias, 0)) - 1 is your neural nonlinearity. With some straightforward simplification you can find that gain*intercepts + bias = 0 and gain + bias = log(max_rates + 1), and from this get the equations below.
import matplotlib.pyplot as plt
import nengo
import numpy as np


class RectifiedExp(nengo.neurons.RectifiedLinear):
    negative = False

    def gain_bias(self, max_rates, intercepts):
        """Determine gain and bias by shifting and scaling the lines."""
        max_rates = np.array(max_rates, dtype=float, copy=False, ndmin=1)
        intercepts = np.array(intercepts, dtype=float, copy=False, ndmin=1)
        gain = np.log1p(max_rates) / (1 - intercepts)
        bias = -intercepts * gain
        return gain, bias

    def max_rates_intercepts(self, gain, bias):
        """Compute the inverse of gain_bias."""
        intercepts = -bias / gain
        max_rates = np.exp(gain + bias) - 1
        return max_rates, intercepts

    def step(self, dt, J, output):
        """Implement the rectification nonlinearity."""
        output[...] = self.amplitude * (np.exp((np.maximum(0.0, J))) - 1)


model = nengo.Network()
with model:
    n_neurons = 100
    encoders = np.ones((n_neurons, 1))
    max_rates = 100 * np.ones(n_neurons)
    intercepts = np.linspace(-0.95, 0.95, n_neurons)

    def make_ensemble(neuron_type):
        return nengo.Ensemble(
            n_neurons=n_neurons,
            dimensions=1,
            max_rates=max_rates,
            encoders=encoders,
            intercepts=intercepts,
            neuron_type=neuron_type,
        )

    a = make_ensemble(RectifiedExp())
    ap = nengo.Probe(a.neurons)

    b = make_ensemble(nengo.PoissonSpiking(RectifiedExp()))
    bp = nengo.Probe(b.neurons)

with nengo.Simulator(model) as sim:
    sim.run(100.0)


plt.plot(-intercepts, sim.data[ap].mean(axis=0))
plt.plot(-intercepts, sim.data[bp].mean(axis=0))
plt.show()