Assume that you have an LTI system with a short memory. The response depends only on the stimulus at the current time point and one time point in the past.
You have data for and for a set of time points , and you can calculate time-averaged values such as . Your goal is to determine the coefficients of the filter, and .
Approach A: minimize the mean squared error between the predicted and actual rate to derive the set of equations you will need to solve to find and . Here are a few hints to get you started:
To minimize E, set the derivatives of E with respect to and to zero:
You should be able to write the system of equations you derive in matrix format as
What are and ? How would you calculate them in numpy?
Approach B: in a linear system, the noise in the prediction () will be uncorrelated with the stimulus . Thus, you can write the following equality:
where stands for convolution and and are the correlation operation. Show how this yields the same result as Approach A. (Hint: the correlation is a function of two time delays, 0 and 1)
In this problem you will estimate a filter response from non-white stimuli. We will use the same two filters from the last problem set:
For both filters, ms and ms.
1. Plot and for ms with a bin size of 1 ms.
2. Your starting point will be a stimulus which is a Gaussian white noise signal 4 s in length (hint: in numpy, use
numpy.random.randn). What is the bandwidth of this signal?
3. Unfortunately, your neuron responds only weakly to this broadband signal. You decide to filter the signal using to obtain a new signal . What is the bandwidth of this signal?
4. is used as a stimulus for a cell that has an impulse response given by . You can obtain “responses” for that cell, , by convolving with . What is the power spectrum of , ? What is the cross-correlation of and ? Hint: in the frequency domain, ; i.e., the product of the Fourier transform () of one signal times the complex conjugate () of the Fourier transform of the other signal.
You’ll need to do some averaging for this to work, as you did in the last problem set. What is the best chunk size and why (i.e., what happens if the chunk is too large or too small)?
5. Plot the cross-correlation in the time domain. Describe what you see and what you expected.
6. Use the equation you derived in Problem 1 (which should be the same as we discussed in class) to compute the proper impulse response of the filter. (Hint: do the normalization in the frequency domain and then return to the time domain for plotting). Describe what you see what what you expected. Did it work?
7. Now try filtering out the high frequencies in the cross-correlation between the stimulus and response before normalizing (hint: what’s filtering in the frequency domain?). Try upper frequency cutoffs of 10, 20, and 30 Hz. Explain your results.
A visual neuron has a temporal impulse response function given by in Problem Set 2, but with ms and ms.
1. Plot at a time resolution of 1 ms from -5 to 150 ms.
2. Generate 100 s of Gaussian white noise stimulus at a time resolution of 1 ms, and plot its power spectrum.
3. Assuming the cell is linear, predict its response (firing rate) to the white noise stimulus generated in 2.
4. Based on the prediction computed in 3, generate a spike train for one trial using a Poisson generator. To do this, you need to rectify the “firing rate” computed in 3 at 0 (make all the negative values 0, because the firing rate of the neuron cannot be negative). The probability of firing a spike at each moment is proportional to this rectified firing rate. Adjust the static nonlinearity so that the mean firing rate of this cell during the white-noise stimulation is 50 Hz.
5. Estimate the temporal receptive field of this cell using spike-triggered averaging. You’ll probably want to write a function to do this. Explore the effect of mean spike rate on your estimated kernel (e.g., 1Hz, 10Hz), assuming you can maintain stable recording from the cell for only 100 sec. Then try keeping the mean firing rate fixed and varying the recording time.
6. Take the kernel computed from one set of simulated data (i.e., 100 sec, 50 Hz mean firing rate) and convolve it with the stimulus – this is your linear prediction of the firing rate. Plot the firing rate of the cell, computed from the spike train, as a function of the linear prediction. Compare the predicted rate against the actual rate used to generate the data to recover the nonlinearity (you’ll need to do some binning). What is the correlation coefficient for your prediction?
A retinal ganglion cell has the following spatiotemporal receptive field (STRF) .
is the neuron’s spatial receptive field, which has an ON/OFF structure. This can be described as a difference of Gaussians, with the following parameters. Center: amplitude 1, standard deviation 1 pixel. Surround: amplitude -0.5, standard deviation 1.5 pixels. For the equation for a 2D Gaussian, see equation 2.45 in Dayan and Abbott - but replace the complicated scaling factors with the amplitudes provided here.
is the neurons’ temporal receptive field, and it’s defined as our good friend the double exponential.
ms and ms.
1. Plot the RF (21 by 21 pixels) at t = 0, 10, 30, 100, 300, and 600 msec. You can either use 3D plots or images, but keep the RF at different delays on the same scale.
2. Generate a random spatiotemporal stimulus with a Gaussian distribution and 0 mean (use
randn in python), 21 by 21 pixels at a frame rate of 100 Hz (one frame every 10 msec), which lasts for 5 min. Predict the linear response of the cell . Plot a 1-sec segment of the predicted .
3. Take the linear response of the cell, and threshold it at 0. (use
maximum in python). Plot r(t) over the same 1-sec period as in (b).
4. Use the whole 5 min of the responses to estimate the STRF of this neuron with the reverse correlation method. Plot the estimated RF as in (a). Compare it to the initial RF in terms of both the shape and the amplitude. What happens if you use rather than ?