Combining Probabilistic Population Codes

1 An Example

Consider the case of a hunter attempting to shoot a pheasant as it flies out of a tree. We'll assume that the hunter uses two cues, a visual cue concerning mo­ tion in the tree and an auditory cue based on rustling of the leaves, to estimate the pheasant's size and velocity. Based on this estimate, he selects a time and place to fire his shotgun.

The combination problem concerns how the two in­ puts should be combined to produce the output. In the simplest version of the combination problem for this ex­ ample, visual motion is confined to one part of the tree, and the auditory signal directly corresponds to this vi­ sual signal. Here these two single-valued inputs (which we will term v and a) give rise to a single output, and the hunter confidently aims his shotgun (to location s).

Evidence exists that the two inputs and the output in­ formation in this example are each represented in neural population codes in some animals. That is, a fixed col­ lection of neurons fire for each of the three variables of interest. The relevant visual input is represented by the activity of a population of motion detectors: in monkeys, a particular cortical area (MT) contains cells that selec­ tively respond to motion of a particular velocity within a small range of visual locations. Similarly, the relevant auditory input is represented in a population of detec­ tors tuned to particular frequencies and spatial locations in owl auditory cortex [Knudsen and Konishi, 1978]; the frequency may contain important information about the bird's size and speed. Directional motor output is also represented in a population code in monkey motor cortex [Georgopoulos et al. , 1986].

Therefore even in the simple version of the problem, the brain does not directly represent the values v, a, and s, but instead represents each in a separate popu­ lation code. The most straightforward way to solve this problem is to perform an intermediate step of extracting separate single values from the input population codes, combine these values, and then encode these into the motor output population code. However, this seems not to be the strategy actually implemented in the brain, where new population codes appear to be generated di­ rectly from old ones.

Another level of complexity is introduced into the problem when we consider that the inputs may be uncer­ tain or ambiguous. For example, if the wind is blowing, then leaves may be moving all over the tree giving rise to multiple plausible motion hypotheses, while at the same time the auditory cues may be too faint to confidently estimate the motion. The experienced hunter may then be able to narrow down the set of candidate motions based on his knowledge of the combinations of auditory and visual cues, but he might not be able to confidently select a single value. Two additional problems are intro­ duced in this more general case. First we must interpret a population code as representing a whole probability distribution over the underlying variable. And then the combination method must preserve the probabilistic in­ formation in the inputs. Thus the aim of a combination network is to infer a population code for the motor ac­ tion that preserves the statistical relationship between the input and output probability distributions.

2 Theory

The basic theory underlying the combination of popula­ tion codes is e:)Ctremely simple. Population codes use the explicit activities r = {r;} of multiple cells (as in area MT) to code information about the value of an implicit underlying variable x (such as the direction and speed of motion of the leaves). We are interested in the case that the activities r code a whole probability distribution over the underlying variable:

Consider the example of the hunter. Activities {ri}

and {rj} represent probability distributions over the mo­ tion position and velocity based on the visual and audi­ tory signals respectively. We will assume that afferent information in the different modalities is independent. The activities {rk} will represent a probability distribu­ tion over the corresponding required position s of the shotgun according to the equivalent of Equation 1.

Two computational operations are required to produce appropriate r 8 : information from the different modal­ ities must be integrated and then expressed in appro­ priate coordinates. These operations have to respect the statistical semantics of rv and ra. We use an un­ derlying analysis-by-synthesis statistical model as in the Helmholtz machine [Hinton et al., 1995]. In such a model, inference is based on the analysis or recognition inverse to a probabilistic synthetic or generative model that specifies probability distributions P[v, a!s] over the visual motion signal v and auditory pattern a given the shotgun location s.

Given true probability distributions P[v!w] and P[a!w] over the visual and auditory information (here w repre­ sents the underlying information available to the hunter) , recognition requires calculating:

where P[s] is the prior distribution over s. Equation 3 establishes the standard by which inferences about the distribution over s should be judged.

