Information-geometric decomposition in spike analysis - CiteSeerX
Information-geometric decomposition in spike analysis Hiroyuki Nakahara, Shun-ichi Amari
Lab. for Mathematical Neuroscience, RIKEN Brain Science Institute 2-1 Hirosawa, Wako, Saitama, 351-0198 Japan fhiro, [email protected]
Abstract We present an information-geometric measure to systematically investigate neuronal ring patterns, taking account not only of the second-order but also of higher-order interactions. We begin with the case of two neurons for illustration and show how to test whether or not any pairwise correlation in one period is signi cantly di erent from that in the other period. In order to test such a hypothesis of di erent ring rates, the correlation term needs to be singled out `orthogonally' to the ring rates, where the null hypothesis might not be of independent ring. This method is also shown to directly associate neural ring with behavior via their mutual information, which is decomposed into two types of information, conveyed by mean ring rate and coincident ring, respectively. Then, we show that these results, using the `orthogonal' decomposition, are naturally extended to the case of three neurons and n neurons in general.
1 Introduction Based on the theory of hierarchical structure and related invariant decomposition of interactions by information geometry [3], the present paper brie y summarizes methods useful for systematically analyzing a population of neural ring [9]. Many researches have shown that the mean ring rate of a single neuron may carry signi cant information on sensory and motion signals. Information conveyed by populational ring, however, may not be only an accumulation of mean ring rates. Other statistical structure, e.g., coincident ring [13, 14], may also carry behavioral information. One obvious step to investigate this issue is to single out a contribution by coincident ring between two neurons, i.e., the pairwise correlation [2, 6]. In general, however, it is not sucient to test a pairwise correlation of neural ring, because there can be triplewise and higher correlations. For example, three variables (neurons) are not independent in general even when they are pairwise independent. We need to establish a systematic method of analysis, including these higher-order
also aliated with Dept. of Knowledge Sci., Japan Advanced Inst. of Sci. & Tech.
correlations [1, 5, 7, 13]. We propose one approach, the information-geometric measure that uses the dual orthogonality of the natural and expectation parameters in exponential family distributions [4]. We represent a neural ring pattern by a binary random vector x. The probability distribution of ring patterns can be expanded by a log linear model, where the set fp(x)g of all the probability distributions forms a (2n ? 1)-dimensional manifold S n . Each p(x) is given by 2n probabilities X pi1 i = 1 pi1 i = Prob fX1 = i1 ; ; Xn = in g ; ik = 0; 1; subject to n
n
i1 ; ;in
and expansion in log p(x) is given by X X log p(x) = i xi + ij xi xj + i<j
X i<j 0; i; j = 0; 1. Among four probabilities, fp00; p01 ; p10 ; p11 g, only three are free. The set of all such distributions of x forms a three-dimensional manifold S 2 . Any three of pij can be used as a coordinate system of S 2 . There are many di erent coordinate systems of S 2 . The coordinates of the expectation parameters, called -coordinates, = (1 ; 2 ; 12 ), is given by i = Prob fxi = 1g = E [xi ] ; i = 1; 2; 3 = 12 = E [x1 x2 ] = p12 ; where E denotes the expectation and i and 12 correspond to the mean ring rates and the mean coincident ring, respectively. As other coordinate systems, we can also use the triplet, (1 ; 2 ; Cov [X1 ; X2 ]), where Cov [X1 ; X2 ] is the covariance,and/or the triplet (1 ; 2 ; ), where is the correlation coecient (COR), = p1 (1?12?1)12(12 ?2 ) , often called N-JPSTH [2].
Which quantity would be convenient to represent the pairwise correlational component? It is desirable to de ne the degree of the correlation independently from the marginals (1 ; 2 ). To this end, we use the `orthogonal' coordinates (1 ; 2 ; ), originating from information geometry of S 2 , so that the coordinate curve of is always orthogonal to those of 1 and 2 . The orthogonality of two directions in S 2 (S n in general) is de ned by the Riemannian metric due to the Fisher information matrix [8, 4]. Denoting any coordinates in S n by = (1 ; :::; n ), the Fisher information matrix G is given by G
= (gij ) ;
gij ( )
=E
@ @i
l (x; )
@ @j
l (x; ) :
(1)
where l (x; ) = log p (x; ). The orthogonality between i and j is de nediby h @ @ gij ( ) = 0. In case of S 2 , we desire to have E @ l (x; 1 ; 2 ; ) @ l (x; 1 ; 2 ; ) = 0 (i = 1; 2): When is orthogonal to (1 ; 2 ), we say that represents pure correlations independently of marginals. Such is given by the following theorem. Theorem 1. The coordinate p11 p00 (2) = log p p i
01 10
is orthogonal to the marginals 1 and 2 . P We have another interpretation of . Let's expand p(x) by log p(x) = p 2i=1 i xi + 12 x1 x2 ? : Simple calculation lets us get the coecients, 1 = log p10 ; 2 = 00 log pp0001 ; = ? log p00 ; and = 12 (as Eq 2). The triplet = (1 ; 2 ; 12 ) forms another coordinate system, called the natural parameters, or ?coordinates. We remark that 12 is 0 when and only when X1 and X2 are independent. The triplet (1 ; 2 ; 12 ) forms an `orthogonal' coordinate system of S 2 , called the mixed coordinates [4]. We use the Kullback-Leibler divergence (KL) to measurePthe discrepancy between two probabilities p(x) and q(x), de ned by D [p : q] = x p(x) log pq((xx)) . In the following, we denote any coordinates of p by p etc (the same for q). Using the orthogonality between - and -coordinates, we have the decomposition in the KL.
Theorem 2.
= D [p : r ] + D [r : q] ; D [q : p] = D [q : r ] + D [r : p] ; (3) q q p p p q r r where r and r are given by = (1 ; 2 ; 3 ) and = (1 ; 2 ; 3 ), respectively. The squared distance ds2 between two nearby distributions p(x; ) and p(x; ; +d) is given by the quadratic form of d, X gij ( )di dj ; ds2 = D [p : q ]
2
i;j (1;2;3)
which is approximately twice the KL, i.e., ds2 2D [p(x; ) : p(x; + d)] : Now suppose is the . Then, the Fisher information matrix 3 2 mixed coordinates is of the form g = 4 ij
g33 (d3 )2 ; ds22
=
P
g11 g12 0 g12 g22 0 0 0 g33
5
i;j (1;2) gij di dj ;
2
and we have ds2 = ds21 + ds22 ; where ds21 =
corresponding to Eq. 3.
This decomposition comes from the choice of the orthogonal coordinates and gives us the merits of simple procedure in statistical inference. First, let us estimate the parameter = (1 ; 2 ) and from N observed data x1 ; :::; xN . The maximum likelihood estimator (mle) ^, which is asymptotically unbiased and ecient, ?^2 +^12 ) , using is easily obtained by ^i = N1 #fxi = 1g and ^ = log ^(^121(1??^12^1)(^ 2 ?^12 ) 1 #fx x = 1g. The covariance of estimation error, and , is given ^12 = N 1 2 1 ?1 ?1 asymptotically by Cov = N G : Since the cross terms of G or G vanish for the orthogonal coordinates, we have Cov [; ] = 0, implying that the estimation error of marginals and that of interaction are mutually independent. Such a property does not hold for other non-orthogonal parameterization such as the COR , the covariance etc. Second, in practice, we often like to compare many spike distributions, q(x(tp)) (i.e, q(t) ) for (t = 1; ; ; T ), with a distribution in the control period p(x), or . Because the orthogonality between and allows us to treat them independently, these comparisons become very simple. These properties bring a simple procedure of testing hypothesis concerning the null hypothesis H0 : = 0 against H1 : 6= 0 ; (4) where 0 is not necessarily zero, whereas 0 = 0 corresponds to the null hypothesis of independent ring, which is often used in literature in di erent setting. Let the log likelihood of the models H0 and H1 be, respectively, l0 = max log p(x1 ; :::; xN ; ; 0 ) and l1 = max log p(x1 ; :::; xN ; ; ): ; The likelihood ratio test uses the test statistics = 2 log ll10 . By the mle with respect to and , which can be performed independently, we have ^ ; 0 ); ^ ; ^12 ); l0 = log p( x; l1 = log p( x; (5) where ^ are the same in both models. A similar situation holds in the case of testing = 0 against 6= 0 for unknown . Under the hypothesis H0 , is approximated for a large N as ^ log p(xi ; ^ ; ^0 ) N g33 ( ? 0 )2 2 (1): (6) ^ ; ) p(xi ; i=1 Thus, we can easily submit our data to a hypothetical testing of signi cant coincident ring against null hypothesis of any correlated ring, independently from the mean ring rate modulation1 . We now turn to relate the above approach with another important issue, which is to relate such a coincident ring with behavior. Let us denote by Y a variable of discrete behavioral choices. The MI between X = (X1 ; X2 ) and Y is written by p(x; y ) = Ep(Y ) [D [p(X jy) : p(X )]] : I (X; Y ) = Ep(X;Y ) log p(x)p(y )
=2
N X
j
Using the mixed coordinates for p(X y) and 0p(X ), we 0 D [ (X y ) : (X )] = D (X y ) : + D : (X ) ; (1 (X y); 2 (X y); 3 (X )) = (1 (X y); 2 (X y); 3(X )):
j
j
j
j
j
j
have D [p(X jy) : p(X )] = where 0 = 0 (X; y) =
1 A more proper formulation in this hypothetical testing can be derived, resulting in using p value from 2 (2) distribution, but we omit it here due to the limited space [9]
Theorem 3.
I (X; Y ) = I1 (X; Y ) + I2 (X; Y ); (7) where I1 (X; Y ); I2 (X; Y ) are given by 0 0 I1 (X; Y ) = Ep(Y ) D (X jy ) : (X; y ) ; I2 (X; Y ) = Ep(Y ) D (X; y ) : (X ) :
Obviously, the similar result holds with respect to p(Y jX ). By this theorem, I is the sum of the two terms: I1 is by modulation of the correlation components of X , while I2 is by modulation of the marginals of X . This observation helps us investigate the behavioral signi cance by modulating either coincident ring or mean ring rates. 0.1
η w/ s1
A
0.05
(a)
0 0
(b)
100
η1
0.05
(c)
η 2 η12
300
(b)
2
(c)
η1 η
12
100
300
500
3 C
D
θ3
2
0.1
1 0 0
100
300
0 0
500
1
0.95
E
MI (bit)
p−value
η (a)
0 0
500
0.2 COR
B
η w/ s2
0.1
0.5 0 0
100
300 time (ms)
500
0.04
100
300
500 I
F
I2
0.02
I
1
0 0
100
300 time (ms)
500
Figure 1: Demonstration of information-geometric measure in two neuron case, using simulated neural data, where two behavioral choices (s1, s2) are assumed. A,B. (1 ; 2 ; 12 ) with respect to s1, s2. C,D. COR,, computed by using = P i p(si ) (si ) with p(si ) = 1=2 (i = 1; 2). E. p-values. F. MI. Fig 1 succinctly demonstrates results in this section. Figs 1 A, B are supposed to show mean ring rates of two neurons and mean coincident ring for two di erent stimuli (s1, s2). The period (a) is assumed as the control period, i.e., where no stimuli is shown yet, whereas the stimulus is shown in the periods (b,c). Fig 1 C, D gives COR, . They look to change similarly over periods, which is reasonable because both COR and represent the same correlational component, but indeed change slightly di erently over periods (e.g., the relative magnitudes between the periods (a) and (c) are di erent for COR and ), which is also reasonable because both represent the correlational component as in di erent coordinate systems. Using in Fig 1 D, Fig 1 E shows p-values derived from 2 (1) (i.e., p > 0:95 in Fig 1 E is `a signi cance with p < 0:05') for two di erent null hypotheses, one of the averaged ring in the control period (by solid line) and the other of independent ring (by dashed line), which is of popular use in literature. In general, it becomes complicated to test the former hypothesis, using COR. This is because the COR, as the coordinate component, is not orthogonal to the mean ring rates so that estimation errors among the COR and mean ring rates are entangled and that the proper metric among them is rather dicult to compute. Once using , this testing becomes simple due to orthogonality between and mean ring rates. Notably, we would draw completely di erent conclusions on signi cant coincident ring given each null hypothesis in Fig 1 E. This di erence may be striking when we are to understand the brain function with these kinds of data. Fig 1 F shows the MI
between ring and behavior, where behavioral event is with respect to stimuli, and its decomposition. There is no behavioral information conveyed by the modulation of coincident ring in the period (b) (i.e., I1 = 0 in the period (b)). The increase in the total MI (i.e., I ) in the period (c), compared with the period (b), is due not to the MI in mean ring (I2 ) but to the MI correlation (I1 ). Thus, with a great ease, we can directly inspect a function of neural correlation component in relation to behavior.
3 Three neuron case With more than two neurons, we need to look not only into a pairwise interaction but also into higher-order interactions. Our results in the two neuron case are naturally extended to n neuron case and here, we focus on three neuron case for illustration. For three neurons X = (X1 ; X2 ; X3 ), we let p(x), x = (x1 ; x2 ; x3 ), be their joint probability distribution and put pijk = Prob fx1 = i; x2 = j; x3 = kg, i; j;P k = 0; 1. The set of all such distributions forms a 7-dimensional manifold S 3 due to pijk = 1. The -coordinates = (1 ; 2 ; 3 ) = (1 ; 2 ; 3 ; 12 ; 23 ; 13 ; 123 ) is de ned by i = E [xi ] (i = 1; 2; 3); ij = E [xi xj ] (i; j = 1; 2; 3; i 6= j ); 123 = E [x1 x2 x3 ] : To single out the purely triplewise correlation, we utilize the dual Porthogonality of - and -coordinates. By using expansion of log p(x) = i xi + P ij xi xj + 123 x1 x2 x3 ? , we obtain -coordinates, = (1 ; 2 ; 3 ) = (1 ; 2 ; 3 ; 12 ; 23p ; 13p ; p123 ).p It's easy to get the expression of these coecients 111 100 010 001 (e.g., 123 = log p110 p101 p011 p000 ). Information geometry gives the following theorem. Theorem 4. 123 represents the pure triplewise interaction in the sense that it is orthogonal to any changes in the single and pairwise marginals, i.e., 1 and 2 . We use the following two mixed coordinates to utilize the dual orthogonality, 1 = ( 1 ; 2 ; 3 ); 2 = ( 1 ; 2 ; 3 ): Here 2 is useful to single out the triplewise interaction (3 = 123 ), while 1 is to single out the pairwise and triplewise interactions together (2 ; 3 ). Note that 123 is not orthogonal to fij g. In other words, except the case of no triplewise interaction (123 = 0), ij do not directly represent the pairwise correlation of two random variables Xi ; Xj . The case of independent ring is given by ij = i j ; 123 = 1 2 3 or equivalently by 2 = 0; 3 = 0. The decomposition in the KL is now given as follows.
Theorem 5. D [p : q ]
= D [p : p] + D [p : q] = D [p : p~] + D [~p : q] = D [p : p] + D [p : p~] + D [~p : q] : (8) where, using the mixed coordinates, we have p2 = ( p1 ; p2 ; q3 ); p1~ = (p1 ; q2 ; q3 ). A hypothetical testing is formulated similarly to the two neuron case. We can exam p ine a signi cance of the triplewise interaction by 2 = 2N D [p : p] N g77 ( p2 )(123 ? q 2 2 123 ) (1): For a signi cance of triplewise and pairwise interactions together, P we have 1 = 2N D [p : p~] N 7i;j=4 gij ( p1 )(ip ? ip~)(jp ? jp~) 2 (4): For the decomposition of the MI between ring X and behavior Y , we have
Theorem 6.
I (X; Y )
= I1 (X; Y ) + I2 (X; Y ) = I3 (X; Y ) + I4 (X; Y )
(9)
where I1 (X; Y ) = Ep(Y ) [D [ 1 (X jy ) : 1 (X; y )]] ; I2 (X; Y ) = Ep(Y ) [D [ 1 (X; y ) : 1 (X )]] ; I3 (X; Y ) = Ep(Y ) [D [ 2 (X jy ) : 2 (X; y )]] ; I4 (X; Y ) = Ep(Y ) [D [ 2 (X; y ) : 2 (X )]] : By the rst equality, I is decomposed into two parts: I1 is conveyed by the pairwise and triplewise interactions of ring, and I2 by the mean ring rate modulation. By the second equality, I is decomposed di erently: I3 , conveyed by the triplewise interaction, and I4 , by the other terms. A
η
(b)
0.02
(c)
ηij η
(d)
100
300
500
700
0 −0.1 0 4
C
0
−2 0
0.1
θ123
θ12, 13, 23
i
ijk
0 0 2
0.2 B
η (a)
COR
0.04
100
300
500
1
500
700
D
0
700
−2 0
0.95
1
100
300
500
700 0.95
F χ2(4)
χ2(1)
300
2
E 0.5 0 0
100
0.5
100
300 time (ms)
500
0 700
0 0
100
300 time (ms)
500
0 700
Figure 2: Demonstration in three neuron case. A = (1 ; 2 ; 3 ) (i ; ij ; ijk ) from top to bottom, since we treated a homogeneous case in this simulation for simplicity. B. COR. C. 12 ; 13 ; 23 . D 123 . E p-value 2 (1). F p-value 2 (4). We emphasize that all the above decompositions come from the choice of the `orthogonal' coordinates. Fig 2 highlights some of the results in this section. Fig 2 A shows the mean ring rates (see legend). The period (a) is assumed as the control period. Fig 2 B indicates that COR changes only in the periods (c,d), while Fig 2 C indicates that 123 changes only in the period (d). Taken together, we observe that the triplewise correlation 123 can be modulated independently from COR. Fig 2 E indicates the p-value from 2 (1) against the null hypothesis of the activity in the control period. The triplewise coincident ring becomes signi cant only in the period (d). Fig 2 F indicates the p-value from 2 (4). The coincident ring, taking the triplewise and pairwise interaction together, becomes signi cant in both periods (c,d). We cannot observe these di erences in modulation of pairwise and triplewise interactions over periods (c, d), when we inspect only COR. Remark: For a general n neuron case, we can use the k-cut mixed coordinates, k = ( 1 ; :::; k ; k+1 ; :::; n ) = ( k? ; k+ ). Using the orthogonality between k? and k+ , the similar results hold. To meet the computational complexity involved in this general case, some practical diculties should be resolved in practice [9].
4 Discussions We presented the information-geometric measures to analyze spike ring patterns, using two and three neuron cases for illustration. The choice of `orthogonal' coordinates provides us with a simple, transparent and systematic procedure to test signi cant ring patterns and to directly relate such a pattern with behavior. We hope that this method simpli es and strengthens experimental data analysis.
Acknowledgments HN thanks M. Tatsuno, K. Siu and K. Kobayashi for their assistance. HN is supported by Grants-in-Aid 13210154 from the Ministry of Edu. Japan.
References
[1] M. Abeles, H. Bergman, E. Margalit, and E. Vaadia. Spatiotemporal ring patterns in the frontal cortex of behaving monkeys. J Neurophysiol, 70(4):1629{ 38., 1993. [2] A. M. H. J. Aertsen, G. L. Gerstein, M. K. Habib, and G. Palm. Dynamics of neuronal ring correlation: Modulation of "e ective connectivity". Journal of Neurophysiology, 61(5):900{917, May 1989. [3] S. Amari. Information geometry on hierarchical decomposition of stochastic interactions. IEEE Transaction on Information Theory, pages 1701{1711, 2001. [4] S. Amari and H. Nagaoka. Methods of Information Geometry. AMS and Oxford University Press, 2000. [5] S. Grun. Unitary joint-events in multiple-neuron spiking activity: detection, signi cance, and interpretation. Verlag Harri Deutsch, Reihe Physik, Band 60. Thun, Frankfurt/Main, 1996. [6] H. Ito and S. Tsuji. Model dependence in quanti cation of spike interdependence by joint peri-stimulus time histogram. Neural Computation, 12:195{217, 2000. [7] L. Martignon, G. Deco, K. Laskey, M. Diamond, W. A. Freiwald, and E. Vaadia. Neural coding: Higher-order temporal patterns in the neurostatistics of cell assemblies. Neural Computation, 12(11):2621{2653, 2000. [8] H. Nagaoka and S. Amari. Di erential geometry of smooth families of probability distributions. Technical report, University of Tokyo, 1982. [9] H. Nakahara and S. Amari. Information geometric measure for neural spikes. in prepration. [10] H. Nakahara, S. Amari, M. Tatsuno, S. Kang, K. Kobayashi, K. Anderson, E. Miller, and T. Poggio. Information geometric measures for spike ring. Society for Neuroscience Abstracts, 27:821.46 (page.2178), 2001. [11] M. W. Oram, N. G. Hatsopoulos, B. J. Richmond, and J. P. Donoghue. Excess synchrony in motor cortical neurons provides redundant direction information with that from coarse temporal measures. J Neurophysiol., 86(4):1700{1716, 2001. [12] S. Panzeri and S. R. Schultz. A uni ed approach to the study of temporal, correlational, and rate coding. Neural Computation, 13(6):1311{49., 2001a. [13] A. Riehle, S. Grun, M. Diesmann, and A. Aertsen. Spike synchronization and rate modulation di erentially involved in motor cortical function. Science, 278:1950{1953, 12 Dec 1997. [14] E. Vaadia, I. Haalman, M. Abeles, H. Bergman, Y. Prut, H. Slovin, and A. Aertsen. Dynamics of neuronal interactions in monkey cortex in relation to behavioural events. Nature, 373:515{518, 9 Feb 1995.
Lab. for Mathematical Neuroscience, RIKEN Brain Science Institute 2-1 Hirosawa, Wako, Saitama, 351-0198 Japan fhiro, [email protected]
Abstract We present an information-geometric measure to systematically investigate neuronal ring patterns, taking account not only of the second-order but also of higher-order interactions. We begin with the case of two neurons for illustration and show how to test whether or not any pairwise correlation in one period is signi cantly di erent from that in the other period. In order to test such a hypothesis of di erent ring rates, the correlation term needs to be singled out `orthogonally' to the ring rates, where the null hypothesis might not be of independent ring. This method is also shown to directly associate neural ring with behavior via their mutual information, which is decomposed into two types of information, conveyed by mean ring rate and coincident ring, respectively. Then, we show that these results, using the `orthogonal' decomposition, are naturally extended to the case of three neurons and n neurons in general.
1 Introduction Based on the theory of hierarchical structure and related invariant decomposition of interactions by information geometry [3], the present paper brie y summarizes methods useful for systematically analyzing a population of neural ring [9]. Many researches have shown that the mean ring rate of a single neuron may carry signi cant information on sensory and motion signals. Information conveyed by populational ring, however, may not be only an accumulation of mean ring rates. Other statistical structure, e.g., coincident ring [13, 14], may also carry behavioral information. One obvious step to investigate this issue is to single out a contribution by coincident ring between two neurons, i.e., the pairwise correlation [2, 6]. In general, however, it is not sucient to test a pairwise correlation of neural ring, because there can be triplewise and higher correlations. For example, three variables (neurons) are not independent in general even when they are pairwise independent. We need to establish a systematic method of analysis, including these higher-order
also aliated with Dept. of Knowledge Sci., Japan Advanced Inst. of Sci. & Tech.
correlations [1, 5, 7, 13]. We propose one approach, the information-geometric measure that uses the dual orthogonality of the natural and expectation parameters in exponential family distributions [4]. We represent a neural ring pattern by a binary random vector x. The probability distribution of ring patterns can be expanded by a log linear model, where the set fp(x)g of all the probability distributions forms a (2n ? 1)-dimensional manifold S n . Each p(x) is given by 2n probabilities X pi1 i = 1 pi1 i = Prob fX1 = i1 ; ; Xn = in g ; ik = 0; 1; subject to n
n
i1 ; ;in
and expansion in log p(x) is given by X X log p(x) = i xi + ij xi xj + i<j
X i<j 0; i; j = 0; 1. Among four probabilities, fp00; p01 ; p10 ; p11 g, only three are free. The set of all such distributions of x forms a three-dimensional manifold S 2 . Any three of pij can be used as a coordinate system of S 2 . There are many di erent coordinate systems of S 2 . The coordinates of the expectation parameters, called -coordinates, = (1 ; 2 ; 12 ), is given by i = Prob fxi = 1g = E [xi ] ; i = 1; 2; 3 = 12 = E [x1 x2 ] = p12 ; where E denotes the expectation and i and 12 correspond to the mean ring rates and the mean coincident ring, respectively. As other coordinate systems, we can also use the triplet, (1 ; 2 ; Cov [X1 ; X2 ]), where Cov [X1 ; X2 ] is the covariance,and/or the triplet (1 ; 2 ; ), where is the correlation coecient (COR), = p1 (1?12?1)12(12 ?2 ) , often called N-JPSTH [2].
Which quantity would be convenient to represent the pairwise correlational component? It is desirable to de ne the degree of the correlation independently from the marginals (1 ; 2 ). To this end, we use the `orthogonal' coordinates (1 ; 2 ; ), originating from information geometry of S 2 , so that the coordinate curve of is always orthogonal to those of 1 and 2 . The orthogonality of two directions in S 2 (S n in general) is de ned by the Riemannian metric due to the Fisher information matrix [8, 4]. Denoting any coordinates in S n by = (1 ; :::; n ), the Fisher information matrix G is given by G
= (gij ) ;
gij ( )
=E
@ @i
l (x; )
@ @j
l (x; ) :
(1)
where l (x; ) = log p (x; ). The orthogonality between i and j is de nediby h @ @ gij ( ) = 0. In case of S 2 , we desire to have E @ l (x; 1 ; 2 ; ) @ l (x; 1 ; 2 ; ) = 0 (i = 1; 2): When is orthogonal to (1 ; 2 ), we say that represents pure correlations independently of marginals. Such is given by the following theorem. Theorem 1. The coordinate p11 p00 (2) = log p p i
01 10
is orthogonal to the marginals 1 and 2 . P We have another interpretation of . Let's expand p(x) by log p(x) = p 2i=1 i xi + 12 x1 x2 ? : Simple calculation lets us get the coecients, 1 = log p10 ; 2 = 00 log pp0001 ; = ? log p00 ; and = 12 (as Eq 2). The triplet = (1 ; 2 ; 12 ) forms another coordinate system, called the natural parameters, or ?coordinates. We remark that 12 is 0 when and only when X1 and X2 are independent. The triplet (1 ; 2 ; 12 ) forms an `orthogonal' coordinate system of S 2 , called the mixed coordinates [4]. We use the Kullback-Leibler divergence (KL) to measurePthe discrepancy between two probabilities p(x) and q(x), de ned by D [p : q] = x p(x) log pq((xx)) . In the following, we denote any coordinates of p by p etc (the same for q). Using the orthogonality between - and -coordinates, we have the decomposition in the KL.
Theorem 2.
= D [p : r ] + D [r : q] ; D [q : p] = D [q : r ] + D [r : p] ; (3) q q p p p q r r where r and r are given by = (1 ; 2 ; 3 ) and = (1 ; 2 ; 3 ), respectively. The squared distance ds2 between two nearby distributions p(x; ) and p(x; ; +d) is given by the quadratic form of d, X gij ( )di dj ; ds2 = D [p : q ]
2
i;j (1;2;3)
which is approximately twice the KL, i.e., ds2 2D [p(x; ) : p(x; + d)] : Now suppose is the . Then, the Fisher information matrix 3 2 mixed coordinates is of the form g = 4 ij
g33 (d3 )2 ; ds22
=
P
g11 g12 0 g12 g22 0 0 0 g33
5
i;j (1;2) gij di dj ;
2
and we have ds2 = ds21 + ds22 ; where ds21 =
corresponding to Eq. 3.
This decomposition comes from the choice of the orthogonal coordinates and gives us the merits of simple procedure in statistical inference. First, let us estimate the parameter = (1 ; 2 ) and from N observed data x1 ; :::; xN . The maximum likelihood estimator (mle) ^, which is asymptotically unbiased and ecient, ?^2 +^12 ) , using is easily obtained by ^i = N1 #fxi = 1g and ^ = log ^(^121(1??^12^1)(^ 2 ?^12 ) 1 #fx x = 1g. The covariance of estimation error, and , is given ^12 = N 1 2 1 ?1 ?1 asymptotically by Cov = N G : Since the cross terms of G or G vanish for the orthogonal coordinates, we have Cov [; ] = 0, implying that the estimation error of marginals and that of interaction are mutually independent. Such a property does not hold for other non-orthogonal parameterization such as the COR , the covariance etc. Second, in practice, we often like to compare many spike distributions, q(x(tp)) (i.e, q(t) ) for (t = 1; ; ; T ), with a distribution in the control period p(x), or . Because the orthogonality between and allows us to treat them independently, these comparisons become very simple. These properties bring a simple procedure of testing hypothesis concerning the null hypothesis H0 : = 0 against H1 : 6= 0 ; (4) where 0 is not necessarily zero, whereas 0 = 0 corresponds to the null hypothesis of independent ring, which is often used in literature in di erent setting. Let the log likelihood of the models H0 and H1 be, respectively, l0 = max log p(x1 ; :::; xN ; ; 0 ) and l1 = max log p(x1 ; :::; xN ; ; ): ; The likelihood ratio test uses the test statistics = 2 log ll10 . By the mle with respect to and , which can be performed independently, we have ^ ; 0 ); ^ ; ^12 ); l0 = log p( x; l1 = log p( x; (5) where ^ are the same in both models. A similar situation holds in the case of testing = 0 against 6= 0 for unknown . Under the hypothesis H0 , is approximated for a large N as ^ log p(xi ; ^ ; ^0 ) N g33 ( ? 0 )2 2 (1): (6) ^ ; ) p(xi ; i=1 Thus, we can easily submit our data to a hypothetical testing of signi cant coincident ring against null hypothesis of any correlated ring, independently from the mean ring rate modulation1 . We now turn to relate the above approach with another important issue, which is to relate such a coincident ring with behavior. Let us denote by Y a variable of discrete behavioral choices. The MI between X = (X1 ; X2 ) and Y is written by p(x; y ) = Ep(Y ) [D [p(X jy) : p(X )]] : I (X; Y ) = Ep(X;Y ) log p(x)p(y )
=2
N X
j
Using the mixed coordinates for p(X y) and 0p(X ), we 0 D [ (X y ) : (X )] = D (X y ) : + D : (X ) ; (1 (X y); 2 (X y); 3 (X )) = (1 (X y); 2 (X y); 3(X )):
j
j
j
j
j
j
have D [p(X jy) : p(X )] = where 0 = 0 (X; y) =
1 A more proper formulation in this hypothetical testing can be derived, resulting in using p value from 2 (2) distribution, but we omit it here due to the limited space [9]
Theorem 3.
I (X; Y ) = I1 (X; Y ) + I2 (X; Y ); (7) where I1 (X; Y ); I2 (X; Y ) are given by 0 0 I1 (X; Y ) = Ep(Y ) D (X jy ) : (X; y ) ; I2 (X; Y ) = Ep(Y ) D (X; y ) : (X ) :
Obviously, the similar result holds with respect to p(Y jX ). By this theorem, I is the sum of the two terms: I1 is by modulation of the correlation components of X , while I2 is by modulation of the marginals of X . This observation helps us investigate the behavioral signi cance by modulating either coincident ring or mean ring rates. 0.1
η w/ s1
A
0.05
(a)
0 0
(b)
100
η1
0.05
(c)
η 2 η12
300
(b)
2
(c)
η1 η
12
100
300
500
3 C
D
θ3
2
0.1
1 0 0
100
300
0 0
500
1
0.95
E
MI (bit)
p−value
η (a)
0 0
500
0.2 COR
B
η w/ s2
0.1
0.5 0 0
100
300 time (ms)
500
0.04
100
300
500 I
F
I2
0.02
I
1
0 0
100
300 time (ms)
500
Figure 1: Demonstration of information-geometric measure in two neuron case, using simulated neural data, where two behavioral choices (s1, s2) are assumed. A,B. (1 ; 2 ; 12 ) with respect to s1, s2. C,D. COR,, computed by using = P i p(si ) (si ) with p(si ) = 1=2 (i = 1; 2). E. p-values. F. MI. Fig 1 succinctly demonstrates results in this section. Figs 1 A, B are supposed to show mean ring rates of two neurons and mean coincident ring for two di erent stimuli (s1, s2). The period (a) is assumed as the control period, i.e., where no stimuli is shown yet, whereas the stimulus is shown in the periods (b,c). Fig 1 C, D gives COR, . They look to change similarly over periods, which is reasonable because both COR and represent the same correlational component, but indeed change slightly di erently over periods (e.g., the relative magnitudes between the periods (a) and (c) are di erent for COR and ), which is also reasonable because both represent the correlational component as in di erent coordinate systems. Using in Fig 1 D, Fig 1 E shows p-values derived from 2 (1) (i.e., p > 0:95 in Fig 1 E is `a signi cance with p < 0:05') for two di erent null hypotheses, one of the averaged ring in the control period (by solid line) and the other of independent ring (by dashed line), which is of popular use in literature. In general, it becomes complicated to test the former hypothesis, using COR. This is because the COR, as the coordinate component, is not orthogonal to the mean ring rates so that estimation errors among the COR and mean ring rates are entangled and that the proper metric among them is rather dicult to compute. Once using , this testing becomes simple due to orthogonality between and mean ring rates. Notably, we would draw completely di erent conclusions on signi cant coincident ring given each null hypothesis in Fig 1 E. This di erence may be striking when we are to understand the brain function with these kinds of data. Fig 1 F shows the MI
between ring and behavior, where behavioral event is with respect to stimuli, and its decomposition. There is no behavioral information conveyed by the modulation of coincident ring in the period (b) (i.e., I1 = 0 in the period (b)). The increase in the total MI (i.e., I ) in the period (c), compared with the period (b), is due not to the MI in mean ring (I2 ) but to the MI correlation (I1 ). Thus, with a great ease, we can directly inspect a function of neural correlation component in relation to behavior.
3 Three neuron case With more than two neurons, we need to look not only into a pairwise interaction but also into higher-order interactions. Our results in the two neuron case are naturally extended to n neuron case and here, we focus on three neuron case for illustration. For three neurons X = (X1 ; X2 ; X3 ), we let p(x), x = (x1 ; x2 ; x3 ), be their joint probability distribution and put pijk = Prob fx1 = i; x2 = j; x3 = kg, i; j;P k = 0; 1. The set of all such distributions forms a 7-dimensional manifold S 3 due to pijk = 1. The -coordinates = (1 ; 2 ; 3 ) = (1 ; 2 ; 3 ; 12 ; 23 ; 13 ; 123 ) is de ned by i = E [xi ] (i = 1; 2; 3); ij = E [xi xj ] (i; j = 1; 2; 3; i 6= j ); 123 = E [x1 x2 x3 ] : To single out the purely triplewise correlation, we utilize the dual Porthogonality of - and -coordinates. By using expansion of log p(x) = i xi + P ij xi xj + 123 x1 x2 x3 ? , we obtain -coordinates, = (1 ; 2 ; 3 ) = (1 ; 2 ; 3 ; 12 ; 23p ; 13p ; p123 ).p It's easy to get the expression of these coecients 111 100 010 001 (e.g., 123 = log p110 p101 p011 p000 ). Information geometry gives the following theorem. Theorem 4. 123 represents the pure triplewise interaction in the sense that it is orthogonal to any changes in the single and pairwise marginals, i.e., 1 and 2 . We use the following two mixed coordinates to utilize the dual orthogonality, 1 = ( 1 ; 2 ; 3 ); 2 = ( 1 ; 2 ; 3 ): Here 2 is useful to single out the triplewise interaction (3 = 123 ), while 1 is to single out the pairwise and triplewise interactions together (2 ; 3 ). Note that 123 is not orthogonal to fij g. In other words, except the case of no triplewise interaction (123 = 0), ij do not directly represent the pairwise correlation of two random variables Xi ; Xj . The case of independent ring is given by ij = i j ; 123 = 1 2 3 or equivalently by 2 = 0; 3 = 0. The decomposition in the KL is now given as follows.
Theorem 5. D [p : q ]
= D [p : p] + D [p : q] = D [p : p~] + D [~p : q] = D [p : p] + D [p : p~] + D [~p : q] : (8) where, using the mixed coordinates, we have p2 = ( p1 ; p2 ; q3 ); p1~ = (p1 ; q2 ; q3 ). A hypothetical testing is formulated similarly to the two neuron case. We can exam p ine a signi cance of the triplewise interaction by 2 = 2N D [p : p] N g77 ( p2 )(123 ? q 2 2 123 ) (1): For a signi cance of triplewise and pairwise interactions together, P we have 1 = 2N D [p : p~] N 7i;j=4 gij ( p1 )(ip ? ip~)(jp ? jp~) 2 (4): For the decomposition of the MI between ring X and behavior Y , we have
Theorem 6.
I (X; Y )
= I1 (X; Y ) + I2 (X; Y ) = I3 (X; Y ) + I4 (X; Y )
(9)
where I1 (X; Y ) = Ep(Y ) [D [ 1 (X jy ) : 1 (X; y )]] ; I2 (X; Y ) = Ep(Y ) [D [ 1 (X; y ) : 1 (X )]] ; I3 (X; Y ) = Ep(Y ) [D [ 2 (X jy ) : 2 (X; y )]] ; I4 (X; Y ) = Ep(Y ) [D [ 2 (X; y ) : 2 (X )]] : By the rst equality, I is decomposed into two parts: I1 is conveyed by the pairwise and triplewise interactions of ring, and I2 by the mean ring rate modulation. By the second equality, I is decomposed di erently: I3 , conveyed by the triplewise interaction, and I4 , by the other terms. A
η
(b)
0.02
(c)
ηij η
(d)
100
300
500
700
0 −0.1 0 4
C
0
−2 0
0.1
θ123
θ12, 13, 23
i
ijk
0 0 2
0.2 B
η (a)
COR
0.04
100
300
500
1
500
700
D
0
700
−2 0
0.95
1
100
300
500
700 0.95
F χ2(4)
χ2(1)
300
2
E 0.5 0 0
100
0.5
100
300 time (ms)
500
0 700
0 0
100
300 time (ms)
500
0 700
Figure 2: Demonstration in three neuron case. A = (1 ; 2 ; 3 ) (i ; ij ; ijk ) from top to bottom, since we treated a homogeneous case in this simulation for simplicity. B. COR. C. 12 ; 13 ; 23 . D 123 . E p-value 2 (1). F p-value 2 (4). We emphasize that all the above decompositions come from the choice of the `orthogonal' coordinates. Fig 2 highlights some of the results in this section. Fig 2 A shows the mean ring rates (see legend). The period (a) is assumed as the control period. Fig 2 B indicates that COR changes only in the periods (c,d), while Fig 2 C indicates that 123 changes only in the period (d). Taken together, we observe that the triplewise correlation 123 can be modulated independently from COR. Fig 2 E indicates the p-value from 2 (1) against the null hypothesis of the activity in the control period. The triplewise coincident ring becomes signi cant only in the period (d). Fig 2 F indicates the p-value from 2 (4). The coincident ring, taking the triplewise and pairwise interaction together, becomes signi cant in both periods (c,d). We cannot observe these di erences in modulation of pairwise and triplewise interactions over periods (c, d), when we inspect only COR. Remark: For a general n neuron case, we can use the k-cut mixed coordinates, k = ( 1 ; :::; k ; k+1 ; :::; n ) = ( k? ; k+ ). Using the orthogonality between k? and k+ , the similar results hold. To meet the computational complexity involved in this general case, some practical diculties should be resolved in practice [9].
4 Discussions We presented the information-geometric measures to analyze spike ring patterns, using two and three neuron cases for illustration. The choice of `orthogonal' coordinates provides us with a simple, transparent and systematic procedure to test signi cant ring patterns and to directly relate such a pattern with behavior. We hope that this method simpli es and strengthens experimental data analysis.
Acknowledgments HN thanks M. Tatsuno, K. Siu and K. Kobayashi for their assistance. HN is supported by Grants-in-Aid 13210154 from the Ministry of Edu. Japan.
References
[1] M. Abeles, H. Bergman, E. Margalit, and E. Vaadia. Spatiotemporal ring patterns in the frontal cortex of behaving monkeys. J Neurophysiol, 70(4):1629{ 38., 1993. [2] A. M. H. J. Aertsen, G. L. Gerstein, M. K. Habib, and G. Palm. Dynamics of neuronal ring correlation: Modulation of "e ective connectivity". Journal of Neurophysiology, 61(5):900{917, May 1989. [3] S. Amari. Information geometry on hierarchical decomposition of stochastic interactions. IEEE Transaction on Information Theory, pages 1701{1711, 2001. [4] S. Amari and H. Nagaoka. Methods of Information Geometry. AMS and Oxford University Press, 2000. [5] S. Grun. Unitary joint-events in multiple-neuron spiking activity: detection, signi cance, and interpretation. Verlag Harri Deutsch, Reihe Physik, Band 60. Thun, Frankfurt/Main, 1996. [6] H. Ito and S. Tsuji. Model dependence in quanti cation of spike interdependence by joint peri-stimulus time histogram. Neural Computation, 12:195{217, 2000. [7] L. Martignon, G. Deco, K. Laskey, M. Diamond, W. A. Freiwald, and E. Vaadia. Neural coding: Higher-order temporal patterns in the neurostatistics of cell assemblies. Neural Computation, 12(11):2621{2653, 2000. [8] H. Nagaoka and S. Amari. Di erential geometry of smooth families of probability distributions. Technical report, University of Tokyo, 1982. [9] H. Nakahara and S. Amari. Information geometric measure for neural spikes. in prepration. [10] H. Nakahara, S. Amari, M. Tatsuno, S. Kang, K. Kobayashi, K. Anderson, E. Miller, and T. Poggio. Information geometric measures for spike ring. Society for Neuroscience Abstracts, 27:821.46 (page.2178), 2001. [11] M. W. Oram, N. G. Hatsopoulos, B. J. Richmond, and J. P. Donoghue. Excess synchrony in motor cortical neurons provides redundant direction information with that from coarse temporal measures. J Neurophysiol., 86(4):1700{1716, 2001. [12] S. Panzeri and S. R. Schultz. A uni ed approach to the study of temporal, correlational, and rate coding. Neural Computation, 13(6):1311{49., 2001a. [13] A. Riehle, S. Grun, M. Diesmann, and A. Aertsen. Spike synchronization and rate modulation di erentially involved in motor cortical function. Science, 278:1950{1953, 12 Dec 1997. [14] E. Vaadia, I. Haalman, M. Abeles, H. Bergman, Y. Prut, H. Slovin, and A. Aertsen. Dynamics of neuronal interactions in monkey cortex in relation to behavioural events. Nature, 373:515{518, 9 Feb 1995.
