Aug '07

Multiply Tuned Neurons

by O'Doherty in Multi-Tuned Neuron Project

Introduction

The goal of this project is to determine what tuning curves look like for neurons that are simultaneously tuned to multiple stimuli. We are using pre-existing data in which monkeys had to focus on two different visual cues at the same time. The neurons recorded in this study are in the pre-frontal cortex and were taken from two different monkeys. The data was acquired by Misha Lebedev at NIH; the paper summarizing that work can be read here. The notes Iyad took during his summer 2007 visit to see Joey and Misha can be read here.

Details of Misha’s Experiment

A single trial of Misha’s experiment proceeds as follows:

Subject presses a button to start the trial
Fixation point appears at center of the screen
Subject fixates for 1-1.5s
Target appears in one of four locations (0°, 90°, 180°, 270°)
After another 1-1.5s, the target revolves around to one of the four locations (could be the same one it started at)
After the target stops, there is a 1-2.5s delay
The target gets either brighter or dimmer for 150ms and then extinguishes. This serves as both the trigger for the subject to saccade, as well as the instruction for which target to saccade to.
Subject saccades to appropriate target; correct responses are rewarded with juice.

Data File Structure

The data come from two monkeys, each of which has had multiple neurons tested. Each neuron has been tested with multiple trials. The list of neurons are stored in the files PF_all_zach.lst and PF_all_hamp_noPMd.lst. A sample line from Zach’s file reads:

za2_1 8 a

The “za2_1″ identifies the trial day number and the electrode number; these are irrelevant details for the purpose of this experiment. The “8″ refers to the unit number and “a” refers to the subunit. There are many other units and subunits for each electrode, but only those listed in these files have been identified by Misha et al as being valid for study. The data for the corresponding electrode would be stored in za2_1_code.mat. That file would contain the results of all the different trials for that electrode, only some of which involve unit 8 and subunit a.

The data file za2_1_code.mat contains data for each of a large number of trials of the experiment. For example, the variable “spike_time” is a cell array of length 1043, meaning that each entry in the cell array is an array with the firing times of all the neurons during one trial. We can determine which channel yielded each spike time by cross-referencing with “spike_channel”. The active channels and subchannels of the “nth” trial are saved in channels{n} and subchannels{n}.

All times must be centered relative to the target blink time by executing the following line: cntr = target_blink;

The code can then march through the trials and determine how many times the target subunit fired during a given time period under different conditions (such as specific combinations of attended and remembered locations).

Determining Qualifying Neurons

According to Misha’s paper, only certain neurons qualify for study – those which are concurrently tuned to both attended and remembered target locations. This is determined by computing two custom statistics: I_Rem and I_Att; if both of these are significantly greater than one, the neuron is said to be concurrently tuned to both stimuli. The formulas for these statistics are given in equations (1) and (2) of Misha’s paper (page 1932). For example, if I_Rem equals unity for a given neuron, then there is little difference between the variance of firing rates within rows (denominator) versus between rows (numerators); such a neuron is not tuned to the remembered target location. If the within-row variance is much smaller than the between-row variance, I_Rem will be greater than one and the neuron will be said to be tuned. In order to test statistical significance of I_Rem for a particular neuron, we repeat the following procedure 1000 times: Reassign the firing rates randomly between trials of with the same Attended location and recompute I_Rem; in theory, this should result in an I_Rem of unity. If the original value of I_Rem is greater than 99% of the randomized I_Rem values, we can be fairly sure that it is significantly greater than one, and that neuron is said to be tuned for remembered target location. The same procedure can be applied for measuring and testing the significance of I_Att. If done correctly, it should be possible to re-create Figure 3D of Misha’s paper (page 1923); the randomization aspect may limit how precisely that figure can be reproduced.

The following table summarizes Iyad’s attempt to recreate Misha’s statistics. Misha’s numbers are taken from page 1923 and Table 2. The differences between the two columns are likely due to the fact that the determination of statistical significance requires randomization, which will have been different between Misha’s and Iyad’s programs. However it is a little surprising that there are as many differences as are observed – it would seem that 1000 repetitions would be enough to virtually eliminate the effects of randomization.

	Misha	Iyad
IAtt	1.84 ± 0.08	1.84 ± 0.08
IRem	1.21 ± 0.02	1.22 ± 0.02

n Att Tuned	186	179
n Rem Tuned	47	53
n Both Tuned	70	77
Total	303	309

Comparing the Fits

The 77 dual-tuned neurons were fitted to three different firing rate models:

1. The SUM of two independent preferred-direction cosine tuned curves

2. The PRODUCT of two independent preferred-direction cosine tuned curves

3. A linear fit of the difference between the two angles

The following procedure was repeated 100 times for each neuron: Use 75% of the trials to fit the model parameters. Then test the fit by using the remaining 25% of the trials to compare predicted versus actual firing rate. The results show very clearly that, for any given neuron, all three models fit the data equally well:

The code producing these traces was backed up on 12/07/2007.