We have therefore reduced the computational prob­ lem to one of mapping activities rv and ra into activities r 8 for which P[s!r 8 ] from Equation 1 is a good approxi­ mation to the integration in Equation 3, where P[v!w] is what rv ,represents (according to Equation 1) and P[a!w] is what ra represents. Figure 1 illustrates the generative and recognition operations, showing the activities, the distributions that they represent, and the various prob­ abilistic relationships.

The remaining questions concern how activities r spec­ ify distributions as in Equation 1, and how rv and ra are actually combined to produce r 8 . We describe two models for Equation 1: a model based on a standard form of function approximation, kernel density estima­ tion (KDE) [Anderson, 1994], and an extension to the conventional statistical population coding model that is designed to handle any form of PDF [Zemel et al., 1997]. Both models form estimates pr (x) of P[x!w] based on r.

2.1 The KDE Model

One way of treating population codes as distributions is in terms of kernel density estimates (KDEs) [Anderson, 1994]. Here, activities r represent distribution pr (x)

through a linear combination of basis functions '1/Ji(x), ie Pr (x) = l:i r:'l/Ji(x) where the {r:} are normalized such that Pr (x) is a probability distribution. The kernel func­ tions '1/Ji(x) are not what are called tuning functions of the cells (the tuning function fi(x) describes the rela­ tionship between the cell's firing rate ri and the under­ lying variable x, where fi(x) could be Gaussian about some preferred value xi) · The kernel functions need have no neural instantiation; instead, they form part of the interpretive structure for the population code. If the '1/Ji(x) are probability distributions, and so are positive, then the range of spatial frequencies in P[xlw] that they can reproduce in Pr (x) is likely to be severely limited.

Evaluating the KDE model requires some method of representing P[xlw] in a firing rate vector r, which we term encoding. One way to encode is to use the Kullback­ Leibler divergence as a measure of the discrepancy be­ tween P[xlw] and l:i r:'l/Ji(x) and use the expectation­ maximization (EM) algorithm to fit the {r:}, treating them as mixing proportions in a mixture model [Demp­ ster et al., 1977]. This relies on { '1/Ji(x)} being probability distributions themselves. The projection method [Ander­ son, 1994] is a one-shot linear filtering based alternative using the £ 2 distance. Here the ri are computed as a projection of P[xlw] onto tuning functions fi(x) that are calculated from '1/Ji(x):

An extremely attractive property of the KDE model is that it makes combination of population codes very A rv simple [Anderson, 1994]. If P (v) = l:i rY'l/Ji(v) is an A ra adequate model of P[vlw] and P (a) = l:i rj'l/Jj(a) is an adequate model of P[alw], then

where

This makes Equation 3 into a simple bilinear sum, even when v and a are not independent given s.

Probabilistically correct combination is not so simple in this model if the underlying multiplication in Equa­ tion 5 is not allowed. In the combination network de­ scribed below, we restricted the combination to the stan­ dard mechanism of linear summation followed by a non­ linear activation function.

2.2 The Extended Poisson Model

Unfortunately, the KDE model has significant difficulties representing probability distributions that involve higher spatial frequencies than the kernel functions . This is a natural, but significant limitation of this method, since the situation of nearly complete certainty in v or a or s is particularly important. An alternative method has been suggested that is based on the same probabilistic analysis that underlies most standard applications of population coding [Zemel et al., 1997].

For these applications the starting point is usually the neurophysiological finding, largely ignored by the KDE model, that the relevant cells have unimodal tuning func­ tions fi(x) . Standard accounts also assume that fir­ ing rates vary, even for a fixed input, eg the Poisson model [Seung and Sompolinsky, 1993] for which

This is an encoding model, since it relates how x is coded in r. The form of Equation 1 is then specified (in an operation referred to as decoding) as the statistical inverse to this encoding model Pr (x) = P[xl {ri }] ex P[x] fL P[rilx], where P[x] is the prior probability dis­ tribution over x. In these terms, the KDE model is spec­ ified by its method of decoding-Equation 4 automati­ cally follows as the method of encoding. Note that a single value of x gives rise to an estimate Pr (x) that is a distribution over x due to the assumed variability.

