Using decoded weights in a direct connection

General Discussion
1

Hello,

I’d like to use a set of weights from a decoded connection in a direct connection; the tuning curves for the ensemble look correct, however a simulation with a direct connection and the decoded weights does not produce correct output. Below is my code, and I would be very grateful if someone could let me know what I am doing wrong. Thanks!

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
import nengo
from nengo.utils.least_squares_solvers import SVD
from scipy.interpolate import pchip



def gaussian(x, mu, sig):
    return np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.)))

# Inputs mimicking receptive fields in a linear space

N = 1000
width = N/20 # receptive field width
input_rates = []
input_space = np.arange(N)
for i in range(3):
    loc = float(i+1)*N/4 # peak of the receptive field of each input
    input_rate = gaussian(input_space, loc, width)
    input_rates.append(input_rate)
input_matrix = np.column_stack(input_rates)

target_loc = N/2 # peak of the receptive field we wish to encode
target_rate = gaussian(input_space, target_loc, width)

trajectory_x = np.arange(250, 650)
trajectory_t = np.linspace(0., 10., len(trajectory_x))
trajectory_inputs = []
for input_rate in input_rates:
    trajectory_rate = np.interp(trajectory_x, input_space, input_rate)
    trajectory_input = pchip(trajectory_t, trajectory_rate)
    trajectory_inputs.append(trajectory_input)
    
# Target function    
target_trajectory_rate = np.interp(trajectory_x, input_space, target_rate)
target_trajectory_input = pchip(trajectory_t, target_trajectory_rate)
 
# Representation of inputs on a trajectory through the input space
def trajectory_function(t):
    result = np.asarray([ 2.*y(t) - 1. for y in trajectory_inputs ])
    return result


X_train = 2.*input_matrix - 1.
y_train = target_rate.reshape((-1,1))

X_test = trajectory_function(trajectory_t).T
y_test = target_trajectory_input(trajectory_t).ravel()

n_features = X_train.shape[1]
n_output = y_train.shape[1]

solver = nengo.solvers.LstsqL2nz(reg=0.01, solver=SVD())
    
ens_params = dict(
        neuron_type=nengo.LIF(),
        intercepts=nengo.dists.Choice([0.0]),
        max_rates=nengo.dists.Choice([20])
        )

# obtain decoder weights for the given inputs and target receptive field
with nengo.Network(seed=19, label="Receptive field tuning") as tuning_model:
    tuning_model.PC = nengo.Ensemble(n_neurons=50, dimensions=n_features, 
                                     eval_points=X_train,
                                     **ens_params)                                
    
    tuning_model.output = nengo.Ensemble(n_neurons=1, dimensions=1,
                                         **ens_params)
                                         
                                   
    PC_to_output = nengo.Connection(tuning_model.PC, 
                                    tuning_model.output.neurons, 
                                    synapse=None,
                                    eval_points=X_train, 
                                    function=y_train, 
                                    solver=solver)
    
tuning_sim = nengo.Simulator(tuning_model)
output_weights = tuning_sim.data[PC_to_output].weights

_, train_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_train)
_, test_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_test)



plt.figure()
plt.plot(y_train)
plt.plot(np.dot(train_acts, output_weights.T))

plt.figure()
plt.plot(y_test)
plt.plot(np.dot(test_acts, output_weights.T))

# use the decoded weights in a direct connection
network_weights = tuning_sim.data[PC_to_output].weights

with nengo.Network(label="Basic network model", seed=19) as model:

        model.Input = nengo.Node(trajectory_function, size_out=X_train.shape[1])
                                     
        model.Exc = nengo.Ensemble(50, dimensions=X_train.shape[1], 
                                   **ens_params)
        model.Output = nengo.Ensemble(1, dimensions=1,
                                   **ens_params)
                                   
            
        nengo.Connection(model.Input, model.Exc)
        nengo.Connection(model.Exc.neurons, model.Output, 
                         transform=network_weights,
                         synapse=0.01)
                         


with model:
    Input_probe = nengo.Probe(model.Input, synapse=0.01)
    Exc_probe = nengo.Probe(model.Exc, synapse=0.01)
    Output_probe = nengo.Probe(model.Output, synapse=0.01)

    
