Table Of ContentMalte J. Rasch, Arthur Gretton, Yusuke Murayama, Wolfgang Maass and Nikos
K. Logothetis
J Neurophysiol 99:1461-1476, 2008. First published Dec 26, 2007; doi:10.1152/jn.00919.2007
You might find this additional information useful...
This article cites 47 articles, 16 of which you can access free at:
http://jn.physiology.org/cgi/content/full/99/3/1461#BIBL
Updated information and services including high-resolution figures, can be found at:
http://jn.physiology.org/cgi/content/full/99/3/1461
Additional material and information about Journal of Neurophysiology can be found at:
http://www.the-aps.org/publications/jn
This information is current as of April 14, 2008 .
D
o
w
n
lo
a
d
e
d
fro
m
jn
.p
h
y
s
io
lo
g
y
.o
rg
o
n
A
p
ril 1
4
, 2
0
0
8
Journal of Neurophysiology publishes original articles on the function of the nervous system. It is published 12 times a year
(monthly) by the American Physiological Society, 9650 Rockville Pike, Bethesda MD 20814-3991. Copyright © 2005 by the
American Physiological Society. ISSN: 0022-3077, ESSN: 1522-1598. Visit our website at http://www.the-aps.org/.
JNeurophysiol99:1461–1476,2008.
FirstpublishedDecember26,2007;doi:10.1152/jn.00919.2007.
Inferring Spike Trains From Local Field Potentials
Malte J. Rasch,1 Arthur Gretton,2 Yusuke Murayama,2 Wolfgang Maass,1 and Nikos K. Logothetis2,3
1InstituteforTheoreticalComputerScience,GrazUniversityofTechnology,Graz,Austria;2MaxPlanckInstituteforBiological
Cybernetics,Tu¨bingen,Germany;and3ImagingScienceandBiomedicalEngineering,UniversityofManchester,Manchester,
UnitedKingdom
Submitted16August2007;acceptedinfinalform24December2007
RaschMJ,GrettonA,MurayamaY,MaassW,LogothetisNK. is used in most recordings to obtain multiple-unit spiking
Inferring spike trains from local field potentials. J Neurophysiol activity(MUA),andalow-passfiltercutoffofabout300Hzto
99: 1461–1476, 2008. First published December 26, 2007; obtain the so-called local field potentials (LFPs). Numerous
doi:10.1152/jn.00919.2007.Weinvestigatedwhetheritispossibleto
experiments have presented data indicating that such a band
infer spike trains solely on the basis of the underlying local field
separation does indeed underlie different neural events (for
potentials (LFPs). Using support vector machines and linear regres-
references see, e.g., Logothetis 2003).
sion models, we found that in the primary visual cortex (V1) of
In summary, depending on the recording site and the elec-
monkeys, spikes can indeed be inferred from LFPs, at least with
moderatesuccess.Althoughthereisaconsiderabledegreeofvariation trode properties, the MUA most likely represents a weighted
across electrodes, the low-frequency structure in spike trains (in the sumoftheextracellularactionpotentialsofallneuronswithin D
100-ms range) can be inferred with reasonable accuracy, whereas a sphere whose radius is about 140 to 300 (cid:2)m, with the ow
exactspikepositionsarenotreliablypredicted.Twokindsoffeatures electrodeatitscenter(Henzeetal.2000).Spikesproducedby n
lo
oftheLFPareexploitedforprediction:thefrequencypowerofbands the synchronous firings of many cells can, in principle, be a
in the high (cid:1)-range (40–90 Hz) and information contained in low- enhanced by summation and thus detected over a larger dis- de
furleaqtiuoennscyaroesciinllfaotrimonasti(v(cid:1)e.10InHfozr)m,wathioenreabnoatlhyspihsasreevaenadlepdowtheartmboodth- tance (Arezzo et al. 1979; Huang and Buchwald 1977). In d fro
general, experiments have shown that large-amplitude signal m
featurescode(mainly)independentaspectsofthespike-to-LFPrela- variations in the MUA range reflect large-amplitude extracel- jn
tcilounsstehriepd, wspitihkitnhge alocwtiv-fitrye.quAenltchyouLgFhPbpohthasefecaotudriensgafnodr tepmrepdoicrtailolyn lular potentials and that small-amplitude fast activity is corre- .ph
y
qualityaresimilarduringseminaturalmoviestimuliandspontaneous latedwithsmallones(BuchwaldandGrover1970;Gasserand sio
activity,predictionperformanceduringspontaneousactivitydegrades Grundfest 1939; Grover and Buchwald 1970; Hunt 1951; lo
much more slowly with increasing electrode distance. The general Nelson 1966). gy
trend of data obtained with anesthetized animals is qualitatively The low-frequency range (i.e., the LFPs) of the mEFP .o
rg
mirrored in that of a more limited data set recorded in V1 of signal, on the other hand, represents mostly slow events re- o
non-anesthetizedmonkeys.Incontrasttothecorticalfieldpotentials, flecting cooperative activity in neural populations. Initially n A
thalamicLFPs(e.g.,LFPsderivedfromrecordingsinthedorsallateral these signals were thought to represent exclusively synaptic p
gacetniivciutyla.tenucleus)holdnousefulinformationforpredictingspiking eavnedntBso(nAdjm19o6n4e-,M1a9r6s7a)n. 1E9v6id5e;nBcuechfowratlhdeeirt aolr.ig1i9n65w;aFsroomftemn ril 14, 2
gatheredfromcurrent-sourcedensity(CSD)analysisandcom- 0
0
bined field potential and intracellular recordings (Mitzdorf 8
INTRODUCTION
1985; Nadasdy et al. 1998). Mitzdorf suggested that LFPs
In a typical electrophysiology experiment, the signal mea- actually reflect a weighted average of synchronized dendroso-
sured by an electrode placed at a neural site represents the matic components of the synaptic signals of a neural popula-
mean extracellular field potential (mEFP) from the weighted tionwithin0.5–3mmoftheelectrodetip(Juergensetal.1999;
sum of all current sinks and sources along multiple cells. If a Mitzdorf 1987). Later studies, however, provided evidence of
microelectrode with a small tip is placed close to the soma or the existence of other types of slow activity unrelated to
axonofaneuron,thenthemeasuredmEFPdirectlyreportsthe synapticevents,includingvoltage-dependentmembraneoscil-
spiketrafficofthatneuronandfrequentlythatofitsimmediate lations (e.g., Kamondi et al. 1998) and spike afterpotentials.
neighbors as well. If the impedance of the microelectrode is Taken together, LFPs represent slow waveforms, including
sufficientlylowanditsexposedtipisabitfartherfromasingle synaptic potentials, afterpotentials of somatodendritic spikes,
large pyramidal cell, so that action potentials do not predom- andvoltage-gatedmembraneoscillations,thatreflecttheinput
inate the neural signal, then the electrode can monitor the of a given cortical area as well as its local intracortical
totality of the potentials in that region. The EFPs recorded processing, including the activity of excitatory and inhibitory
undertheseconditionsarerelatedbothtointegrativeprocesses interneurons.
(dendritic events) and to spikes generated by several hundred Given the different natures of LFPs and MUA, we felt that
neurons. it would be interesting to address the question of whether one
Thetwodifferentsignaltypescanbesegregatedbyfrequency- can infer spiking of neurons from the locally measured field
bandseparation.Ahigh-passfiltercutoffofabout300–500Hz potentials.Hereinweaddressthisquestioninastraightforward
Address for reprint requests and other correspondence: N. K. Logothetis, Thecostsofpublicationofthisarticleweredefrayedinpartbythepayment
Max Planck Institute for Biological Cybernetics, Spemannstrasse 38, 72076 ofpagecharges.Thearticlemustthereforebeherebymarked“advertisement”
Tu¨bingen,Germany(E-mail:[email protected]). inaccordancewith18U.S.C.Section1734solelytoindicatethisfact.
www.jn.org 0022-3077/08$8.00Copyright©2008TheAmericanPhysiologicalSociety 1461
1462 RASCH, GRETTON, MURAYAMA, MAASS, AND LOGOTHETIS
manner. We use methods derived from the field of machine cluding drives, holder, and preamplifier, was custom designed to
learning in an attempt to infer exact spike timings from the minimizecrosstalkofsignalsbetweenelectrodes(typically(cid:5)1ppm).
underlying LFPs. We compare the accuracy of spikes trains Thesignalswereamplifiedandfilteredintoabandof1–8kHz(Alpha
predictedbysupervisedlearningalgorithmsonawiderangeof Omega Engineering, Nazareth, Israel) and then digitized at 21 kHz
recordingsfromtheprimaryvisualcortex(V1)aswellasfrom with 16-bit resolution (National Instruments, Austin, TX), ensuring
the lateral geniculate nucleus (LGN) of anesthetized and non- enoughresolutionforbothlocalfieldandspikingactivities.Binocular
visual stimulation was provided through a two-fiber optic system
anesthetized macaques and investigate what kinds of features
(Avotec, Stuart, FL) after fine alignment to each of the animal’s
of the LFP are important for inferring spikes from LFP.
foveas by a modified retinoscope coupled with a stimulus projector
holder.
METHODS Inthecaseoftheanesthetizedanimalswedifferentiatebetweentwo
different conditions: spontaneous activity (“spo”) and movie-driven
Data acquisition
data(“stm”).Intheformertheinputscreenisblankforabout5min;
Electrophysiologicaldatarecordedfromnineanesthetizedandtwo inthelattera4-to6-minsegmentofacommerciallyavailablemovie
non-anesthetized monkeys (Macaca mulatta) are included in the isshown.Movieframesweresynchronizedwiththerefreshrateofthe
present study. All animal experiments were approved by the local monitor(60Hz,twosynchronizationspermovieframe)andcovered
authorities (Regierungspraesidium) and are in full compliance with 7–12° of the visual field. Most of the electrodes were confirmed to
theguidelinesoftheEuropeanCommunity(EUVD86/609/EEC)for haveareceptivefieldwithinthemoviepresentationarea(seeFig.1B
the care and use of laboratory animals. Surgical procedures are foranexample).Multipletrialsofmoviepresentationsandspontane-
described elsewhere (Logothetis et al. 2002) together with hardware ous activity are run within one recording session (intermingled with
detailsoftherecordingsetup. recordingsofotherstimulinotconsideredhere).Fromthesedatawe
Toperformtheneurophysiologicalrecordingsinanesthetizedmon- include1,304recordedtimeseriesinthepresentstudy,whichwecall D
o
keys, the animals were anesthetized [remifentanil (typically 1 trialsthroughoutthisarticle.Thedatasetconsistedofrecordingsfrom w
(cid:2)g(cid:1)kg(cid:2)1(cid:1)min(cid:2)1)], intubated, and ventilated. Muscle relaxation was nineanimalscollectedin12recordingsessions.Fromeachsessionwe nlo
achievedwithmivacurium(5mg(cid:1)kg(cid:2)1(cid:1)h(cid:2)1).Bodytemperaturewas takefiverepeatsofmoviepresentationandfiverepeatsofspontaneous a
d
keptconstant,andlactatedRingersolutionwasgivenatarateof10 activity trials (with the exception of j97nm1, where only one movie e
d
mmol(cid:1)nkkge(cid:2)y1(cid:1)ahn(cid:2)d1.thDeudreinpgththoef aennetisrteheesxiapewriemreenct,onthtienuvoiutasllysimgnosniotofrtehde. terlieacltirsodaevacihlaabnlnee).lsTpoearvsoeisdsiaonnyasruebijneccltuivdeeds.elTehctiisornesbuialtssailnlm67e0asturiraelds from
Drops of 1% ophthalmic solution of anticholinergic cyclopentolate for spontaneous activity and 634 for movie stimulus recorded using jn
hydrochloridewereinstilledintoeacheyetoachievecycloplegiaand 134 electrode placements. Movies are identical within a session but .p
h
mydriasis. Refractive errors were measured and contact lenses [hard maydifferbetweensessions. y
s
poly(methylmethacrylate)lensesbyWo¨hlkContactlinsen,Karlsruhe, In three sessions (two anesthetized monkeys) no more than four io
Germany]withtheappropriatedioptricpowerwereusedtobringthe electrodes were simultaneously placed in the LGN. Thus stimuli lo
g
animal’seyeintofocusonthestimulusplane.Simultaneousrecording reflected in these data (100 trials) are exactly identical to those in y
.o
of neural activities were made from the primary visual cortex (V1) correspondingrecordingsoftheV1data.Thisdataisseparatedagain rg
using8–16electrodesconfiguredin4(cid:3)4or2(cid:3)8matricesinagrid according to the conditions spontaneous activity (“spo(L)”) and o
of1–2mm.Electrodetipsweretypically(butnotalways)positioned movie-driven activity (“stm(L)”). Data from non-anesthetized mon- n A
intheupperormiddlecorticallayers.Theimpedanceoftheelectrode keysaremorelimitedandwereincludedinthepresentstudyonlyto p
varied from 300 to 800 k(cid:4). In the case of simultaneous LGN corroborate results for the data described earlier. They are recorded ril 1
recordinganadditionalsetofdrives,usuallyconsistingoftwotofour usingonetofourtetrodesfromchronicimplantspenetratingV1(for 4
electrodes, was additionally positioned. The electronic interface, in- a detailed description of surgical methods and recording setup see , 2
0
0
8
A B Receptive fields
600
z]
H
pike rate [ 240000 °ual field [] −05
s s
st. Vi
n 0 −10
I
5 −5 0 5
Visual field [°]
s]
nit
u
D 0 V1 electrode
S
P [ LGN electrode
F Stimulus region
L
−5
25 30 35 40 45
Time [sec]
FIG.1. A:representativeelectroderecordingfromSessiona98nm5ofananesthetizedmonkey.Topplot:theinstantaneousfiringrateofanexperimentduring
moviepresentationfromaprimaryvisualcortex(V1)electrodeinasmalltimeregion.Moviepresentationstartsafterablankperiodat30s.Recordingtime
of170-sdurationstarting5saftermovieonsetisusedforpredictionperformanceevaluationandiscalledatrial(seeMETHODS).Bottomplot:thelocalfield
potential(LFP)tracecorrespondingtotheV1electrodeabove.B:thearrangementofreceptivefieldsrelativetothemoviepresentationareaforSessiona98nm5,
wheresimultaneousV1andlateralgeniculatenucleus(LGN)recordingsareavailable.Othersessionshavesimilarelectrodearrangements.
JNeurophysiol•VOL99•MARCH2008•www.jn.org
INFERRING SPIKE TRAINS FROM LFPS 1463
Tolias et al. 2007) in a total of 56 trials. Unlike the data from The sampling interval (cid:8) is 5 ms, in accordance with the sampling
anesthetized animals the stimulus conditions here are mixed, with frequencyoftheLFPsignal(200Hz).
spontaneousactivity(notask)andafixationtaskshowinggratingsof
differentorientations.Thislastdataislabeled“awake”inthefollow- LEARNING ALGORITHMS. A support vector machine (SVM) was
ing. All data are processed in the same way as outlined in the next used (Vapnik 1999) as our learning algorithm. For a more detailed
section. introduction to SVMs see, for example, Bishop (2006), Burges
(1998), and Scho¨lkopf and Smola (2002). SVMs perform binary
PROCESSING. Thedatapreprocessingstepsareasfollows.Electrode classificationinasupervisedfashion.
signals were decimated to 7 kHz. Spiking activity is inferred from Briefly,themodelcanbestatedasfollows(fordetailsseeBishop
high frequencies of the resulting signal (see following text). The 2006)
recordinghardwareintroducesahigh-passfilterwithacutoff(cid:1)1Hz;
1Hzisthusthelowestfrequencyconsideredhere. h(cid:10)x(cid:11)(cid:4)wT(cid:12)(cid:10)x(cid:11)(cid:5)b (1)
The7-kHzsignalislow-passfilteredwithacutofffrequencyof250
Hzandresampledfirstto500Hzforcomputationalconvenience.The whereonelooksforthedecisionboundaryorweightvectorw;bisa
resulting signal is low-pass filtered at 90 Hz to derive local field biastermand(cid:12)isaprojectionintoaspaceoffeatures.SVMschoose
potentials (LFPs). For low-pass filtering we use a custom finite thehyperplanethathasthewidestmarginbetweenbothclassesrather
impulseresponsefilter(FIR):aKaiserwindowFIRfilterwith60-dB thananarbitraryseparatinghyperplane.Thisisachievedbyenforcing
attenuationinthestopband,a0.01-dBpassbandripple,andatransi- appropriate constraints in the optimization. For nonseparable prob-
tion band of 1 Hz. To eliminate possible phase shifts, signals were lems,suchasourreal-worlddata,oneintroducestheconceptofsoft
filteredforwardandbackward(usingtheMatLabfiltfiltfunction).The margins, i.e., in the optimization one now allows for incorrectly
signalisthenresampledtoafinalsamplingrateof200Hz. classified examples, where an additional parameter C regulates the
GoodpropertiesofFIRfiltersarewonattheexpenseoflargefilter penalty.
sizes(afewseconds).However,sincewediscardleadingandtrailing SVMs have the power to do nonlinear separation (seen from the D
pcoonrtcieornnshoefr(cid:6)e.15sofeachtrial,filteron-andoffsetartifactsareofno pimerpslpiceicttliyvedeofifniensputhtespfeaacteu)rebymcahpo(cid:12)os.inHgeraeninapnpornolpinreiaatrerkaedrinalelbtahsaist own
function(RBF)kernelsareused. lo
a
SPIKEEXTRACTION. Spiketimesaredetectedbyapplyingathresh- AsasimplealternativetoSVMsweusedstandardlinearregression d
e
old to the high-pass filtered 7 kHz signal described earlier (fourth- (withaconstantbiasterm)onthelabelvectorandthesamples(see, d
orderButterworth,cutofffrequency500Hz).SincethisMUAsignal e.g.,Bishop2006).Briefly,usingalinearmodelh (x)(cid:7)wTx (cid:9)b fro
isusuallyasymmetric,thedetectionthresholdisautomaticallyapplied wecalculatedtheoptimalweightvectorw*bymirneigmiizingthemi ean m
to that side where spike waveforms exhibit greater deflection. To squarederror(cid:13)[h (x)(cid:2)y]2(cid:14)onthetrainingsamples.Classlabelson jn
avoid dependence of the size of spikes the threshold is applied at thetestsetwereroebgtaiinedbiythresholdingwiththesignfunction,i.e., .ph
3.5SD ((cid:3)) of the “noise component” of the MUA signal. (cid:3) is y (cid:7)sign[h* (x)]. ys
estimated by calculating the SD of the signal, neglecting the 4.55% j reg j io
((cid:6)2(cid:3))absolutehighestvaluesdividedbythepercentageofvariance EXTRACTIONOFLFPFEATURES. AnLFPfeaturecouldbeanyaspect log
y
thatiskeptingeneral,whensettingtheprobabilityofabsolutevalues oftheLFPthatonemightdeemhelpfulforinferringwhetherthereis .o
(cid:6)2V(cid:3)istuoazleirnos.pection confirms that spikes are well detected. If the aressppeikcet atot tti.i,Initsouproawnearlyastisd,iwffeeruesnetdfrtheqeuLeFncPieast,dainffdertehnet lpahgassweiothf onrg
assumption of a Gaussian “noise component” is correct, then the oscillationsatparticularfrequencies(alsoatdifferentlags).Multiple A
rNfwareootermteeodtfmhowanutelrtotiwhnpeigltehlryenssdieuneulgttreiloencn-gttseipsdp(essilekpeeeickttrDerosIadSiieCnssUsaSwwSbeIoiOluNdlto)m.1no.B6osettHculaizsukese(eflaoynmryin(cid:3)okcsit(cid:7)lnudrde3eco.o5fsrpdsSipiknDikge)ess. fmieneaaetcotahuntredeasisnmtdaaernunednsiisttoirmSnaDipnol.ifynFgtehcaoestenurtcrseea.sstuealntriaentegedxstairnamctptheldeespsar{imxoirp}tloeisdvnievocirtdmoirnaglxiitz.heeNdsotatomezptehlreaost pril 14, 20
sorting. feaItfugre(tim)raeyprbeesednetfisnthede(ansoTrm(atl)iz:e(cid:7)d)gv(tolta)g,ewahtesraem(cid:6)pl(cid:7)ingk(cid:8)binretpi,raesteimntes 08
k i i(cid:9)k
the time lag (we neglect boundaries to simplify the description).
Learning to infer spike trains
FeaturesT (t)simplyrepresenttheLFPtimecourserelativetosample
k i
The learning algorithm has to learn to map from LFP waveforms timet.
i
(or other LFP features) to spikes. Ideally, the learning algorithm Another type of feature, which we denote P (t), is the estimated
f,k i
should output all predicted real-valued spike timings at once if the powerat(center)frequencyfoftheLFPtimecourseattimet .To
i(cid:9)k
LFP time course is given as input, although this task requires too obtain an estimate for the power at a given frequency and time, we
muchdata.Instead,wesimplifythetaskbyassumingthatspikesare calculatedthespectrogram,usingthemultitaperapproachintroduced
independent and that the spike-to-LFP relationship remains constant by Thomson (Jarvis and Mitra 2001; Percival and Walden 2002;
over time. With these assumptions one can use a binary classifier, Thomson1982).Becausespikesaresingleeventsonthetimescaleof
which yields the prediction of a spike (or no spike) at time t. (cid:8), we chose a high temporal resolution at the expense of frequency
Concatenating the prediction for each t results in a predicted spike resolution. We set the moving window to 150 ms and the time–
trainforagivenLFP.Notethattheindependenceassumptiondoesnot bandwidthproductto1.6,whichimplicitlysetsthehalf-bandwidthto
implythatpredictedspiketrainsarenecessarilyuncorrelatedbecause W (cid:7) 10.67 Hz. Spectral estimation was averaged over two Slepian
temporalcorrelationcanbeinducedbytheunderlyingLFPs. tapers. Because this window setting does not allow accurate power
In supervised fashion the binary classifier is trained on a set of estimation(cid:1)20Hzweusedlargerwindowsforbands(cid:1)20Hz(500
trainingexamplesandtestedonadistincttestset.Wetrainabinary ms)and(cid:1)6Hz(2s).Thisreducedthehalf-bandwidthto3.2and0.8
classifier on LFP features summarized in the sample vectors x, i (cid:7) Hzforfrequencies(cid:1)20and6Hz,respectively.WealsotriedMorlet
i
1,...,L,topredictthelabely (cid:1){1,(cid:2)1}.Theindexiistheithpoint waveletswithvariablebandwidthperfrequency,butthisdidnotalter
i
of the discrete LFP time series with sampling frequency 1/(cid:8) at predictionperformance.
recordingtimet (cid:7)i(cid:8)(cid:9)t .Thusy (cid:7)1statesthatthereoccurs(at To access phase information at particular frequencies of the LFP,
i 0 i
least)onespikewithintimebint andy (cid:7)(cid:2)1indicatesnospike.In we first band-pass filtered the recorded signal with the FIR Kaiser
i i
this framework prediction is temporally restricted to the sampling filter(describedearlier)withabandwidthof2(Fig.6)or4Hz(Fig.
resolutionoftheLFPs,makingitnecessarytobinthespiketimings. 7), and then used the Hilbert transform to extract an instantaneous
JNeurophysiol•VOL99•MARCH2008•www.jn.org
1464 RASCH, GRETTON, MURAYAMA, MAASS, AND LOGOTHETIS
pfehaatsueres(cid:7)f((cid:8)tif),k(atit)f:r(cid:7)equceonsc[y(cid:7)ff(.ti(cid:9)Frko)]mhathveinseg saiglnaaglsofwe(cid:6)d(cid:7)efikn(cid:8)ed. Tphheassee syn(cid:10)F1,F2(cid:2)S(cid:11):(cid:7)I(cid:10)S;F1,F2(cid:11)I(cid:10)(cid:11)S;IF(cid:10)S;,FF1(cid:11)(cid:11)(cid:11)I(cid:10)S;F2(cid:11) (4)
features have phase information identical to that of the band-passed 1 2
signalsbutaredevoidofanyamplitudemodulation(AM).Addition- This measure ranges from (cid:2)1 in the case of completely redundant
ally,weused(cid:8)ˆf,k:(cid:7)sin[(cid:7)f(ti(cid:9)k)]inthefeatureanalysisofFig.6to information to 1 for completely synergic information. The measure
help the classifier linearly extract phase locking at phases where the syn(F ,F (cid:2)S)iszeroifbothfeaturesF andF conveyindependent
1 2 1 2
cosinewouldbenearthezerocrossing. informationaboutthespikesS.
To analyze prediction accuracy on different timescales, we used
PERFORMANCEMEASURES. Thekappameasurewasusedasamea- spectral coherence (Jarvis and Mitra 2001). Spectra were again esti-
sureofperformance(Cohen1960).Letp bethefractionofsamples
l,r matedviaamultitaperapproachdesignedforpointevents(Jarvisand
havingtargetlabell(cid:1){(cid:2)1,1}andpredictedlabelr(cid:1){(cid:2)1,1}and Mitra 2001). Here the time–bandwidth product was set to TW (cid:7) 3
letqlandq˜rbethefractionofsamplesinthetestsetthathavethelabel usingtheaverageofK(cid:7)5tapers,yieldingahalf-bandwidthofW(cid:7)
l in the target or the label r in the prediction, respectively. Then the 0.001HzforT(cid:7)17s.
chancelevelforclassificationisgivenby(cid:9) (cid:7)q q˜ (cid:9)q q˜ .Ifwe
c (cid:2)1 (cid:2)1 1 1
define (cid:9)0 (cid:7) p(cid:2)1,(cid:2)1 (cid:9) p1,1 to be the overall fraction of correctly PERFORMANCE EVALUATION. We evaluated the prediction perfor-
classifiedsamples(bothpositiveandnegative),then(cid:10)isgivenby mance for each trial separately, using 10-fold cross-validation. We
analyzed a 170-s region, avoiding the on- and offset of the movie
(cid:9) (cid:11)(cid:9)
(cid:10):(cid:7) 0 c (2) stimulus. Spontaneous activity trials are also restricted to 170-s
1(cid:11)(cid:9) duration. In the case of tetrode recordings, performance is estimated
c
astheaverageperformanceofthefourwiresofthetetrode.
Thismeasureisanormalizedabove-chanceclassificationrate.Itcan HyperparametersfortheSVMalgorithmwereestimatedasfollows.
beeasilyseenthat(cid:10)equalszeroifpredictionisatchancelevel(i.e., The RBF kernel width (cid:9)was selected by a heuristic procedure. We Do
(cid:9)(cid:9)00(cid:7)(cid:7)(cid:9)1c)). andequalsoneifthepredictedclassificationisperfect(i.e., tdoiostkan(cid:9)cteosbine1th.7e7tr(aoirn3in.5g4s)etti.mFeosrtheeacmhetdriiaalnwdiestcahnoceseotfhaaltlCEu(calniddea(cid:9)n) wnloa
AnotherperformancemeasureistheSpearmanrankcorrelationr(cid:3) showingthebestperformance(averagedover10cross-validationruns d
e
betweensmoothedpredictedspiketrainsandtargetspiketrains.This on a logarithmic grid of 25 values from 0.25 to 400). We visually d
might give a more intuitive picture of the prediction quality. If not confirmedthatthisrangeisappropriateforourdata(notshown).We fro
stated otherwise, spike trains are smoothed by a Gaussian kernel of used the libSVM library (http://www.csie.ntu.edu.tw/(cid:1)cjlin/libsvm/) m
width25ms. forallSVMcalculations. jn
Yet another measure for prediction quality is the mutual informa- Since the sample sizes were heavily biased toward the negative .ph
tionbetweenclasslabels.Themutualinformation(MI)betweentarget (nonspiking)classwerandomlypickedapproximatelythesamenum- ys
spikes S and the prediction outcome of a classifier C(F) using LFP ber of samples of both classes from the training region. This effec- iolo
featuresFandlabelsL:(cid:7){(cid:2)1,1}is tivelychangesthelossfunctionfromequallosstohigherimportance g
y
I(cid:15)S;C(cid:10)F(cid:11)(cid:16)(cid:4) (cid:1)l(cid:1)L r(cid:1)(cid:1)L PI,rlog2qpllq,rr (3) frccoaloatrness)ss.tepasWinkte(ecoslrua(ssatesbhdoerau1ttm,i0oaa0)x0fiaimnvadenufdmeoml1d,pa2ivi0nrai0ccirlaaeslbaalymsleefp,oliduenesnpdtfheotenrhdistisrpnatigoikniobinnneggatahrngeedogmoinodeonancnsowpmfiiiktprhiirnnogga- on Ap.org
wherewetaketheprobabilitiesdefinedearlier.ThisestimationofMI m(cid:6)i1s,e00b0etswameepnlepsreodnilcytiomnarqguinaaliltlyy aimndprocovmedputhtaetiroensaulltssp(eneodt sbheocwauns)e. ril 1
is different from nonparametric approaches in that it can access Trainingtheclassifieronallpossiblesampleswasprohibitivedueto 4, 2
dependencethatisonlyinreachoftheclassifier;thusonehastomake the enormous sample size. We tried to use class biasing in the C 0
0
sure that the classifier captures the main aspects of its dependence. parameters(Musicantetal.2003),butthisonlyincreasedcomputation 8
Notethatweusethenaiveestimatorformutualinformation[without timewithlittlegaininpredictionquality.
biascorrection(Panzerietal.2007)].SinceallMIvaluecalculations The test set was always a temporally contiguous region to avoid
involve an identical number of bins—that is, two, one for each featurecorrelationoftrainedandtestedsamplesthatmightlienearby
class—we can nevertheless safely compare results even for classifi- intimeifarandomizedsetofsampleswereused.
cationswithdifferentnumbersoffeatures.However,theabsoluteMI
FEATURESELECTIONALGORITHM. Wenowdescribehowtodeter-
valuesmightbebiased.
mine the usefulness of different features for spike prediction. Al-
To access redundancy, synergy, and independence of information
thoughfeaturesimportantfortheSVMclassifierarehardtointerpret,
(Polaetal.2003;Schneidmanetal.2003)conveyedbytwofeatures
inthecaseoflinearregression(withsquaredloss)onecanderivean
F andF aboutthespikingactivityS,weestimatemutualinformation
1 2 analytical expression for the prediction error of a set of features
using two classifiers trained on features F and F individually,
1 2 involving only the STAs and correlation among features. Based on
yieldingI(S;F ):(cid:7)I[S;C(F )]andI(S;F ):(cid:7)I[S;C(F )].Thenathird
1 1 2 2 thispredictionerror,wederivedanalgorithmthatforwardlyselectsa
classifier is trained on both features jointly, yielding I(S; F1, F2) :(cid:7) small subset of features out of a much larger pool of features. As
I[S;C(F1,F2)].Ifbothfeaturescarriedindependentinformationfrom subsequently explained, the selected subset will show minimal pre-
theperspectiveoftheclassifier,bothfeaturestogetherwouldconvey diction error compared with other subsets with the same number of
identical information as individual features; i.e., I(S; F , F ) (cid:7) I(S; features.Inthatsensethesubsetoffeaturesselectedbyouralgorithm
1 2
F )(cid:9)I(S;F ).Ifbothfeatureswererelatedbyaone-to-onemapping represents the most useful features from a given pool. Because this
1 2
(completely redundant information about the spikes), then all terms feature selection can be efficiently done for huge feature pools, we
wouldbeequal:I(S;F )(cid:7)I(S;F )(cid:7)I(S;F ,F ).Ifthetwofeatures restrictedfeatureanalysestolinearclassification,ratherthanusingthe
1 2 1 2
didnotcarryinformationindividually,i.e.,I(S;F )(cid:7)I(S;F )(cid:7)0,but SVM classifier. This is not too restrictive in our case because linear
1 2
carriedinformationtogether,I(S;F ,F )(cid:6)0,theywouldbetermed classification achieves almost 90% of the performance of an SVM
1 2
(completely) synergistic. Thus we define a normalized degree of classifier(seeRESULTS).
synergy of information about the spikes (as measured by the classi- In linear regression we look for the weight vector w that has
ficationalgorithm)as(Schneidmanetal.2003) minimalerrorinameansquaredsense
JNeurophysiol•VOL99•MARCH2008•www.jn.org
INFERRING SPIKE TRAINS FROM LFPS 1465
(cid:12)(cid:10)w(cid:11)(cid:4)(cid:13)(cid:10)y (cid:11)wTx(cid:11)2(cid:14) (5) investigatewhichLFPfeaturesareimportantfortheprediction
i i
task and which aspect of the spikes they code for.
wherethebracketsindicateaveragingoverallsamplesi.Minimizing
theerrorisstraightforwardandresultsintheoptimalweightvector1
Average spike to LFP-power and -phase relation
w*of
w*(cid:4)(cid:13)xxT(cid:14)(cid:2)1(cid:13)yx(cid:14) (6) Figure1showsspikingactivityand(normalized)LFPsofa
i i i i representative electrode. Relationships between spiking activ-
provided that the estimated correlation matrix A: (cid:7) (cid:13)xxT(cid:14) has full ity and underlying LFPs are visualized in Fig. 2.
i i
rank.Wenotethatwehaveabinaryclassificationandthusyi(cid:1){(cid:2)1, InFig.2Athespike-triggeredaverage(STA)oftheLFPsof
1}.Thus anexampleelectrodeduringmoviestimulusisplotted.Clearly,
thereisalinearrelationbetweenspikesandLFPs.Onenotesa
c c
(cid:13)yx(cid:14)(cid:4) (cid:11) (cid:2)1 m (cid:5) 1 m sharp negativity at spike position at zero time lag and a
i i c (cid:5)c (cid:2)1 c (cid:5)c 1
1 (cid:2)1 1 (cid:2)1 prominentupswingforpositivetimelags,i.e.,afterspikinghas
wherem andm aretheclassmeansforthenonspikingandspiking occurred. Likewise, in the STA of the spectrogram (Fig. 2B)
(cid:2)1 1
class,respectively,andc andc arethenumberofsamplesineach power is enhanced in the high-frequency range (40–90 Hz)
(cid:2)1 1
class. With the condition that each feature is normalized to variance during spiking activity. Enhancement of power at high LFP
oneandzeromean,Eq.6reducesto frequencies as a response to spikes is common among elec-
trodes,stimulusconditions,andmonkeys,aswewillseeinthe
w*(cid:4)cA(cid:2)1m (7)
1 next sections.
Figure 2C shows the probability of spiking activity at the
with
oscillation phase of a particular LFP frequency for the same D
o
2c example electrode, but averaged (cid:6)30 repeats of the movie pre- w
c:(cid:7)c (cid:5)1c sentation((cid:5)120minofrecordingtime).Onenotesthatthephases nlo
1 (cid:2)1 a
of all LFP frequencies are at least weakly related to spiking d
e
Theminimalerroristhengivenby activity(Raleightestofnonuniformangulardistribution).Most d
(cid:12)(cid:10)w*(cid:11)(cid:4)1(cid:11)c2mTA(cid:2)1m (8) strikingly,spikesarerelativelytightlylockedtophasesoflow from
1 1 frequencies ((cid:13)10 Hz). The generality of this behavior is jn
Inotherwords,theerrorofthelinearregressionisdependentonlyon illustrated in Fig. 2D, where the phase preferred by spikes is .p
h
the STA m1 and the correlation among features. We note that if we plotted as an average across all data from V1 (anesthetized ys
restrictourselvestotheuseofnfeatures(thedimensionsofxi)outof animals). The average preferred phase shifts with frequency io
apoolofN(cid:6)nfeatures,theprecedingequationsremainvalidifthe fromtheonsetofapositivehalf-wavetothevalleyoftheLFP log
correlation matrix and the mean are also restricted to these features y
only(thatis,Ahassizen(cid:3)nwithranknandm1isn-dimensional). ocoscnislilsatteionnt, (rceogmarpdalerseswofithwhFeigth.er2Ea)c.tivTihtyisisbeshpaovniotarneisouvseroyr .org
aalrggTomoritashxemle.(cid:2)cWt(maes)setta(cid:2)ro(tftwhneitihmvatphroieartnfaecnaettuforefeatetuharacethshwfaeseattuhusereehtihbgeehiefnosgtlloSnwToArinm,gai.lietiez.,reafd1tiv(cid:7)toe mquoevnicey-ddrievpeenn.deFnocresiosmsleigehltelyctrdoidffeesretnhte(apsreinferthreedexpahmaspel-efroef- on Ap
osenaer)c.hAtshjsruomugeh1njaolwlNth(cid:2)atnn(cid:9)(cid:2)11rfeematauirneisngarfeeaatlurereasdyansdelcehcoteods.eTohneenthwaet Ftoi-gs.p2ikCefroerlahtiiognh fsreeeqmusentcoiebse).mFiorrroarefdewate(cid:14)lec(tnrootdeshsotwhen)p.hase- ril 14
minimizes the error (Eq. 8), where the restriction of the correlation The gray-shaded area in Fig. 2D shows the (average) phase , 2
matrix is now enforced on the n features rather than on n (cid:2) 1 (and rangewithinwhich50%ofthespikesfall.Itwouldbezerofor 00
8
analogouslywithm ).Westopthisiterationwhenthedesirednumber perfect phase locking and (cid:14)for no phase–spike relation. One
1
offeaturesm(cid:1)(cid:1)Nisselected. notes that this range is somewhat smaller for low frequencies
Thisalgorithmishighlyefficientinfindingagoodsetoffeatures, (0.85(cid:14)), but approaches 0.98(cid:14)for frequencies (cid:6)20 Hz, indi-
since we need to calculate the correlation matrix only between the
catingthatthephaselockingisfarfromperfectatallfrequen-
selectedfeaturesandallotherfeatures(whichcostsmuchlesseffort
cies and is especially weak for high frequencies.
thancalculatingitforallpairs).
In summary, we have seen that there is indeed a consistent
relationship between LFPs and spiking activity on average. In
RESULTS thenextsectionsweasktowhatdegreeitispossibletoexploit
these relations (and maybe other information available in the
This section is organized as follows. After showing the
LFP) in a systematic way to infer spikes from the LFP.
general spike-to-LFP relationships present in our data, we
report the population performance for the task of predicting
spike trains from LFP, focusing first on data from V1 of Population prediction performance
anesthetized monkeys (nine monkeys) collected during the
presentationof5minofcommercialmoviestimuliandequally In Fig. 3 a typical example of a predicted spike train is
long periods of spontaneous activity. We then compare these depicted together with the used LFP features. Figure 3, A and
resultswithamorelimitedsetofsessionsrecordedfromV1of Bshows8sofLFPspectrogramandthetimecourseinthetest
awake monkeys (two monkeys) and with results on data from region, respectively. Small vertical lines in Fig. 3B indicate
LGN of anesthetized monkeys (two monkeys). Finally, we spike times before binning to 5-ms resolution. Several inter-
estingpointscanbenoted.AsexpectedfromtheLFP-to-spike
relations discussed in the last section, spikes preferentially
1Note that the optimal weight vector w* has only minimal error for the
occurintheupswingandvalleysofverylowandmediumLFP
regression,andthattheremaybeabetterweightvectorforclassification.We
neglectthishereforthesakeofsimplicity. oscillations, as seen for instance at times 171 and 173.4 s in
JNeurophysiol•VOL99•MARCH2008•www.jn.org
1466 RASCH, GRETTON, MURAYAMA, MAASS, AND LOGOTHETIS
A B
0.4
LFP−STA 80 N
0.6
o
U] Shuffled spikes SD rm
0.3
d LFP [SD 00..24 ency [Hz] 4600 0.2 alized LFP
Normalize −0.02 Frequ 20 00.1 power [S
D
−0.4 −0.1 U
]
0
−1 −0.5 0 0.5 1 −2 0 2
Time lag [sec] Time lag [sec]
C D E
2.3
6 6
5 2.2 S 5
d] p d] D
P [ra 4 2.1 ike p P [ra 4 own
Phase of LF 23 2 robability [% Phase of LF 23 loaded from
1 1.9 ] 1 movie jn
spont .p
h
0 1.8 0 y
s
0 20 40 60 80 100 101 102−1 0 1 io
Frequency [Hz] Frequency [Hz] log
y
aneFsIGth.e2t.izedSpmikoen–kLeFyPs(rDelaantidonEs)h.iAps:sfpoirkoe-nterigegleecrterdodaeveirnagVe1(SoTfAan)LanFePs.thFeotrizseigdnmifiocnaknecyeldeuvreilnsginmteorvspieiksetiimntuelravtiaolsn((IAS–IsC))araendshaucfflroesdsaanlldrtehceoSrdDinogfsthfreormesuVl1tinign .org
o
STAiscalculated.B:STAoftheLFPspectrogram(seeMETHODS),withpowerseriesnormalizedtozeromeanandunitSD.Powerathighfrequenciesisclearly n
modulatedbyspikingactivity,whereaspoweratlowerfrequenciesshowsonlydiffusedependenceonspikes.C:probabilitydistributionofLFPphasesata A
p(saarmticeuelalercstproikdiengaspboesiftoioren)..LNFoPtepthhaastetshaerceocloomrmpuatpedshvoiawHsoilnbleyrtatrnaanrsrfoowrmra(n1g-Hezofbpanrodbsa).bHilietireesaallndsptihkaetsvoavleures30abroevpeeaotsrobfelmowovtihee-dlriimveitnsaacretivtrituyncaareteidn.cBluldacekd pril 1
dotsindicatethepreferred(i.e.,mean)phase.Nophaselockingofspikeswouldresultinauniformdistributionat2%perbin.Althoughlockingtolowfrequencies 4
isstrong,lockingtohighfrequenciesisonlyweak(butpresent).D:theaveragepreferredphasesforallelectrodesacrossallanesthetizeddataindividuallyfor , 2
0
movie-drivenactivity(“movie”)andspontaneousactivity(“spont”).BarsindicateSEs.Thephaserangecontaininghalfofthespikesaroundthepreferredphase 0
8
isindicatedbytheshadedarea.E:theinterpretationofphase,showingthatspikesarelockedatverylowfrequenciestotheonsetofapositivehalf-waveand
athighfrequenciestothevalley.
Fig. 3B. Additionally, the power of multiple frequencies is closelypredicted(172.5s),althoughnoclearmarkintheLFP
enhancedwhenaburstofspikingactivityoccurs,assuggested time course or spectrogram can be seen with the naked eye.
byFig.2B,butthefrequencyresponsetoburstsisdiffuseand However, their length (176.6 s) and exact position (176.3 s)
variable (compare the burst at 173.4 s to that at 174.5 s). The sometimes seem inaccurate. There are also occasions where
clustering of spiking activity on a timescale of a few hundred spikes are simply missed (172.9 s) or fabricated (173.9 s).
millisecondsinthisexampleisactuallyquitetypicalinourV1 Prediction performance is evaluated in different ways. One
data (see following text). For the single spikes in between the measure is the (cid:10) performance, a measure defined on the
clusters no feature of the LFP is immediately predictive. samplesinthetestset,whichispositiveforabove-chanceclassi-
Figure 3C shows target and predicted spiking activities. fication; it equals one for perfect classification (see METHODS
Prediction of spikes is made for individual sample times at a for definitions). In contrast, the correlation measure r (see
(cid:3)
resolution of 5 ms using a set of LFP time course and fre- METHODS) is defined as a local average in the time domain and
quencyfeatures(seefollowingtext).Concatenatingthepredic- is therefore less sensitive to small temporal inaccuracies. The
tion over time yields a predicted spike train that is compared performancemeasure(cid:10)ofthepredictedspiketraininFig.3Chas
with the target spike train. One notes that the prediction a value of (cid:10)(cid:7) 0.40, which is relatively good (but not the best
captures, at least approximately, the overall structure of the possible;seefollowingtext).Rankcorrelationisr (cid:7)0.60.
25ms
spiketrain.Theoccurrenceofburstsofspikingactivity,which IntheexampleofFig.3,thepredictedspiketrainresembles
are associated with easily seen traces in the LFP time course theoriginaltoacertaindegree.Weaskwhetherthisprediction
and spectrogram, is well predicted. Nevertheless, the exact qualitycarriesovertoLFPsrecordedfromdifferentmonkeys,
onsetsandoffsetsoftheburstsaresomewhatinaccurateinthe electrodes, and stimulus conditions. For that we estimated
prediction. Even some smaller bursts and single spikes are prediction performance using a large data set (see METHODS).
JNeurophysiol•VOL99•MARCH2008•www.jn.org
INFERRING SPIKE TRAINS FROM LFPS 1467
A Normalized power
−2 0 2 4 6 8 10
80
z]
H 60
y [
c
n
e 40
u
q
e
Fr 20
0
B
4
U] LFP and spikes
D 2
S
e [ 0
g
a
olt −2
V
D
−4 o
w
C n
lo
300 Target SDF (σ=0.01s) ad
kes/s] 200 Prediction (κ=0.40) ed fro
Rate [spi 1000 jn.phym
s
io
171 172 173 174 175 176 177 log
Time [s] y.o
FIG.3. ExampleofspikepredictionfromLFP(anesthetizedmonkey,Sessiona98nm5,spontaneousactivity).A:the(normalized)spectrogramofthe8sof org
LFPactivity.B:correspondingLFPtimecourseandspikingactivity.SpikesareindicatedbymarksbeforebinningtotheLFPresolution(5ms).C:thebinned n
targetspikesandtheirspikedensityfunction(blue)togetherwiththepredictedspikesandtheirspikedensityfunction(red).Thepredictionisrelativelygood A
(th(cid:10)a(cid:7)tre0g.i4o0n,sro2f5hmisgh(cid:7)a0ct.i6v0i)tyoanrethwiseltlriparle,dbiuctteodt,hwerhetrrieaalsstshheolwoceavtieonnboeftsteinrgpleersfpoirkmeasniscele(scsoamccpuarraetew.iCthlars2s5ifimcsavtiaolnueissdoofnoethweirthtrtihaelssuinppFoigrt.v4eBc,to“rspmoa”c)h.iOnene(SnVoMtes) pril 1
radialbasisfunction(RBF)classifiertrainedontheregion35–160susingthesamefeaturesasforthepopulationanalyses(Fig.4). 4
, 2
0
InferencesaremadeonthebasisofasetofLFPfeatures,with predictspikingactivity(allconditionshighlysignificant;t-test, 0
8
whichweobservedadependencebetweenspikesandLFPsin P (cid:1) 10(cid:2)6, Wilcoxon signed-rank test for zero median, P (cid:1)
the previous section. In the population analysis we include as 10(cid:2)6). Second, prediction quality for the stimulus condition
features the time course around each sample position (in a and for spontaneous activity differs only slightly: indeed, one
windowof100msbeforeand300msafterspikeposition)and cannot reject the hypothesis that the underlying distributions
an estimate of the frequency content of LFPs at zero time lag have identical means (two-sided unpaired t-test, P (cid:7) 0.21;
[Pf,0(ti);seeMETHODS],resultinginatotalof116features.This linear, P (cid:7) 0.18). However, if one compares pairwise record-
feature set generally produced good performance (with a rea- ings during spontaneous activity and stimulus presentation
sonable computational speed) over a wide range of data. For done with identical electrodes, mean and median prediction
thepredictionitselfanonlinearsupportvectormachineisused performances on spontaneous activity are significantly better
with radial basis functions such as kernel (SVM–RBF) and thanthoseonstimulus-drivenactivity(one-sidedpairedt-test:
a linear classification is employed (for details see METHODS). P(cid:1)10(cid:2)4;linear,P(cid:1)10(cid:2)4;forthedistribution-freeWilcoxon
In Fig. 4 the prediction performance over all trials is eval- matched-pairssigned-rankstest:P(cid:1)10(cid:2)4;linear,P(cid:1)10(cid:2)3).
uated(on10cross-validationruns)andaveraged.Theanesthe- Average prediction performance for spontaneous activity is
tizedV1datasetislabeled“spo”forspontaneousactivityand (cid:10) (cid:7) 0.211 (cid:17) 0.006 (linear (cid:10) (cid:7) 0.185 (cid:17) 0.005) and (cid:10) (cid:7)
“stm” for movie stimulus-driven activity. We shall focus on 0.201(cid:17)0.005(linear(cid:10)(cid:7)0.175(cid:17)0.005)forstimulus-driven
thesedataforthemoment.Theremainingconditionsshownin activity.
this plot are discussed in the following text. Third, nonlinear margin classification is consistently better
PlotAshowstheaverageperformance(cid:10)fortheSVM–RBF than linear classification (one-sided paired t-test: “stm” P (cid:1)
classifierandforlinearclassification.Fromtheresultswedraw 10(cid:2)6, “spo” P (cid:1) 10(cid:2)6; Wilcoxon matched-pairs signed-ranks
the following insights. First, since performance measure (cid:10)is test,P(cid:1)10(cid:2)6).Itamountstoanincrease,onaverage,ofabout
(cid:6)0 for above-chance prediction it can be said that both 12%inperformance.ThissuggeststhatthemappingfromLFP
classifiers can exploit information in the LFP time course to featuresisnonlinear.However,sinceasimplelinearregression
JNeurophysiol•VOL99•MARCH2008•www.jn.org
1468 RASCH, GRETTON, MURAYAMA, MAASS, AND LOGOTHETIS
A B
0.35 1
SVM RBF a98nm5
0.3 Linear regression ms] ac9988nnmm16
0.8 c98nm2
e 5 d04nm1
nc 0.25 =2 d04nm2
a σ d98at
rm 0.2 n [ 0.6 gg0927nnmm11
o o georgios3
κg. perf 00.1.15 correlati 0.4 jlrs99907772nnnnmmmm1111
Av k
n 0.2
a
0.05 R
0 0
stm spo awake stm (L) spo (L) stm spo awake stm (L) spo (L)
Condition Condition
C D
[%]mm 80 sstpmo ±± SSEE −−30..32%%mmmm−−11++4443%% 0.25 SCaromses eelleeccttrrooddees ( s(stmtm)) Dow
κκce /0 60 mance 0.2 SCSairgomnsiesfi ceealleencct ttprrooeddrfeeos r( ms(spapono)c)e nloaded
man 40 erfor 0.15 Not significant from
ative perfor 20 κAvg. p 0.00.51 jn.physiolog
el y
R .o
0 0 rg
0 2 4 6 8 10 LFP: V1 LGN V1 LGN o
n
Electrode distance [mm] Spikes: V1 V1 LGN LGN A
p
Electrode configuration ril 1
FIG.4. Population performance for spike prediction from LFP. A: average prediction performance (cid:10)for SVM and the linear regression classifier across 4, 2
conditions[anesthetizedmonkeyV1,moviestimulus(“stm”),anesthetizedmonkeyV1,spontaneousactivity(“spo”),non-anesthetizedmonkeyV1withmixed 0
stimuli (“awake”), and spontaneous activity or movie-stimulus driven activity in anesthetized monkeys from LGN, “spo (L)” and “stm (L),” respectively]. 08
Prediction is above chance level for all conditions (see RESULTS for significance tests). B: prediction accuracy of the nonlinear classifier evaluated by rank
correlationbetweentargetandpredictedspikestrainsmoothedbyaGaussiankernelofwidth(cid:3)(cid:7)25ms.Redhorizontallinesindicatetheaverageperformance
withineachconditionanditsSE.Smallblacklinesshowthequalityonindividualtrials.Insomecasespredictionyieldsveryaccurateresults,withcorrelations
ashighas0.8–0.9.Blackcurvesonthesidesofeachconditionaresmoothedhistogramsovertrials.Symbolsindicatetheaverageperformanceofindividual
sessions(i.e.,1dayofrecording).Averagesessionperformanceclustersneartheoverallmean,butvarianceforindividualelectrodesishigh.C:cross-electrode
prediction.LFPsaretakenfromoneelectrodeandspikesfromanother.Relativepredictiondecreasesmuchfasterwithincreasingelectrodedistanceinthe“stm”
condition than for the “spo” condition. A linear fit of the data points is also shown. Vertical lines indicate SEs. D: average prediction performance for
simultaneousLGNandV1recordings(3sessions,2monkeys).Performanceiscomparedforavailablecross-electrodepredictionfromdifferentareasandwith
averageperformancewhenusingthesameelectrodeforspikesandforLFPs.X-axislabelsindicatewhereLFPandspikesoriginatefrom.
classifier already achieves almost 90% of the accuracy of the Thereisnot muchvariabilityinperformance overtime: the
nonlinearclassifier,onecouldstatethattheLFPfeaturespace average SD for the (cid:10)performance of five repeats of 170-ms
exploited here seems expressive enough for this task. recordings for the same electrodes is 0.023 (cid:17) 0.002 for
We found that for individual trials performance varies stimulus-drivenactivity,0.026(cid:17)0.002forspontaneousactiv-
widely. For selected trials prediction performance can reach ity,and0.045(cid:17)0.002forbothtogether.Thisisincontrastto
(cid:10)(cid:7)0.65.PlotBofFig.4showstherank-correlationmeasure the variance across electrodes recorded simultaneously. Here
r oftheSVM–RBFprediction.Eachthinshortlinerepre- theaverageSD(in(cid:10))is0.113(cid:17)0.004forstimulus-drivenand
25ms
sents performance for an individual trial. Whereas the corre- 0.130 (cid:17) 0.005 for spontaneous activity. The roughly 25-fold
lation for some trials is as high as 0.8–0.9 on this moderately increase in variance across electrodes compared with within-
small timescale (25 ms), it is almost zero in others. There are electrode variance suggests that prediction performance is a
some trials where prediction fails completely in each of the matter of which electrode is being observed, rather than stim-
conditions.Thefailingtrialsarenotallfromthesamesessions ulusconditionortime.Electrodetipsmightbepositionedina
since the session means (markers) tend to cluster around the region where the arrangement of current sources and sinks
overall mean. mightdiffer(e.g.,indeeporsuperficiallayers),orwhereactive
JNeurophysiol•VOL99•MARCH2008•www.jn.org
INFERRING SPIKE TRAINS FROM LFPS 1469
neuronsmightbelesswellcorrelatedwiththebulkactivityof comparespontaneousandstimulus-drivenactivitybecausethe
the cortex. Since we cannot distinguish the layers from which electrode placements do not change with the condition.
electrodes record, we cannot pursue this further. Because LGN data were collected while other electrodes
Up to now we have presented results for recordings only simultaneously recorded from V1, we can investigate whether
from V1 of anesthetized monkeys. We also have a limited the LFPs of V1 can be predicted on the basis of spikes from
amount of data available from V1 where monkeys were not LGN and vice versa. This is shown in Fig. 4D averaged over
anesthetized and free to behave. The stimuli for this data set data from the three sessions recording simultaneous measure-
(labeled “awake”) are mixed and include spontaneous activity ments from V1 und LGN (see METHODS). Performance is
averaged either across electrode predictions (regardless of
andfixationtasks.Anotherpoolofdataconsistsofrecordings
distances) or over all predictions using the same electrode for
from LGN of anesthetized monkeys (labeled “stm(L)” and
bothsignals.Althoughresultsaredifficulttointerpretbecause
“spo(L)”; see METHODS).
of the limited size of the data set, one notes that using LFPs
WeseefromFig.4Athatpredictiondiffersquitedrastically
from LGN and spikes from V1 results in performance above
for the different data types. Spike prediction for the anesthe-
chance, whereas LFPs from V1 seem to hold no information
tized monkey data from V1 is more than fivefold better than
about spikes in LGN (unpaired Wilcoxon signed-rank test for
thatintheLGN,whereperformanceishardlyabovechance:on
medianperformancedifferentfromzero,significancelevel0.05).
average, (cid:10)(cid:7) 0.035 (cid:17) 0.005 (linear (cid:10)(cid:7) 0.033 (cid:17) 0.005) for
movie-driven activity and (cid:10) (cid:7) 0.027 (cid:17) 0.003 (linear (cid:10) (cid:7)
0.022 (cid:17) 0.003) for spontaneous activity. Temporal accuracy of predicted spike trains
AsinV1,thereislittledifferencebetweenspontaneousand We found an average (cid:10)value of about 0.2, which is well
movie-driven activity in LGN, although there is a reversed above chance but nevertheless far from perfect prediction at D
tendency for spike prediction to be easier on movie-driven (cid:10)(cid:7)1.Ontheotherhand,inexampleFig.3Csomefeaturesof ow
activity than that on spontaneous activity. This tendency is the target spike trains seem to be well captured by the predic- nlo
barely significant (one-sided paired t-test: P (cid:7) 0.02; linear, tion, especially regions of high and low activity, which alter- ad
P (cid:7) 0.01; Wilcoxon matched-pairs signed-ranks test: P (cid:7) nate on a timescale of about 0.5 s in this example. Thus one ed
0.05; linear, P (cid:7) 0.08). mightaskatwhattimescalethepredictedspikingactivitymost fro
Wefindthataveragepredictionperformanceonawakedata closely resembles the target spiking activity or at what timing m
is(cid:10)(cid:7)0.063(cid:17)0.005(linear(cid:10)(cid:7)0.046(cid:17)0.005).Thisismuch accuracy the prediction fails. jn
worse than that on anesthetized V1 data (unpaired t-test, P (cid:1) To answer this we evaluated the coherence between target .ph
10(cid:2)5), but still significantly better than that on LGN data (all and predicted spike train (Fig. 5). Coherence is a correlation ys
unpaired one-sided t-test, P (cid:1) 0.05). Figure 4B reveals that measure in the frequency domain. Coherence at a particular iolo
g
individual trials have a correlation of target and prediction frequency makes a statement about the exactness of the pre- y
similar to that in anesthetized monkeys. There are trials with diction on a timescale of one over that frequency. We also .org
correlationuptor (cid:7)0.6,whereasinthecaseoftheLGN, estimated the temporal accuracy directly in the time domain o
25ms n
no trial exceeds 0.3 correlation. (by varying the correlation kernel width) where one arrives at A
similar conclusions (not shown). pril 1
Cross-electrode predictions 4
1 , 2
ThevolumeofcortexthatcontributestothegenerationofLFPs 500ms 100ms 25ms 10ms 5ms sstpmo 00
8
isdifferentfromthatproducingourspikingsignal(seeINTRODUCTION). 0.8 awake
stm (L)
Thus it might be interesting to see how the relationship be- spo (L)
tween the two signals changes with distance. Because record- ce 0.6
n
ingsweredonesimultaneouslywithmultipleelectrodes(inthe e
data set from anesthetized animal), we tried to infer spikes her
o 0.4
from LFPs collected with two different electrodes. In Fig. 4C C
the average performance is plotted against the (three-dimen-
sional)distanceoftheelectrodetips.Tofacilitatecomparison, 0.2
performance is evaluated relative to the average performance
achieved using the spiking signal from the electrode from 0
which the LFPs were taken. 1 10 100
One notes that prediction performance drops to about 40% Frequency [Hz]
when electrodes are 1 mm apart (the minimal distance in our FIG.5. Coherencelevelsofpredictedandtargetspiketrains.Coherenceis
recording setup). Interestingly, for stimulus-driven activity plottedagainstfrequencyandisaveragedoveralltrialsforeachcondition.For
performance degrades significantly with distance (rank corre- comparisonweincludedsurrogatedata,inwhichspiketrainsweregenerated
withthesameISIdistributionasinthemonkeydata.Graydashedlinesshow
lation between distance and relative kappa performance using
coherencebetweensurrogatespiketrainanditsjitteredversionwithGaussian
allmeasurements:(cid:2)0.20,P(cid:1)10(cid:2)4),whereasforspontaneous noiseofdifferentSDs((cid:3)from5to500ms,aslistedintheplot).Coherence
activity no significant correlation with distance can be found dropsforhigherfrequencies,suggestingthat,onaverage,predictionisreason-
for distances (cid:13)1 cm (rank correlation 0.015, P (cid:7) 0.2). Note ably good only for slow structure in the spike trains. Note that the chance
coherencelevelisatabout0.15here,asshownbythesurrogatedata.Chance
that the number of samples becomes relatively small for
levels do not tend to zero because coherence is estimated on 10 cross-
distances(cid:6)6mmsincerectangularelectrodegridswith1-mm
validationregions(each17-sduration)andonlysubsequentlyaveragedover
spacing are used for most sessions. However, we can safely alltrials.ColoredareasindicateSEs.
JNeurophysiol•VOL99•MARCH2008•www.jn.org
