Questions about the article"A Biologically Plausible Spiking Neural Model of Eyeblink Conditioning in the Cerebellum"

There’s an easier way to construct the Legendre polynomials using scipy.special.legendre. You can use this to create a vector of polynomials (evaluated at the points you want to reconstruct), and then apply them with a dot product, like so:

def decoder(m, points_to_reconstruct=2):
    output = np.zeros(number_of_points)
    width = int(theta / dt)
    window = np.linspace(0, 1, points_to_reconstruct)
    D = np.asarray([legendre(i)(2 * window - 1) for i in range(q)])
    return m.dot(D)

This decodes a two-dimensional array, where the first element is an approximation of u(t) and the second is an approximation of u(t - theta). Increasing the value of points_to_reconstruct gives you evenly spaced points in between. You can also replace np.linspace(0, 1, points_to_reconstruct) with whatever points you want within the window (normalized to [0, 1]).

However, one thing to note is that the reconstruction is approximating the window using polynomials up to order $q - 1$. If you have a step function, such as the one in your example, then it is impossible for a finite number of polynomials to fit such a curve, as it has discontinuous derivatives wherever there are discontinuous changes in the input (whereas, in contrast, polynomials have continuous derivatives). You do better if you increase $q$, though, in the same way that higher-order Taylor series form better approximations as you add more polynomials to the series.

In general, the reconstruction will work best for low-frequency inputs, where the frequency scales linearly with $q$. That is, if you double the frequency of your input signal, then you need twice as many polynomials to fit it as accurately.

Here’s an example taking your previous code, increasing q to 24, using ZOH discretization, and the above decoder function:

import numpy as np
import matplotlib.pyplot as plt
import scipy.signal
from scipy.special import legendre

q = 24
theta = 0.4
number_of_points = 10000
dt = 1e-4



def get_a_and_b():
# A:
    a = np.zeros([q, q])
    for i in range(q):
        for j in range(q):
            if i < j:
                a[i, j] = (2 * i + 1) * -1
            else:
                a[i, j] = (2 * i + 1) * ((-1) ** (i - j + 1))
# B:
    b = np.zeros(q)
    for i in range(q):
        b[i] = (2 * i + 1) * ((-1) ** i)
    return a/theta, b/theta


def encoder():
    Acont, Bcont = get_a_and_b()

    Adisc, Bdisc, _, _, _ = scipy.signal.cont2discrete(
        (Acont, Bcont[:, None], np.zeros((1, q)), 0),
        dt=dt,
        method='zoh',
    )
    Bdisc = Bdisc.squeeze(axis=1)
    
    m = np.zeros(q)
    output = np.zeros([number_of_points, q])

    for i in range(number_of_points):
        if 1000 < i < 2000:
            u = 1
        else:
            u = 0
        # B_u = u * B
        # m = m + dt * (np.dot(A, m) + B_u)
        m = Adisc.dot(m) + Bdisc * u
        output[i, :] = m

    return output


def decoder(m, points_to_reconstruct=2):
    output = np.zeros(number_of_points)
    width = int(theta / dt)
    window = np.linspace(0, 1, points_to_reconstruct)
    D = np.asarray([legendre(i)(2 * window - 1) for i in range(q)])
    return m.dot(D)


output = encoder()
output2 = decoder(output)

plt.plot(range(number_of_points), output2)

plt.show()

delayed_step

Note the ringing around the edges occurs from the higher-order polynomials trying to fit the instantaneous changes in the input. It would look much better with a filtered or low-frequency input.

Thank you very much for your answer! I am very interested in this article, I found its poster and saw the corresponding Github project. But unfortunately, this website cannot be accessed. May I ask, is the website I visited wrong?
https://github.com/ctn-waterloo/cogsci2020-cerebellum

Just tagging the poster’s author @astoecke to rectify the repository link. :slight_smile:

Thank you for pointing this out! I actually forgot to create the repository, but I have done that now.

The code you are looking for is in model/granule_golgi_circuit.py.

Edit: Oh, and the code for the Delay Network decoders is in https://github.com/ctn-waterloo/cogsci2020-cerebellum/blob/master/notebooks/Delay%20Network%20Plots.ipynb

1 Like

Thank you very much for your answer! But I encountered some problems when running this program. What is the problem?I would be very grateful if you can answer my question
Operating environment:
Ubuntu 18.04
python 3.6.6
nengo 3.1.0

Hm, did you make sure that you are running a version of NengoBio from the end of January 2020 as stated in the readme? The code is as it was back when we submitted the paper, and I haven’t updated it to work with newer versions of NengoBio, which no longer feature the max_n_post_synapses parameter.

I think that the latest version of NengoBio that should be compatible with this code is Git commit 8ff2274f8f5e09877ed0a554ca268f9fa396e3a5. You can simply checkout that commit from the NengoBio Git repository and install it from there by running pip3 install -e .. Make sure to uninstall any currently installed NengoBio version before you do that.

Thanks for your answer, I have solved this problem!

When I was running the model exploration file in the notebook, I found that the running result was different from the expected result. Could it be my mistake? I have been thinking about this issue for some time, I would be very grateful if I can get your advice.

I think that your results are exactly what they should be. The committed version of the notebook was executed on Jan 22, and the default parameters model changed quite significantly since then. If you revert back to the Jan 22 revision 79b4fb5d48a5feb91eef15e02bc356d093160a59 (and revision 37a5e17d4b62a495e46fb21af0d55b07d7bf5e15 of nengo-bio) you will get the results that were stored in the notebook.