with nengo.Simulator(model) as sim:
    sim.run(np.max(trajectory_t))



plt.figure(figsize=(15,8))

#plt.plot(sim.trange(), sim.data[Input_probe], label='Input')
plt.plot(sim.trange(), sim.data[Exc_probe], label='Exc')
plt.plot(sim.trange(), sim.data[Output_probe], label='Output')
plt.legend();
2

What do you mean when you say “does not produce correct output”? Can you describe what you expect to happen, and what actually is happening?

3

In the code I attached, the tuning curve output closely resembles the function being computed across the connection:

image

(where X_test are the three inputs, y_test is the function we want to compute, and output is the dot product of the tuning curve “activities” vector and the decoded weights)

Then, the decoded weights are extracted and used in a new connection:

network_weights = tuning_sim.data[PC_to_output].weights
...  
model.Exc = nengo.Ensemble(50, dimensions=X_test.shape[1], ...)
model.Output = nengo.Ensemble(1, dimensions=1, ...)
nengo.Connection(model.Exc.neurons, model.Output,  transform=network_weights)

When this network is provided the input defined in X_test, I would expect an output that is similar to the one obtained by the dot product of the tuning curve “activities” vector with the decoded weights. However, the actual output is nothing like this:

image
image873×466 132 KB

4

I edited your post to add in the code so it’s a bit easier to see.

One thing to keep in mind is that when you’re using the nengo.utils.ensemble.tuning_curves function, you’re plotting the theoretical responses of the neurons with idealized firing rates, it’s not actually simulating the spiking behaviour. So when you compare that to actually running the model and plotting the output, the latter is going to look a lot “spikier” (which is what you’re observing). You could increase the synaptic filtering to smooth this out more. Or you could use a rate neuron model, like LIFRate, where the actual output will be the same as the theoretical tuning curves, in order to debug things to start off.

Another thing to keep in mind is that decoded weights are tied to the encoders of the Ensemble as well. So if you use the same decoding weights but different encoding weights, the decoding won’t work properly. So you will also need to make sure that model.Exc has the same encoders as tuning_model.PC

Haven’t looked through the code in enough detail to say that that’s everything, but those are the two issues that jumped out at me! If it still doesn’t give the output you expect after that we can dig into it some more.

5

Hi Daniel,

Thank you very much for the detailed response! First, your point about the encoder is well taken: the example code uses the same intercepts and maximum rate in both models. Is this sufficient to sure the encoding weights are the same? Even though the output of model.Exc is noisy, as you say, the decoded activity from the Exc ensemble appears to follow qualitatively the same trajectory as the tuning curve output. But if you have a suggestion for a more rigorous comparison of encoders, I am happy to implement it.

Second, I followed your suggestion to use LIFRate, and the result was that the output of model.Exc now follows the tuning curve even better, as the noisiness is all but gone, but the activity of the output ensemble is still nonsensical. I will paste those results in a separate post below, as I slightly changed the input parameters to make the decoded function somewhat less trivial. Thanks again for all your help.

6

Here is a version of the same model when using LIFRate:

First is the model used to compute tuning curves and its output below:

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
import nengo
from scipy.interpolate import pchip

def gaussian(x, mu, sig):
    return np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.)))

N = 1000
width = N/20
input_rates = []
x = np.arange(N)
for i in range(3):
    loc = float(i+1)*N/5
    input_rate = gaussian(x, loc, width)
    input_rates.append(input_rate)
input_matrix = np.column_stack(input_rates)
target_rate = gaussian(x, N/2, width)

trajectory_x = np.arange(250, 650)
trajectory_t = np.linspace(0., 10., len(trajectory_x))
trajectory_inputs = []
for input_rate in input_rates:
    trajectory_rate = np.interp(trajectory_x, x, input_rate)
    trajectory_input = pchip(trajectory_t, trajectory_rate)
    trajectory_inputs.append(trajectory_input)
    
    
target_trajectory_rate = np.interp(trajectory_x, x, target_rate)
target_trajectory_input = pchip(trajectory_t, target_trajectory_rate)
 
def trajectory_function(t):
    result = np.asarray([ 2.*y(t) - 1. for y in trajectory_inputs ])
    return result