Upon further review, it appears that these results are predicted by theory. This was determined by the following experiment. I created a neuron whose firing rate was the sum of the firing rates predicted cosine tuning curves for two independent inputs plus some zero-mean Gaussian random noise. A series of randomized trials were generated for this neuron. The resulting firing rates were reverse fitted to a different model, one where the firing rate is the product of the cosine tuning curves. Finally, the percent error between the “actual” firing rate (from the additive model) and the “predicted” firing rate from the fitted multiplicative model was computed. It was found that the mean percent error was roughly zero, with neither mean nor standard deviation changing dramatically as the amount of Gaussian noise summed into the additive neuron model was increased. This verifies that both the additive and multiplicative neuron models are flexible enough to fit the data reasonably well. The code generating this trial is found here.

A more interesting question is whether the competing models predict different preferred directions for attention and memory. This was tested as follows: for each of the 77 dual-tuned neurons, the two preferred directions were computed twice, once according to the “summation” model and once according to the “multiplication” model. For each neuron, the preferred directions were recalculated 100 different times and these were averaged respectively; each iteration used different trials for training versus testing the model coefficients. The results are shown below. Each plot has “preferred direction according to summation model” on the x-axis and “preferred direction according to multiplication model” on the y-axis. The left plot is predicted preferred directions for “attention” whereas the right plot is for “memory”. As is evident, they are linearly related with slope of 0.6 in both cases. The code used to produce this plot is backed up on 12/10/2007.

Background on Tuning Functions

Cosine tuning describes how the firing rate of a given neuron varies as a function of the direction of a particular motor task. It was first described by Georgopoulos in 1982 for 2-D movement and again in 1986 for 3-D movement.

Single Variable Tuning

The basic equation is given as follows:

$FR = a + b*cos(\theta - \theta_p)$

where a represents the mean firing rate and b represents the maximum deviation from a. The maximum firing rate is therefore a + b and the minimum firing rate is a-b. The cosine tuning equation can equivalently be represented as a function of the sine and cosine of the direction of movement:

$FR = a + \beta_1sin(\theta) + \beta_2cos(\theta)$

where

$b = \sqrt{\beta_1^2 + \beta_2^2}$

and

$\theta_p=tan^{-1}\left(\frac{\beta_1}{\beta_2}\right)$

The concept of a cosine tuned neuron with a preferred direction is therefore equivalent to a linear regression of firing rate against the sine and cosine of the movement direction.

Multiply Tuned Curves

Based on this original concept, there are a number of models to consider for multiply tuned neurons. In the following lines, I will assume that there are two tuning variables, $\theta_{att}$ and $\theta_{mem}$ :

Additive

$FR = \left(a_{att}+\beta_{1,att}sin(\theta_{att})+\beta_{2,att}cos(\theta_{att})\right)+\left(a_{mem}+\beta_{1,mem}sin(\theta_{mem})+\beta_{2,mem}cos(\theta_{mem})\right)$

which is equivalent to:

$FR=c_1 + c_2sin(\theta_{att}) + c_3cos(\theta_{att}) + c_4sin(\theta_{mem}) + c_5cos(\theta_{mem})$

The Matlab code developed for this project computes the coefficients $c_1 - c_5$ using a least squares approximation. From those coefficients, we can work backwards to compute the estimated preferred direction for the two tuning variables as follows:

$\hat{\theta}_{att,preferred}=tan^{-1}\left(\frac{c2}{c3}\right)$

$\hat{\theta}_{mem,preferred}=tan^{-1}\left(\frac{c5}{c6}\right)$

Multiplicative

$FR = \left(a_{att}+\beta_{1,att}sin(\theta_{att})+\beta_{2,att}cos(\theta_{att})\right)\times\left(a_{mem}+\beta_{1,mem}sin(\theta_{mem})+\beta_{2,mem}cos(\theta_{mem})\right)$

which is equivalent to:

$FR=c_1 + c_2sin(\theta_{att}) + c_3cos(\theta_{att}) + c_4sin(\theta_{mem}) + c_5cos(\theta{mem})$

$+ c_6sin(\theta_{att})sin(\theta_{mem}) + c_7sin(\theta{att})cos(\theta_{mem}) + c_8cos(\theta_{att})sin(\theta_{mem}) + c_9cos(\theta_{att})cos(\theta_{mem})$

The Matlab code developed for this project computes the coefficients $c_1 - c_9$ using a least squares approximation. Working backwards from those coefficients to compute the estimated preferred directions requires solving a system of non-linear equations. There is no closed-form expression for this. Rather, a Levenberg-Marquardt approximation can be executed in Matlab to estimate the beta values. From there, we can work backwards to compute the estimated preferred directions:

$\hat{\theta}_{att,preferred}=tan^{-1}\left(\frac{\hat{\beta}_{1,att}}{\hat{\beta}_{2,att}}\right)$

$\hat{\theta}_{mem,preferred}=tan^{-1}\left(\frac{\hat{\beta}_{1,mem}}{\hat{\beta}_{1,mem}}\right)$

I attempted to evaluate how reliable this inverse method is for estimating the preferred directions from the estimated beta values. I created a neuron that was tuned according to the product of two cosine tuned parameters. I added noise with standard deviation σ_n to both of the cosine tunes before multiplying them to get the estimated final firing rate. This represented one trial. The beta values were computed from a set of 100 trials and the preferred directions were computed as stated above. This procedure was again repeated 100 times. The histograms of the preferred direction prediction errors are shown below: left column is attenuation, right column is memory, top row is σ_n = 0.1, bottom row is σ_n = 0.3. As is evident from the plots, the predicted preferred direction is occasionally 180 degrees off, which biases the error measures. The code used to generate the figure below can be found here.