Identifying features in spike trains using binless similarity measures
IIdentifying features in spike trains using binlesssimilarity measures
Shubhanshu Shekhar and Kaushik Majumdar25th November, 2012
Abstract
Neurons in the central nervous system communicate with each otherwith the help of series of Action Potentials, or spike trains. Various studieshave shown that neurons encode information in different features of spiketrains, such as the fine temporal structure, mean firing rate, synchronyetc. An important step in understanding the encoding of information byneurons, is to obtain a reliable measure of correlation between differentspike trains. In this paper , two new binless similarity measures for spiketrains are proposed. The performance of the new measures are comparedwith some existing measures in their ability to detect important featuresof spike trains , such as their firing rate, sensitivity to bursts and commonperiods of silence and detecting synchronous activity.
The human brain contains around neurons , which form a very large net-work interconnected with the help of nearly synapses. The neurons arecapable of generating an all-or-none impulse of voltage called the Action Poten-tial(AP) , which is the fundamental unit of information processing in the nervoussystem.Series of these action potentials are used by the neurons to communicatewith each other. It has been observed that the amplitude and the shape of theaction potentials do not vary too much over different trials . Thus it can be as-sumed that the information is contained mainly in the temporal structure of thespike trains. Single electrode and multi-electrode recordings can be performedto measure these spike trains from single neurons or a population of neurons.For the analysis of these spike trains, it is often important to quantify thesimilarity or dissimilarity between two spike trains. This can be required in thestudy of synchroniztion of the activity of a population of neurons[1] , for study-ing the reliability of neuronal response when repeatedly presented the samestimulus[2] , for testing the discrimination ability of auditory neurons[3], forbenchmarking quantitative neuron models[4] etc. Various measures have beenproposed for this purpose , such as the Victor-Purpura family of distance met-rics[5], the van Rossum distance metric [6] , the ISI distance measure [7], the1 a r X i v : . [ q - b i o . N C ] O c t orrelation based similarity measure [8], the Hunter-Milton similarity measure[9]and many more. Comparative studies of these measures have been performedin some recent works such as [10],[11] etc.From the results and observations of works comparing various simlaritymeasures , such as [10], it can be concluded that there exists no such simi-larity/distance measure which gives an optimum performance consistently in allthe benchmark tests. In particular, taking the examples of the correlation basedsimilarity measure given by Schreiber et al (SC) and the van Rossum distancemetric (VR), it was reported that the SC measure performed better than theVR metric in synchrony detection task ,but it wasn’t able to give the desiredresponse in mean firing rate detection tests , in which the VR distance metricperformed better[10] . Similar observations can be made in the case of othermeasures also . The performance of a measure is actually governed by the as-sumptions on which it has been constructed. In this paper we propose two newsimilarity measures . Our aim is to capture most of the relevant aspects of theneural activity with the help of these two measures . A similar approach wastaken in [12] , in which the authors proposed two different similarity measureswhich were sensetive to bursts and to common periods of silence respectively,and then suggested that a linear combination of the two terms be used as thesimilarity measure.The paper is organized as follows: In Sec.2 , the relation between similaritymeasures and neural code is briefly explained. The existing similarity measuresstudied in this paper are described in Sec.3 and the new similairty measures areintroduced in Sec.4 . A thorough comparison of the performance of differentsimilarity measures is presented in Sec.5 and Sec.6 contains a discussion of theresults obtained. Any idea of spike train similarity or dissimilarity is very strongly associated withthe concept of neural code. The representation of information in the temporalstructure of spike trains of a single neuron , or a population of neurons is calledthe neural code[13]. In [14] , it has been alternatively defined as the minimumset of neural symbols capable of representing all the biologically relevant infor-mation. Neural symbols could be the mean firing rate or some other statisticalfeatures of the patterns of spike trains. The term encoding time window is alsodefined in [14] as the duration of the spike train which corresponds to a singlesymbol of the neural code. Its limiting value can be obtained by taking theinverse of the maximum frequency with which the neural code is updated torepresent dynamic changes in the stimuli . If we take the encoding time windowto be N ms long, sampled at the rate of 1kHz , we can represent that segmentof a spike train as a binary vector of length N , with ’1’ present at the pointsof occurence of spikes. We can have two limiting cases of the encoding processfor the given segment. If we only count the number of spikes occuring in theencoding time window , then all the possible spike train segments can be re-2uced to a single value called the rate of firing. This coding scheme is calledthe rate coding scheme. On the other extreme if the position of each and everyspike is important , we will have a set of N different values corresponding to allthe different possible spike train segments. Neural coding requires that the N different possible binary vectors be reduced to a smaller set of neural symbols,or aspects of neural activity which carry the relevant information. In temporalcoding it is assumed that not just the number of spikes occuring within theencoding window , but also the pattern of spikes in the window carry significantinformation.In this paper we focus mainly on four different aspects of neural activitywhich might carry important information , namely rate coding, temporal cod-ing and interval coding , and information coding by common periods of silence.In rate coding , it is assumed that within a particular encoding time window, it is only the number of spikes occuring which corresponds to the informa-tion about the stimulus. Rate coding is particularly likely in case of neuronswhere the integration time is larger than the mean Inter Spike Intervals(ISI).[13]As explained in [14] , for stimuli with a single time scale, rate coding is onlyconcerned with the number of spikes in the encoding window. If however, thestimulus is dynamic in time, rate coding requires correlation between same fre-qency components of the stimulus and the corresponding spike train , in thefrequency domain . Various methods for estimating the firing rates of spiketrains have been discussed in [[15]Temporal coding hypothesis assumes that the relevant information is con-tained in the precise patterns of the spikes within the encoding window, notjust their count. Analogous to the definiton of rate coding in frequency do-main for dynamic stimuli , requirement for temporal coding as defined in [14],is that there exists some correlation between certain frequency components ofthe stimuli with higher frequency components of the spike trains. Thus , withinthe framework of temporal coding, two spike trains which have different timesof occurence of spikes , even though both have equal number of spikes, will beconsidered different. Alternatively, two neurons which fire synchronously, willbe considered to be similar.Another aspect of information encoding by neurons has been reported in [16], called the interval coding scheme. In that paper , it is reported that the shortInter Spike Intervals(ISIs) that occur during bursts are a distinct feature of theneural code. They have also given an account of some experimental resultson pyramidal cells, which when presented with a broadband current injection,responded with bursts , in which the burst ISIs were correlated with the intensityof stimulus upstrokes. It has furthur been stated that the number of spikes in aburst is also correlated with the slope of the stimulus upstrokes. This encodingof stimuls features in the burst ISIs was termed as interval coding. The role ofbursts in information coding has also been studied in some other previous works.In [17] , it is reported that bursts encode information in their timing and theirdurations and that the bursts are much more reliable than isolated spikes, i.e thetiming of the bursts are more reproducible across various trials than the timingof isolated spikes. Also the reliability was in direct proportion to the length of3he bursts. In [18] , the authors report that although the signalling in manycentral synapses are very unreliable to individual spikes , the reliability increasesconsiderably to the bursts because of facilitation. Thus those synapses act as afilter , which rejects single spikes while allowing bursts to pass through. Burstshave also been shown to encode information about the direction of movementmore reliably than isolated spikes in electrosensory midbrain neurons [19].Another important feature of neural activity, observed mainly in the Cere-bellar Purkinje cells (PC) , are the long periods of silence, which are oftensynchronized over different PCs.[20] The presence of long pauses in the activityof these neurons hints at an inherent bistability in their activity in which themembrane potential transitions between states of continuous firing and periodsof quiescence [21]. It has been reported that approximately half of the cerebellarPCs exhibit this behaviour of alternating between long pauses and firing simplespikes [22]. Some experimental studies have shown that the synchronous pausesin firing of converging PCs can induce Deep Cerebellar Nuclei to fire in a reliablefashion [23]. Thus synchronous periods of silence also seem to be a significantaspect of neural code.In this paper we propose two new measures of spike train similarity . Withthese two measures we intend to capture all the four different firing propertiesof neurons discussed above. We also compare the ability of some existing spiketrain similarity/distance measures in detecting these firing properties. A large number of measures for spike train similarity and dissimilarity have beenproposed , for the purpose of spike train analysis. In this paper , the followingspike train similarity / distance measures were studied:-
Schreiber et al [8] proposed a new corrleation based measure for quantifyingthe spike timing reliability . This measure is dependent on a single parameter,which is related to the timescale of the precision of spike timing. For calculationof this measure , first the convolution of the two spike trains, represented as asum of delayed dirac delta functions , is performed with a gaussian kernel,having a standard deviation σ c . Then the inner product of the two waveformsobtained after smoothing the spike trains is calculated , and normalized withthe norms of the individual waveforms, to obtain the measure of similarity.Basically in the discrete time case, if we consider the smoothed spike trainsto be vectors in an m dimensional vector space, m being the number of timepoints, the similarity measure is the cosine of the angle between the two vectors.Thus if the spike trains are represented by s1 and s2 and the gaussian kernelwith standard deviation σ c is represented by G( σ c ), the the similarity measurebetween s1 and s2 can be calculated as4 = s ∗ G ( σ c ) and x = s ∗ G ( σ c ) (1) S corr ( s , s
2) = x .x | x || x | (2) This distance metric was introduced by van Rossum in [6] . Like the SC mea-sure, this distance also depends on one parameter which defines the timescaleconsidered. For the calculation of this distance measure, first the spike trainsare convolved with a decaying exponential signal . As a result of this convolu-tion we get a smooth curve representing the spike train. From the two smoothsignals, the distance is calculated as their L norm. The calculations involvedare as follows f ( t ) = n s (cid:88) i =1 u ( t − t i ) e ( t − ti ) /τ s (3) D V R ( s , s
2) = 1 τ s (cid:90) ∞ t =0 [ f − f dt (4)The parameter τ s determines the time scale over which the influence of a singlespike is extended. When the value of this paramter is taken to be very small,the range of influence of each spike is very limited and thus the measure acts asa coincidence detector. In the limiting case of τ s − > , the distance measuregives the count of non coincident spikes. On the other extreme , with τ tendingtowards ∞ , the distance measure returns the difference in the total spike countof the two spike trains. Another spike train distance measure which was proposed in [7] , is the ISIdistance measure. This distance measure is different from the other measures,in that it takes into account the inter spike intervals instead of the actual spikeoccurence times for the calculation of the distance value. For calculating thisdistance measure , first a spike train is represented as the set of time of occurenceof spikes, i.e S = { t , t ...t n }. To each such spike train, two additional spiketimes corresponding to the start of the time window , and the end of the timewindow are added. That is, if the entire duration of the spike train is T ms,then S is modified to include 0 and T as the first and last terms. From thismodified spike train , a function f is calculated , such that f(t) = t i +1 − t i ,for t i < t < = t i +1 . Thus for two spike trains S1 and S2, after obtaining thecorresponding f functions , f1 and f2, another function I is calculated in thefollowing way 5 ( t ) = 1 − min ( f t ) , f t )) max ( f t ) , f t )) (5)From the function I , the distance between the two spike trains is calculatedas D ISI ( s , s
2) = 1 T (cid:90) Tt =0 I ( t ) dt (6)An advantage of the ISI distance measure is that it does not invlove anyfree parameter which has to be selected by the user. However , as a result ofthis very property , this measure cannot be used to study one particular featureof the encoding of information in the spike train ,as it cannot be adjusted toconcentrate on a particular aspect of the spike train . Another measure of similarity was proposed in [9] , by Hunter and Milton. Inthis scheme also the spike trains are represented as a set of their occurencetimes. For each element in the first spike train ( t i ), its nearest neighbour inthe second spike train ( t i ) is obtained. After getting the nearest neighbour ,degree of synchronization or coincidence is obtained by calculating e − | t i − t i | τH .These values are then calculated for all the spikes in the first spike train, andtheir mean is taken t= to obtain r . Similarly r is calculated , by taking themean of the coincidence values over all the spikes occuring in the second spiketrain. Finally the similarity measure is the arithmetic mean of r and r . Likethe SC and the VR measures, the HM measure also has a free parameter τ H ,which determines the range of coincidence of the spikes. In [12] , the authors proposed two similarity measures , one sensetive to burstsand the other sensetive to common periods of silence. A convex combinationof these two measures could then accordingly be used as a similarity measure ,depending upon the situation.
This measure is a modification of the SC measure. In this measure, the firststep involves convolving the spike train with a gaussian kernel of width σ toobtain f(t). A piecewise linear transformation N(x) was applied to each f(t) ,where N(x) was defined as N ( f ( t )) = H ( f ( t ) − ηT )( f ( t ) − ηT ) (7)Here T is a threshold function which decides which segment of the spike trainis considered to be a burst. It depends upon n( minimum number of spikes in a6urst) , b( maximum inter burst ISI) and σ width of the gaussian kernel. Theparameter η takes values in the range (0,1) and can be used to set the amountof emphasis to be given to bursts over isolated spikes. T can be calculated bythe following expression. T ( n, b, σ ) = k = n (cid:88) k =1 e ( p − kb )2 σ (8) p b ( n + 1) / p b ( n + 2) / T = max ( T ( p , T ( p (9)The similarity measure is then calculated by taking the standard correlationmeasure of the N(f(t) ) of the two spike trains. For calculating this measure, the spike trains are mapped to a function(g(t))which is set to zero at all the spikes, and then rises linearly in the inter spikeinterval . The linearly rising part begins after a time delay τ , in order to ignorepauses of smaller lengths. τ can be set to the mean value of all the ISIs , sothat only ISI of large durations are considered in this similarity measure. Afterobtaining g(t) for both the spike trains, the standard correlation measure iscalculated for the similarity measure. In this section , we define the following two similarity measures ( SM1 and SM2):- The first similarity measure (SM1) is a modified version of the SC measure.In this measure , we first construct a function f from the given spike train .In the case of the SC measure , the function f is obtained by performing theconvolution of the spike train with a gaussian kernel of a certain width. Hereinsted of using a gaussian kernel , we obtain the smoothed signal using thefollowing differential equations. df ( t ) dt = − f ( t ) τ f + S ( t ) .u ( t ) (10) du ( t ) dt = − u ( t ) − u τ u + ∆ u.S ( t ) (11)The function S(t) represents the spike train as a sum of delayed dirac deltafunctions . The initial value of the variable f is taken to be zero, i.e f(0) = 0.If a spike occurs at time t i , the value of variable f jumps by an amount u( t i ).This is because of the following property of the delta function. (cid:90) t i + (cid:15)t i − (cid:15) δ ( t − t i ) dt = 1 (12)7etween two spikes , the value of the variable f decays exponentially with atime constant τ f . This parameter is similar to the corresponding parametersin the VR and SC measures, in the sense that it defines the timescale overwhich the coincidence of two spikes is considered. The amount by which thevariable f jumps with each spike is not constant, it is represented by anothervariable u. The variable u also jumps at the occurence of every spike by afixed value u , and decays exponentially with a time constant τ u between twospikes. The value of the time constant τ u is taken to be much smaller than τ f .This method of generating the variables f and u is quite similar to some simplephenomenological models of facilitating synapses[24]. A similar approach wastaken in [25] , in which the van Rossum distance measure was studied with asynaptic filter . However , in that paper ,the authors reported that they didnot get any significant gain in the performance with facilitating synapses. Inthis paper , our motivation for using a model analogous to that of a facilitatingsynapse, is to exploit its properties to make the measure more sensetive tobursts. Whenever a spike occurs , the variable f jumps by a value u( t − ), where t − is the time just before the arrival of the spike. With the arrival of the spike, the variable u(t) also jumps by a value u , i.e u ( t + ) = u ( t − ) + ∆ u . Both thevariables then start decaying with their respective time constants. If ,however,another spike arrives before u has decayed to a value close to its base valueof u , the next jump in f would be higher than the previous jump. This kindof situation occurs in the case of bursts , where a group of spikes occur withvery small inter spike intervals . Thus if the value of the time constant τ u is comparable to the interspike intervals , the variable f will rise considerablyduring the occurence of the bursts , and if the bursts occur simultaneously intwo spike trains, it will be reflected in the higher value of the similarity measure.An important distance measure used in information theory is the Hammingdistance , which calculates the distance between two vectors as the number ofpositions at which they differ. Alternatively it counts the minimum number ofchanges to be made to transform one vector to the other. It is not feasible todirectly apply the Hamming distance to the binary representation of two spiketrains, as it would consider two coinciding spikes as no different from two pointsfrom inter spike intervals. For the second distance measure(SM2), we use thebasic idea of the Hamming distance and modify it to make it more suitable forapplication to spike trains.Given two spike trains S1 and S2, we first obtain their smoothed versions byconvolving them with a suitable smoothing kernel. We have used the decayingexponential kernel in this paper, but other kernels can also be applied. Let thesmoothed spike trains be r1(t) and r2(t). Then the similarity measure (SM2) iscalculated as follows: d S , S
2) = 1 T (cid:90) Tt =0 x ( t ) dt (13) SM S , S
2) = 1 − d (14)where x is defined as 8 ( t ) = (cid:26) , min ( r t ) , r t )) ≤ k ∗ max ( r t ) , r t ))0 , otherwise Here two points of r1(t) and r2(t) are taken to be similar if the smaller ofthe two values is greater than k times the larger. The parameter k ∈ (0,1) , andit denotes the amount of tolerance allowed for two points of r1(t) and r2(t) tobe considered similar .A typical value we used in our simulations was 0.7. Usingthis definition , the distance d2 is calculated for r1 and r2 , and 1 - d2 is definedas the similarity measure. In this section , we compare the performance of all the different similarity mea-sures discussed so far in the paper. We will be comparing the performance ofthe measures in the following • Firing rate discrimination • Burst Sensitivity • Sensitivity to common periods of silence • Syncronous firing detectionThese tests correspond to the different aspects of neural coding that we discussedin Section 2.
To study the firing rate discrimination ability of the measures we first used atest similar to the one described in [10]. For the test, we generated artificialspike trains , each 5s long, using a homogenous poisson process[26]. The meanfiring rate of one spike train was kept fixed at 20 Hz, while the firing rate ofthe second spike train was varied from 2Hz to 40 Hz in steps of 2. For eachcombination of firing rates, 100 pairs of spike trains were generated , and theentire process was repeated for 1000 times. The parameters used are shown inTable 1.Figure 1 shows the performance of the similarity measures. If a similaritymeasure has the ability to detect differences in the mean firing rate of two spiketrains, then the similarity values should be maximum for the rate of 20Hz , andit should decrease as we move away from 20 Hz on both sides. For plotting allthe measures on the same graph , we had to modify some of the measures. TheISI was converted into a similarity measure by taking S ISI = 1 − D ISI . Themaximum mean value of the van Rossum distance measure over all the trialswas around 364. So , for showing it on the same plot, we converted it into asimilarity measure by defining S V R = 1 − D V R / . The performance of SM1,9 imilarity Measure Parameter Values used SM1 τ f =100ms , τ u =5ms , u =0.3 , ∆ u =0.2SM2 τ =100 ms , k = 0.7SC σ c = 100 msVR τ s = 100 msHM τ H =100 msISI -Table 1: Various parameter values used in the firing rate test mean firing rate s i m il a r i t y m eas u r e va l u es Firing rate discrimination
SM2SM1SCISIHMVR
Figure 1: The figure shows the performance of the six similarity measures , inthe test for firing rate difference discrimination. The VR and SM2 were the twobest performers in this test.SC and HM measures are almost identical. They show the correct behaviourwhen the firing rate of the second spike train is below the reference firing rate(20 Hz) . However for higher firing rates, the similarity measure values eitherstay almost constant or go up instead of decreasing. The two best performingmeasures were the VR and SM2 measures.To furthur study the ability of the similarity measures to classify spike trainsbased on their mean firing rates , we generated 8 sets of spike trains of 5s lengthwith mean firing rates varying from 5Hz to 40 Hz , with each set containing 50spike trains. We used the classification scheme suggested in [5] . A confusionmatrix N( r i , r j ) was constructed using the various similarity/distance measures.N( r i , r j ) represents the number of times a spike train of mean firing rate r i hasbeen classified in the set corresponding to mean firing rate r j . The confusionmatrix is initialized as a null matrix of dimension N s XN s ( N s =number ofsets=8) . For each spike train s, its mean distance from the spike trains of the10istance/Similarity Measure H/H max
SC 0.045VR 0.426ISI 0.334SM1 0.035SM2 0.375HM 0.0267LFB 0.021LFS 0.142Table 2: The performance of the different measures in clustering the spike trainson the basis of their mean firing rates are shown in the table. The H values havebeen normalized by H max , which corresponds to the case of perfect clustering.different sets is calculated according to the following formula. d ( s, r i ) = (cid:20)(cid:68)(cid:0) d ( s, s r i ) (cid:1) z (cid:69)(cid:21) z (15)Here d ( s, s r i ) is the distance between the spike train under consideration s,and a spike train from the set with mean firing rate r i . The given spike trains is alloted to the set , with which it gives the minimum value of the meandistance. Using this process the confusion matrix is completed. For an idealclassifier the confusion matrix N would be a diagonal matrix. To quantify theperformance of the classification , we used the information theoretic measure ,transmitted information (H) [5] . We normalized the values of H by its maximumvalue H max , which is given by H max = log ( N s ) corresponding to the case ofperfect classification , where N s is the number of sets in which the data is to beclassified. The values of normalized H are listed in the Table 2From the Table 2, we can see that the VR measure performed the best inclassifying spike trains based on their mean firing rates , followed by the SM2and ISI measures. All the similarity measures which involved calculating innerproducts performed poorly in this test. Bursts are a very important component of neural signalling , as they are knownto increase the reliability of synapses in the central nervous system. Manysynapses transmit bursts but do not respond to single spikes , thus acting asa filter [18]. To test the sensetivity of similarity measures to bursts , we willfollow a procedure similar to the one used in [12]. We first generated twopoisson spike trains , 5 seconds long , with mean firing rates of 20 spks/sec and30 spks/sec respectively. The similarity/distance measures were calculated forthese two spike trains. Then , the two spike trains were modified by adding somebursts at same positions and deleting an equivalent number of isolated spikes.11 bursts
ISI VR SC HM LFS LFB SM2 SM13 1.019 0.998 1.309 1.146 1.014 1.9544( η =0.5) 1.132 1.5236 1.038 0.996 1.541 1.274 1.019 2.940( η =0.5) 1.266 1.881Table 3: The S B values obtained for the different measures in the two cases areshown in the tableThe similarity /distance measures were again calculated for the modified spiketrains. This process was repeated for 1000 times. For a similarity ( distance )measure sensitive to bursts, the value should be higher (lower ) in the case ofthe modified spike trains. To quantify the sensitivity of the measures to burstsof action potential , the following term was calculated S B = 1 n i = n (cid:88) i =1 Sim ( s ∗ , s ∗ ) Sim ( s , s (16)Here Sim denotes any similarity measure and n is the total number of tri-als. s1 and s2 are the original spike trains, whereas s ∗ and s ∗ represent themodified spike trains. In the case of VR measure , the ratio was inverted , whilethe ISI measure was changed into a similarity measure, by taking its differencefrom unity. Thus for a measure which emphasizes bursts of action potentialsmore than isolated spikes, the value of S B should be greater than one. Table 3shows the values of S B obtained for various measures.The LFB measure and the SM1 measures were the most versatile measures, with regards to burst sensitivity. The parameter η in LFB , and u and ∆ u in SM1 can be tuned to adjust the emphasis to be given to the bursts. The S B value for SM1 in the table correspond to the values of ( u , ∆ u ) = (0.2, 0.4) . Byincreasing ∆ u and reducing u ( or keeping it constant ) , the weightage givento bursts can be increased.So when ( u , ∆ u ) were changed to (0.1,0.5) , the S B values of SM1 increased to 1.887 and 2.43 with the addition of 3 and 6 burstsrespectively. Other measures which performed well in this test were SC , HMand SM2 . In the case of LFS, ISI and VR , the introduction of bursts did notseem to have much effect on the similarity/distance measure values. Common periods of quiescence have also been suggested as an important aspectof information coding , especially in the case of cerebellar micorcircuits. Totest the sensitivity of measures to common periods of silence , we followedthe following procedure . Two poisson spike trains of length 5ms each, weregenerated and the similarity/distance values between them were obtained. Thethe two spike trains were modified to include a period of silence of length L s atthe same position , and the similarity/distance values were calculated for thesetwo trains. The value of L_s was varied from 100 ms to 500 ms, and the entireprocess was repeated 1000 times. The sensitivity was quantified using a term S S , defined in a similar way to S B . 12
00 150 200 250 300 350 400 450 50011.522.533.54
Variation of LFS with length of silent period
Length of silent period
100 150 200 250 300 350 400 450 5000.9511.051.11.151.21.251.3
Length of silent period
Variation of SM2 with length of silent period
Figure 2: LFS and SM2 the only measures which displayed any significantvariation with the introduction of a silent segment. LFS has been speciallyconstructed to be responsive to silent periods , and it showed very large variationwith the length of segment inserted. SM2 varied almost linearly with the lengthof silent period, although the slope of variation was quite small as compared toLFSIt was observed that the insertion of common silent periods had negligibleeffect in the case of all the measures except LFS and SM2 . LFS was themost responsive to the introduction of common periods of silence , with its S S values increasing rapidly with increasing length of the segment introduced. SM2showed an almost linear variation with increasing length of the silent segmentembedded, although the rate of increase was much smaller than in the case ofLFS. Synchronous neural activity is an important component of information process-ing in the brain. It has been hypothesized to be involved in information encod-ing, influencing transmission of activity from one group of neurons to another,facilitating a group of neurons with common post synaptic targets to depolarizethem more effectively.[13]. Here , by synchronous spike trains , we mean spiketrains in which the timings of spikes are correlated.For studying the ability of measures to detect synchronous firing , we gener-ated the dataset using a procedure resembling those used in [8] and [10]. First areference poisson spike train with mean firing rate of 20 spks/sec was generated. From this reference spike train , new spike trains were generated using thefollowing procedure:- • Each spike of the reference spike train was considered individually. Anyspike of the reference spike train was retained in the new spike train with13
Jitter Standard Deviation
Variation of SM1, SC , HM and SM2
SM1SCHMSM2
Jitter Standard Deviation
Variation of VR
VR0 5 10 15 20 25 3000.20.40.60.81
Jitter Standard Deviation
Variation of ISI
ISI
Jitter Standard Deviation
Variation of LFS and LFB
LFSLFB
Figure 3: The figure shows how the different measures responded to changein the values of jitter in spike timing. All the measures except LFS and LFBfollowed the expected trend of decreasing similarity with increasing jitter inspike timing.a probability p r . The term p r reflects the reliability of the data . A highvalue of p r (i.e close to 1) corresponds to greater reliability. • The time of occurence of the spikes of the reference spike train was carriedover into the new spike train with some jitter . The amount of jitter wastaken randomly from a gaussian distribution with a standard deviation of σ j . The term σ j represents the precision in the timing of the spikes overdifferent trials. Smaller σ j means greater precision in spike timing. • To represent background activity , the spike train obtained after the pre-vious two steps was superimposed with a background spike train, withmean firing rate r b . The value of r b was taken to be much smaller thanthe mean firing rate of the reference spike train.The values of σ j was varied in the range (1, 30) ms in steps of 1 , while p r ranged from 0.40 to 0.98 in steps of 0.02. r b was kept fixed at 5 spks/sec. Foreach combination of ( σ j , p r ) , ns = 50 spike trains were generated . In this waya large number ( 30*30*50) of spike trains were generated , and the behaviourof the different measures were studied on this data set.In Figure 3, the variation of the different measure values for a fixed p r (0.8)and varying σ j is shown. As the σ j value increases, the precision of spike timingover different trials decreases. So it is expected that the similarity(distance )values should decrease ( increase) with increasing σ j . This general trend isfollowed by all the different measures except LFS and LFB .Figure 4 shows how the different measures vary with σ j kept fixed (10 ms), and varying p r . As the value of p r is increased , the reliability of the spikes14 .4 0.5 0.6 0.7 0.8 0.9 100.20.40.60.81 Probably of retaining a spike p r Variation of SM1, SC , HM and SM2
SM1SCHMSM2
Probably of retaining a spike p r Variation of VR VR Probably of retaining a spike p r Variation of ISI
ISI
Probably of retaining a spike p r Variation of ISI
LFSLFB
Figure 4: The figure shows the variation of the different measures , with in-creasing reliability of spikes . All the measures except LFS and LFB followedthe expected trend.Set No. σ j (ms) p r r b (spks/sec)Set 1 (12,20) (0.85,1) (6,9)Set 2 (6,12) (0.70,0.85) (3,6)Set 3 (0,6) (0.55,0.70) (0,3)Table 4: The table shows the range of values of the three parameters used forcreating three different sets of spike trains.over different trials increases . Thus the values of similarity (distance) measuresshould increase ( decrease) with increasing probability of retention of spikes . SC, SM1 ,SM3 , HM and ISI measures showed the desired variation with increasingvalues of p r . Although the VR measure generally followed the desired trend,it also showed large fluctuations. In this case also the LFS and LFB failed torespond in the expected manner.To furthur study the ability of the different measures to detect synchrony ,they were used in the following classification problem . First a reference poissonspike train with a mean firing rate of r was generated. From this referencespike train , three different sets of spike trains with different ranges of the threeparameters σ j , p r and r b were obtained. The values of these parameters usedare given in Table 3 . 100 spike trains for each set were genrated by randomlyselecting the different parameters from the ranges shown in the table. The firstset corresponds to spike trains with high reliability , but low precision in spiketimes and high background activity . The third set contains spike trains withlow reliability , but high precision and small amout of background activity. Inthe second set , moderate values of reliability , precision and background noiseare used. The same clustering scheme , which was used in classifying spike trains15able 5: The performance of the different measures in the clustering problem , atlow mean firing rate ( r =15 spks/sec) . The SM1 showed the best performancein this case.Similarity/Distance Measure H/ H max SM1 0.455VR 0.118SC 0.392SM2 0.147HM 0.386ISI 0.352LFB 0.246LFS 0.064Similarity/Distance Measure H/ H max SM1 0.516VR 0.327SC 0.557SM2 0.598HM 0.4536ISI 0.1325LFB 0.367LFS 0.148Table 6: The table lists the performance index for the different measures inclustering the spike trains, when the mean firing rate of the reference spiketrain was high( 50 spks/sec). In this case the best performing measure was SM2, closely followed by SC and SM1.based on their mean firing rates , is used here too. Two different values of themean firing rate of the reference spike train ( 15 spks/sec and 50 spks/sec) wereused corresponding to low and high values of mean firing rate. The results ofthe two clustering problems are shown in Table4 and Table 5. Some measureslike SM1 , SC and HM performed well in the classification task at both, highas well as low mean firing rates , while SM2 , LFB and VR worked better inthe case of spike trains with higher firing rate. ISI measure on the other hand ,performed better with r = 15. In this paper, we have introduced two new measures for quantifying similaritybetween a pair of spike trains. The first of these measures (SM1) is a modifi-cation of an existing correlation based measure (SC) proposed in [8] . Here ,instead of using a gaussian kernel for smoothing the spike trains, we have em-16loyed a process which is similar to that used in obtaining the post synaptictrace in phenomenological models of facilitating synapses. Thus, in general thismeasure exhibits behaviour quite similar to that of the SC measure, but has anadded feature of being burst sensitive. The burst sensitivity of this new measurecan be tuned according to the requirements, by varying the parameters u and ∆ u . Like the SC measure, this new measure also did not perform well in classi-fying spike trains based on their mean firing rates. However , using this measure, very good classification of spike trains based on synchrony was obtained.The second similarity measure(SM2) introduced in this paper , has beenmotivated from the well known distance measure used in information theory ,called the Hamming distance. Hamming distance between two vectors is thenumber of positions at which the two vectors differ . In case of spike trains,two coincident action potentials are a much more significant event than twopoints on the silent parts of the spike trains. However , these two occurenceswill be treated equally in calculating Hamming distance , and hence it cannotbe directly applied to spike trains. So , for calculating SM2, we first smoothedout the spike trains with an exponential kernel. Then two corresponding pointsof the smoothed spike trains were defined to be equal , if the absolute valueof their difference is less than (1-k) times the greater of the two values. Usingthis definition , the modified Hamming distance (d2) was calculated , and SM2was obtained by subtracting d2 from 1. This measure performed well in theclassification of spike trains based on their mean firing rate , as well as in theclassification of spike trains based on synchronous firing, when the mean firignrate of the reference spike train was high(50 spks/sec). SM2 was also the onlymeasure , other than LFS , which was sensitive to common pauses introducedin the spike trains . However , unlike the LFS measure, it also performed wellin other tests.We compared the performance of the two new measures , with some existingmeasures on their ability to detect important features of spike trains , suchas their firing rate, sensitivity to bursts and common periods of silence anddetecting synchronous activity. SM1 showed the ability to emphasize burstsover isolated spikes , and also performed well in detecting synchronous firing .SM2 , on the other hand, performed well in classification based on firing rate, and also exhibited sensitivity to common periods of silence. SM2 was alsoable to outperform other measures in classifying the spike trains at higher meanfiring rate. In this way , these two measures can be used to identify all the fourimportant characteristics of spike trains mentioned above.In future , work can be done to develop new spike train clustering techniqueswhich utilizes both the new similarity measures proposed in this paper. Sucha clustering technique could benefit from the complementary advantages of thetwo measures , and provide optimum classification of spike trains.17 eferenceseferences