X_train = 2.*input_matrix - 1.
y_train = target_rate.reshape((-1,1))

X_test = trajectory_function(trajectory_t).T
y_test = target_trajectory_input(trajectory_t).ravel()

n_features = X_train.shape[1]
n_output = y_train.shape[1]

from nengo.utils.least_squares_solvers import SVD
solver = nengo.solvers.LstsqL2nz(reg=0.01, solver=SVD())
    
ens_params = dict(
        neuron_type=nengo.LIFRate(),
        intercepts=nengo.dists.Choice([0.0]),
        max_rates=nengo.dists.Choice([20])
        )
    
with nengo.Network(seed=19, label="Receptive field tuning") as tuning_model:
       
    
    tuning_model.PC = nengo.Ensemble(n_neurons=1000, dimensions=n_features, 
                                     eval_points=X_train,
                                     **ens_params)                                
    
    tuning_model.output = nengo.Ensemble(n_neurons=1, dimensions=1,
                                         **ens_params)
                                         
                                   
    PC_to_output = nengo.Connection(tuning_model.PC, 
                                    tuning_model.output.neurons, 
                                    synapse=0.01,
                                    eval_points=X_train, 
                                    function=y_train, 
                                    solver=solver)
    
tuning_sim = nengo.Simulator(tuning_model)
output_weights = tuning_sim.data[PC_to_output].weights

_, train_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_train)
_, test_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_test)

plt.figure()
plt.plot(np.dot(test_acts, output_weights.T), label="output")
plt.plot(y_test, label="y_test")
plt.plot(X_test, label="X_test")
plt.legend()

image

As we can see, there is reasonable agreement between the target function y_test and the decoded weights output.

Next is the model in which I wish to use the decoded weights:


network_weights = tuning_sim.data[PC_to_output].weights

with nengo.Network(label="Basic network model", seed=19) as model:

        model.Input = nengo.Node(trajectory_function, size_out=X_test.shape[1])
                                     
        model.Exc = nengo.Ensemble(1000, dimensions=n_features, 
                                   neuron_type=nengo.LIFRate(), 
                                   intercepts=nengo.dists.Choice([0.0]),
                                   max_rates=nengo.dists.Choice([20])
                                  )
                                   
        model.Out = nengo.Ensemble(1, dimensions=1, 
                                      neuron_type=nengo.LIFRate(),
                                      intercepts=nengo.dists.Choice([0.0]),
                                      max_rates=nengo.dists.Choice([20]))
            
        nengo.Connection(model.Input, model.Exc)
        nengo.Connection(model.Exc.neurons, 
                         model.Out.neurons, synapse=0.01,
                         transform=network_weights)
                         
with model:
    Input_probe = nengo.Probe(model.Input)
    Exc_probe = nengo.Probe(model.Exc)
    Output_probe = nengo.Probe(model.Out, 'input')

    
with nengo.Simulator(model) as sim:
    sim.run(np.max(trajectory_t))
                         
plt.figure(figsize=(15,8))

plt.plot(sim.trange(), sim.data[Input_probe], label='Input')
plt.legend()

plt.figure(figsize=(15,8))
plt.plot(sim.trange(), sim.data[Exc_probe], label='Exc')
plt.plot(sim.trange(), sim.data[Output_probe], label='Output')
plt.legend();

The output of the model.Exc ensemble looks very close to the tuning curve output, but now the model.Output ensemble appears to have no activity:

image
image889×466 28.5 KB

image
image882×466 35.5 KB

7

No, the encoders are distinct from the intercepts/max rates. You can set them when you create the Ensemble, like

ens = nengo.Ensemble(
    n_neurons, dimensions, max_rates=..., intercepts=..., encoders=x)

where x is an array with shape (n_neurons, dimensions). And you can see the built encoders for an Ensemble in sim.data[ens].encoders (similar to how you access the connection weights).

8

Awesome, that worked! Thanks a lot, below is the final working code for posterity.

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
import nengo
from scipy.interpolate import pchip

def gaussian(x, mu, sig):
    return np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.)))

N = 1000
width = N/20
input_rates = []
x = np.arange(N)
for i in range(3):
    loc = float(i+1)*N/5
    input_rate = gaussian(x, loc, width)
    input_rates.append(input_rate)