Whereas the KDE model has problems representing peaked distributions, the conventional Poisson model has problems representing broad distributions, since Equa­ tion 6 assumes that there is a single underlying value of x. If the information provided to the hunter is somewhat unspecific with respect to the visual information, then it does not make sense to assume that r 11 are all deter­ mined (stochastically) according to some single particu­ lar value. In the extended Poisson model, the recorded activities rare allowed to depend on general P[xlw], hav­ ing Poisson statistics with mean:

This equation is identical to that for the KDE model {Equation 4), except that variability is built into the Poisson statistics, and decoding is now required to be the Bayesian inverse of encoding. Note that since ri de­ pends stochastically on P[xlw], the full Bayesian inverse will specify a distribution P[P[xlw]lr] over possible dis­ tributions. The distribution in Equation 1 can then be generated as the maximum likelihood (or rather maxi­ mum a posteriori with respect to a smoothness prior) distribution that P[P[xlw]ir] implies.

Decoding in this model may be performed by approxi­ mating P[xlw] as a piece-wise constant histogram which takes the value ¢i in (xi, Xj+I], and h(x) by a piece-wise constant histogram that take the values hi in (xi, xi+ll· The maximum a posteriori estimate for { ¢i} can be de­ rived by maximizing:

where E is the variance of an assumed smoothness prior. A form of EM may be used to maximize the likelihood L( { ~i} . By comparison with the linear decoding of the KDE method, the extended Poisson model offers a non­ linear way of combining a set of activities {ri} to give a probability distribution pr (x) over the underlying vari­ able x. Figure 2 describes the full extended Poisson model and illustrates the underlying probabilistic frame­ work for population codes.

Unfortunately, there are no short cuts like Equation 5 for combining these extended Poisson population codes. Instead, we have adopted some particular functional form for the combination and optimized its parameters in order to satisfy Equation 3 as best as possible.

3 Experiments

The central question for both coding methods is whether the combination network can form a set of activities {rU · r · that makes P (s) a close approximation to the implicit distribution P[siw]. For the purpose of evaluating net­ work performance, and adapting its weights, several as­ pects need to be determined: 1). the implicit distribu­ tion P[siw] . In the hunting example, this distribution describes the likely shotgun motions based on all infor­ mation available to the hunter, so multiple peaks corre­ spond to different possible motions and entropy corre­ sponds to uncertainty about these motions. 2). the gen­ erative model P[v,als]. The implicit distributions for the network inputs are produced by applying P[v, ais] to P[siw]. In these simulations, we made the simplifying assumption that the visual and auditory signals are in­ dependent givens. 3). the encoding model. The inputs r 11 and ra are obtained from the input implicit distri­ butions via appropriate encoding model (Equation 4 for KDE; Equation 7 for the extended Poisson method) . 4). a combination function. The network inputs produce an output r 8 based on a weighted combination of r 11 and ra . In these simulations we had both excitatory W and in­ hibitory weights U between each input and output unit, and the combination function was:

Note that this is not quite general enough to imple­ ment Equation 5 exactly.

We evaluate the networks' performances by comparing · r· the P (s) obtained by decoding the explicit representation s in the network to the true implicit distribution P[slw]. The error function is the Kullback-Leibler diver­ gence:

The network weights (excitatory and inhibitory) are adapted to minimize this error function using the delta­ rule.

We examined the combination network performance under several target distributions and generative mod­ els. Here we focus on a difficult, important target dis­ tribution, a bimodal Gaussian distribution 1/2N(2, r) + 1/2N( -2, r) that contains both uncertainty in a given position (r) and multiple possible values (x = {2, -2}). The generative model is simple: P[vls] and P[ais] are Gaussians with variance 0.5 and means at 1.0 and -1.0 respectively.

We compared the two methods of probabilistic in­ terpretation on this problem. For the kernel func­ tions (KDE) and the tuning functions (extended Pois­ son model), we used Gaussians N(xi , 0.3). To reduce the number of units in the network, we considered a simplified situation where all signals are constrained to one dimension. This reduction is only a matter of conve­ nience; the network can readily be extended to include higher-dimensional implicit spaces. We used 50 units in each population code, and the Xi were spaced evenly in the range x = [-10, 10] .