Edit: As the name of the notebook suggests, we only used this notebook for playing around with the model. The data shown there are not what we reported in the paper. We used the notebook Run Model.ipynb for that, which I apparently didn’t commit until Apr 2. It thus wasn’t included in the repository that I uploaded to GitHub. I have now added Run Model.ipynb.

Thank you very much for your detailed guidance, my problem has been solved. But there is another question that has always puzzled me. Is there any basis for setting q=6 in the paper? I tried q=4 and q=10, and the results they achieved are as follows:
4 10
I would be very grateful if you can answer my doubts!At the same time any suggestions or ideas are welcome!

Excellent question! We don’t really discuss this in the paper. There are two reasons, really.

  1. The dynamical system becomes unstable for larger $q$ when building the Granule-Golgi circuit, as you have noticed. $q = 6$ is kind of a compromise between being stable and having a sufficiently good temporal basis. We are not exactly sure where the instabilities come from, but we think that it has to do with the higher frequencies of some of the state dimensions.

  2. As Voelker and Eliasmith show in their 2018 Neural Computation paper [1] the Delay Network matches time cell data from the brain for $q = 6$, but not for other dimensionalities.

[1] https://www.mitpressjournals.org/doi/full/10.1162/neco_a_01046

Thank you for your detailed and patient answer! :smiley:
In the paper I saw that you used the data from Heiney et al. (2014), but I did not find his experimental data in his paper or appendix. Is his data public? Would you share his data or the way to obtain the data with me please?In addition, how does your paper model the 6 mouse models? There are ANTAGONIST and AGONIST in the code. Is this used to model different mouse models? Is there a basis for setting the values of ANTAGONIST and AGONIST?

That data is in the paper in the form of diagrams. We manually digitized the diagrams to get the data (our extracted data points should be in one of the notebooks). We (very roughly) fitted the motor model to the eyeblink trajectories reported in one of the papers.

Thank you for your answer. I searched the notebook carefully and did not find the data points you extracted based on the paper (Heiney et al.(2014)). Could you tell me the exact path of the data in the notebook?

In your paper, the eyeblink trajectories of 6 rats is simulated. Did you fit the data based on the eyeblink trajectories of different rats?

Ah, the traced diagrams were in the experimental_data folder, which I originally didn’t upload to GitHub because it contained tons of copyrighted material. I’ve cleaned that folder up and uploaded it to GitHub. The data are loaded at the beginning of Plot - Results - Basic.ipynb.

Regarding to the motor system parameters, I’ve added my original notebook Eyeblink Motor System.ipynb to the repository; I’m noting therein which paper and which figure I fitted the data to.
We’re using the same set of motor system parameters for all experimental trials.

Thanks a million. :grinning:Now I have plotted the neurobiological data, but in the process of running plot-result-basic, I regret that I did not reproduce the model data in the paper. Is it because I am using the wrong software package or code version?
nengo-bio: commit 8ff2274f8f5e09877ed0a554ca268f9fa396e3a5
cogsci2020-cerebellum-master: f0f14a66a4d6210c6e513e67982acbaa06532bbe
System:Ubuntu18.04


At the same time, I also tried code( 79b4fb5d48a5feb91eef15e02bc356d093160a59), the result is the same as above.
Any suggestions or ideas are welcome. :grinning:

Interesting; I just used this code last week and I’m getting exactly the same results we got last year. How exactly did you generate the data that you’re visualising here, i.e., can you provide me with a step-by-step list of things that you do to produce these figures (starting with a fresh copy of the repository)?

Certainly.
step :
1.Install nengo using command:pip install nengo
2.Install nengo-bio using command: pip install nengo-bio
3.reset nengo-bio version to commit 8ff2274f8f5e09877ed0a554ca268f9fa396e3a5
4.Download the cerebellum model code
5.Run the plot-result-basic.py program with mode=two_populations_dales_principle
6.Then the above result appeared
I have been struggling with this problem for a long time, if you can help me solve this problem, I would be very grateful! :smile:

I looked into this some more and there are multiple potential issues:

  • You have to use pip install -e . if you want changes made to the source repository to be visible. In other words, you cannot use pip install nengo-bio and then make changes to the source repository (e.g., check out an older revision).
  • The code for running the model in Plot - Results - Basic.ipynb was just for debugging; as I wrote above [1] you have to use Run Model.ipynb and then load the results in the appropriate places in Plot - Results - Basic.ipynb. You’re likely using the debug data to plot your results (at least your plots are very similar to what I get when I plot the debug data).
  • The Plot - Results - Basic.ipynb in the repository was not exactly the version we used for the paper; that’s why the overall layout was not what’s in the paper. I’ve restored a version from Feb 3 from our backups and replaced the file on GitHub.

To make things easier for you, I’ve prepared a Docker container that sets everything up. See the Dockerfile [2] in the repository for the exact sequence of commands required in a fresh Debian 10.8 installation. See the updated README [3] for how to use the Docker container. Note that my knowledge of Docker boils down to ten minutes of using Google; so I cannot help you if you run into any issues with Docker.

Sorry for making this entire process so complicated; but (unfortunately) this is how things typically are with research code.

[1] Questions about the article"A Biologically Plausible Spiking Neural Model of Eyeblink Conditioning in the Cerebellum" - #12 by astoecke

[2] cogsci2020-cerebellum/Dockerfile at master · ctn-waterloo/cogsci2020-cerebellum · GitHub

[3] cogsci2020-cerebellum/README.md at master · ctn-waterloo/cogsci2020-cerebellum · GitHub

Thank you very much for your patience and meticulous help! I have solved this problem. This paper deeply inspired me and made me very interested in this research direction.