input_matrix = np.column_stack(input_rates)
target_rate = gaussian(x, N/2, width)

trajectory_x = np.arange(250, 650)
trajectory_t = np.linspace(0., 10., len(trajectory_x))
trajectory_inputs = []
for input_rate in input_rates:
    trajectory_rate = np.interp(trajectory_x, x, input_rate)
    trajectory_input = pchip(trajectory_t, trajectory_rate)
    trajectory_inputs.append(trajectory_input)
    
    
target_trajectory_rate = np.interp(trajectory_x, x, target_rate)
target_trajectory_input = pchip(trajectory_t, target_trajectory_rate)

def trajectory_function(t):
    result = np.asarray([ 2.*y(t) - 1. for y in trajectory_inputs ])
    return result


X_train = 2.*input_matrix - 1.
y_train = target_rate.reshape((-1,1))

X_test = trajectory_function(trajectory_t).T
y_test = target_trajectory_input(trajectory_t).ravel()

n_features = X_train.shape[1]
n_output = y_train.shape[1]

from nengo.utils.least_squares_solvers import SVD
solver = nengo.solvers.LstsqL2nz(reg=0.01, solver=SVD())
    
ens_params = dict(
        neuron_type=nengo.LIFRate(),
        intercepts=nengo.dists.Choice([0.0]),
        max_rates=nengo.dists.Choice([20])
        )
    
with nengo.Network(seed=19, label="Receptive field tuning") as tuning_model:
       
    
    tuning_model.PC = nengo.Ensemble(n_neurons=1000, dimensions=n_features, 
                                     eval_points=X_train,
                                     **ens_params)                                
    
    tuning_model.output = nengo.Ensemble(n_neurons=1, dimensions=1,
                                         **ens_params)
                                         
                                   
    PC_to_output = nengo.Connection(tuning_model.PC, 
                                    tuning_model.output.neurons, 
                                    synapse=0.01,
                                    eval_points=X_train, 
                                    function=y_train, 
                                    solver=solver)
    
tuning_sim = nengo.Simulator(tuning_model)
output_weights = tuning_sim.data[PC_to_output].weights

_, train_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_train)
_, test_acts = nengo.utils.ensemble.tuning_curves(tuning_model.PC, tuning_sim, inputs=X_test)

plt.figure()
plt.plot(y_train)
plt.plot(np.dot(train_acts, output_weights.T))

plt.figure()
plt.plot(np.dot(test_acts, output_weights.T), label="output")
plt.plot(y_test, label="y_test")
plt.plot(X_test, label="X_test")
plt.legend()

network_weights = tuning_sim.data[PC_to_output].weights
network_encoders = tuning_sim.data[tuning_model.PC].encoders


with nengo.Network(label="Basic network model", seed=19) as model:

        model.Input = nengo.Node(trajectory_function, size_out=X_test.shape[1])
                                     
        model.Exc = nengo.Ensemble(1000, dimensions=n_features, 
                                   neuron_type=nengo.LIFRate(), 
                                   intercepts=nengo.dists.Choice([0.0]),
                                   max_rates=nengo.dists.Choice([20]),
                                   encoders=network_encoders
                                  )
                                   
        model.Out = nengo.Ensemble(1, dimensions=1, 
                                      neuron_type=nengo.LIFRate(),
                                      intercepts=nengo.dists.Choice([0.0]),
                                      max_rates=nengo.dists.Choice([20]))
            
        nengo.Connection(model.Input, model.Exc)
        nengo.Connection(model.Exc.neurons, 
                         model.Out, synapse=0.01,
                         transform=network_weights)
                  
with model:
    Input_probe = nengo.Probe(model.Input)
    Exc_probe = nengo.Probe(model.Exc)
    Output_probe = nengo.Probe(model.Out, 'input')

    
with nengo.Simulator(model) as sim:
    sim.run(np.max(trajectory_t))

plt.figure(figsize=(15,8))
plt.plot(sim.trange(), sim.data[Input_probe], label='Input')
plt.legend()

plt.figure(figsize=(15,8))
plt.plot(sim.trange(), sim.data[Exc_probe], label='Exc')
plt.plot(sim.trange(), sim.data[Output_probe], label='Output')
plt.legend();