Results are shown in Figures 3 and 4. At low values of r, the target distribution approximates two spikes, at x = 2 and x = -2. The KDE method is not as accurate here, since it is unable to retain the high frequency in­ formation required for precise recovery of the these two location. As the target distribution contains more un­ certainty, both methods are able to recover the implicit distribution with high fidelity. Note that an error rate of 0. 7 bits for this target distribution would be obtained ' r · if P (s) has the correct peaks and is off by a factor of 2 in T (see Equation 9).

We have also conducted a number of other experi­ ments with this combination method. In one set of ex­ periments, we modeled the task of combining monocular and stereo cues to estimate depth in a particular visual il­ lusion. In the double-nail illusion, the task is to estimate the depth of a nail aligned directly behind another nail in the observer's line of sight. Here computational vision systems based on binocular stereo produce a PDF for depth estimates with two peaks, one at the correct value and another at the illusory frontoparallel interpretation (both nails side-by-side). A PDF based on monocular cues will not have the same ambiguity, but it is typically a much broader distribution [Yuille and Biilthoff, 1994]. These two PDFs must be combined multiplicatively to produce the correct peak.

We simulated this problem by training a combina­ tion network identical to the network described above except in the generative model. Here P[blt] is a mul­ timodal Gaussian 1/3N[t, 1/2] + 2/3N[t + 2, 1/2] (with a frontoparallel bias) and P[mlt] is a broader unimodal Gaussian N[t, 1] , where b, m and t are the binocular, monocular and true depth estimates, respectively. Af­ ter training on 300 cases in which the target distribution was a narrow Gaussian N[t, .01], the network produced output distributions on novel inputs that were within .1 bits of the true distributions.

Other experiments have examined the combination network's ability to recover PDFs in which the certainty as to the presence of the output ( ie the integral under the PDF) is < 1. Good performance on this task sug­ gests that the method can be useful for recognition ( eg recognizing an instance of an object based on the spatial locations of its features).

4 Discussion

We have presented a general framework for mapping be­ tween population codes that approximates statistically . correct inference. The framework applies and extends two recent methods for the probabilistic interpretation of population codes to the problem of combining these codes. This framework has a wide variety of applica­ tions, including any context in which probabilistic in­ formation from several sources, each represented in a distributed manner, must be combined. The simulation results demonstrate that a feedforward network can cap­ ture the appropriate probabilistic relationships between some simple population-coded PDFs. Generally, several population-coded inputs should be multiplied (to com­ pute a full joint PDF), but we found empirically that they can be combined reasonably using a non-linearity.

A straightforward alternative to the proposed frame­ work would extract single values from the input pop­ ulation codes, combine these values, and then form a new population code at the output. Aside from biolog­ ical realism, the computational advantage of construct­ ing direct mappings between population codes without requiring an intermediate step of extracting single val­ ues is that information about whole distributions can be brought to bear-including the ambiguity and uncer­ tainty in the underlying variables.

Integral to the framework is an interpretation of a pop­ ulation code as encoding a probability distribution over the underlying quantity. The framework can thus be seen as a generalization of [Salinas and Abbott, 1995], in which a network is trained to map one population code to another, where each code is interpreted as represent­ ing a single value. Our method extends this mapping to probabilistic interpretations while maintaining the bio­ logically realistic representations.

There are many open issues, .Particularly understand­ ing the nature of encoding and decoding. Both op­ erations are only implicit in the system so some free­ dom exists in choosing ones appropriate for particular tasks. Based on neurobiological and engineering consid­ erations, one expects a consistent interpretation across levels; maintaining this interpretation should lead to a simple learning rule. Noise is a second key issue. H constructing one population code from others introduces substantial extra noise, the system will be unable to con­ vey information accurately. Here the restriction of the network to feedforward connections might be relaxed in order to allow lateral connections between units within a population, which may be useful in cleaning up the codes.