Single neuron sequential interspike intervals: Measurement, teleprocessing, and analysis

COMPUTERS

Single

AND

BIOMEDICAL

Neuron

Laboratoire

RESEARCH

13, 21l-224 (1980)

Sequential lnterspike Intervals: Teleprocessing, and Analysis G. SANDNERAND

E. DREYFUS

de Neurophysiologie, Centre 67000 Strasbourg,

de Neurochimie France

Measurement,

du C.N.R.S.,

Received May 2, 1979

Successive neuronal interspike intervals (ISI) were measured sequentially. Code numbers, each one characteristic of an experimental manipulation, were added to the measured ISI. The resulting file was sent automatically through MODEM and telephone to an UNIVAC 1110 computer by way of a microprocessor system. The UNIVAC computer was programmed to extract each code number from the IS1 file and to create a descriptive index by use of a set of questions to the experimenter, each one depending on the code number encountered in the IS1 file. Each element of the descriptive index contained, as pointer, the address of the very IS1 where the manipulation started or ended. From these data, a classical spike rate versus time diagram as well as a sequential probability density function were computed and some spike train characteristics were automatically estimated from an autocorrelogram and a serial correlogram. An example is given to illustrate the results obtained with the system.

INTRODUCTION

According to a definition taken up by Zadeh, any system “is an aggregation or assemblage of objects united by some form of interaction or interdependence” (I). With regard to the central nervous system, our specific interest lies in the study of the functional interactions occurring in neuronal assemblages thought to play an essential role in the elicitation of escape behavior. The stimulation of several brain sites such as the central gray of the mesencephalon or the medial hypothalamus is well known to elicit escape. Rats learn easily to press a lever or to cross a photobeam to switch off such a stimulation (2). In order to get direct information on the neuronal activities correlated with an escape response induced by stimulating the mesencephalic central gray or the medial hypothalamus, a study was undertaken that consisted of two experiments, namely a chronic one and an acute one, carried out in each of a number of rats (3). More concretely, a stimulation was applied 211 OOlO-4809/80/03021 I -14$02.00/O Copyri@t @ 1980 by Academic Ress. Inc. Ail rights of reproduction in my form reserved.

212

SANDNER

AND

DREYFUS

to the central gray or to the medial hypothalamus and both the escape latencies (in the chronic experiment) and the neuronal spike activity (in the acute one) were studied as functions of the stimulation parameters. Since this neuronal activity is an all-or-none phenomenon, it is easily digitized. The interspike intervals (ISI), namely the time intervals between the spike potentials, are considered to be the biologically relevant data (4). Neurophysiological findings are commonly expressed in terms of mean IS1 or mean spike rate analysis through the latter brings out only a rather partial information, since Eckhom et al. have shown the neuronal firing rate to be a rather unreliable code higher up in the brain (5). Therefore, we have developed a tool for spike train analysis based on the following principles: -Storage of each successive IS1 so that no information gets lost and that data are easily and quickly accessible to further computation. -Clear-cut specification of the characteristics of any kind of spike train, through the use of a set of computations as little redundant as possible. The initially used methodology must be briefly described before presenting the more recent steps undertaken in order to meet the above-stated principles. The hardware and the software will then be described in some detail. This methodology makes it possible to use at any time any FORTRAN written program to get information on the neuronal discharge patterns as illustrated by a short description of a set of such programs. INITIALLY

USED METHODOLOGY

We usually recorded the spike potentials of a single neuron in the brain of the curarized rat before, during and after each of several manipulations such as electrical brain stimulations, electrophoretic drug applications, or peripheral skin touches or pinches. Four channels of an analogic tape recorder (MP 5521, Schlumberger) were used. The amplified (x 1000) and filtered (500-10000 Hz F 3db) neuronal activity was recorded by means of frequency modulation on channel 3 (MA 1217/1128). The pulses of the stimulation train were recorded by means of amplitude modulation on channel 2 (MA 2219/2114). When studying the effects of local electrophoretic drug applications, the following signals were recorded in order to precisely indicate the time of injection: a + 0.5 V DC signal during the injection of glutamate and a - 1.O V DC signal during the injection of GABA were recorded by means of frequency modulation on channel 1 (MA 1217/1128). Channel 4 was used to record a set of standardized vocal comments on the stimulation parameters, the iontophoretic ejection current intensity and the handling of the animal (MA 4201/4101) (Figs. 1 and 2). The neuronal spike activity was digitized and a DIDAC 800 (INTERTECHNIQUE) computer was used to estimate the mean IS1 or the mean spike rate during the manipulations and between them. These estimations were read by the experimenter from the printer output of the DIDAC and sent manually to an UNIVAC 1110 computer through an ALPHA 20 (SINTRA) terminal for further computation.

W

Laboratoire

de

I

~S~NTRA

FIG.

Neurophysiologie

20

1

de

Calcul

1. Block diagram of the data acquisition, transmission,

ALPHA

Cmtre

c

and

de

,I,O

computation

&F

INlVAC

STRASBOIJRG

system.

COMPUTER

CRONENEIOURG

214

SANDNER

AND

DREYFUS

RECENT DEVELOPMENTS General

Methodology

In fact, the way in which the information is coded by successive IS1 is still a controverted question. If the information is being processed by means of an interval code, any ISI may be useful in attempt to describe this code and to describe thereby functional relations between neural structures (6). To this end, we have completed the above-described setup by adding a special circuitry to the DIDAC and by using a Motorola 6800 microprocessor system that would allow us to send automatically each ISI to the UNIVAC computer. The successive ISI were measured by the DIDAC computer by means of clock pulses (2 x lo4 pulses . set-l) accumulating in a memory, its address being increased by one step by the electrical pulse corresponding to the occurrence of a spike potential; a window discriminator being used to transform the spike potentials into electrical pulses. Each memory location thus contained the duration of one IS1 (maximum: 800 successive intervals up to 50 set each). The lowest possible IS1 has a duration of 5 x lo-’ see, i.e., 10 clock pulses, as there exists a “refractory period” following the occurrence of any spike potential (7). Thence, the numbers between 0 and 9 can never be ISI. In fact, these numbers were used as code symbols, each one corresponding to a specific experimental manipulation. These code symbols were automatically included in the memory file of the DIDAC at the very moment when the manipulation started or ended. This was achieved by adding to the data a pair of impulses separated from each other by various specific time intervals shorter than 5 x 10v4 sec. A special circuitry, fully described in Fig. 2 was added to the DIDAC for this purpose. The resulting data file is first stored in the DIDAC and later transmitted to the memory of the microprocessor system which is able to store up to five contents of the DIDAC. This memory was used as buffer for the automatic transmission of the data by telephone to the UNIVAC computer. The UNIVAC computer was programmed to read the data and to create both a formatted file containing the sequential ISI (in seconds) and a descriptive index based on a set of questions to the experimenter, depending on the code number uncountered in the IS1 file. Microcomputer

Hardware

Our basic material consisted of a kit microprocessor system (MEK 6800 D 2, Motorola) available on two printed circuitry cards, the first one supporting the microprocessor itself (MC 6800, Motorola), its clock generator, two parallel interface adaptors (PIA), one asynchronous serial interface (ACIA), 1 kbyte ROM, 1 kbyte EPROM, and 384 bytes of RAM memories, the other card being made of an hexadecimal keyboard and an hexadecimal LED display. The EPROM memories contained the program for the data transmission. Because of the low RAM memory capacity, we added a 16 kbyte memory card (MMS 68103, Motorola). Moreover, two apparatuses, namely the ALPHA 20 and the DIDAC, had to be connected to this system. The ALPHA 20 terminal received

NEURONAL

INTERSPIKE

INTERVAL

ANALYSIS

T -I-g x

215

216

SANDNER

AND

DREYFUS

&bit parallel ASCII symbols from the PIA of the kit card. The transmission control symbols were fed into the microprocessor system from the transmission MODEM (TRADAN 1200, LTT) through the ACIA and an impedance adaptor (LM 741 integrated circuitry). The parallel output of the DIDAC was connected to the microprocessor system by a 38-bit wiring, i.e. 36 bits for data and addresses and 2 command bits, one of them signalling the “free access” to the memory of the DIDAC, the other being used to increase this address step by step. This parallel wiring connection required the addition of a double PIA card to our microprocessor system (MEX 6820, Motorola). Microcomputer

Software

(Fig. 3)

The transmission program, stored in a I-kbyte EPROM (MCM 68708), is made of an endless loop enquiring about the DIDAC “free access” signal and about the following symbol: “f”, coming from the Centre de Calcul when the UNIVAC computer is ready to receive the data. When one of these two conditions is fulfilled, then the program branches to a subroutine. There are two subroutines for two specific tasks, namely, read the DIDAC memory content or send the Motorola memory content to the UNIVAC computer. At the end of each subroutine, the program returns to the initial waiting loop. The output of the DIDAC starts automatically at a preset memory address. The latter is then increased step by step by a command signal sent by the DIDAC reading subroutine. Each address is compared with the just preceding one, and the corresponding memory content is read only if the difference between both addresses does not differ from one. Each data made of six DCB symbols is stored in 3 bytes of the microprocessor system memories. The DIDAC reading subroutine ends when the address of the DIDAC reaches a preset value, since the DIDAC “free access” signal is then reset. The user is informed of whether the memory of the microprocessor is empty or full by a specific green or red LED accessible to the program through the control outputs of a PIA. One reading of the 800 x 36 bits lasts 65 msec. For the transmission itself, the data have first to be gathered in blocks of 72 ASCII symbols, i.e., the maximum possible with regard to the computer’s input format, since we wanted to optimize the data versus control symbol ratio, the number of control symbols being the same whether each block contains a small or a great total number of symbols. In other words, the transmission subroutine has to find the DCB numbers, two per byte, in order to convert them into ASCII symbols and to gather them into blocks. Moreover, this subroutine has to recognize and to send command symbols. A “>” symbol allows the transmission of a block of data which has to be followed by the emission of a “CR” symbol. The end of the complete transmission is signalled by an emission of six “*” symbols. The transmission of one DIDAC memory content through a full duplex 300/1200 baud telephone channel lasts about 5 min. The duplex procedure allows visual control of the transmitted data on the screen of the SINTRA terminal.

NEURONAL

EeEJ

Yes

INTERSPIKE

INTERVAL

217

No

FIG. 3. Flow diagram of the microprocessor

Software for the UNIVAC

ANALYSIS

software.

Computer

A program was written to collect automatically clear-cut information on each manipulation signalled by a code symbol in the file transmitted by telephone. This is done by branching to a subroutine programmed to put questions to the experimenter about the manipulation, when a stimulation onset symbol has just been read by the main program. The subroutine includes the resetting of an accumulator. Since the latter is used to sum the successive ISI, it measures the actual duration of the manipulation (precision & 2 ISI), its content having to be stored at the very moment when the stimulation offset symbol is read by the main program. The characteristics of the manipulation as well as its actual

218

SANDNER

AND

DREYFUS

TABLE ALGOL

PIDGIN procedure begin

ALGORITHM

I

OF THE

DATA

STORAGE

DEMAND INPUT J + 1: I + 0; K + 0; COMP t 0: wr,te ‘1 ST TIME?‘: REP = ‘YES’ then call subroutine DEBUT end if

if read

+ input

while

PROGRAM

read REP -

answer.

FIFO file:

I
= = = = =

I

6 then I then (2 OR 4) then (3 OR 5) rhea 0 then

call call call call call

subroutine subroutine subroutine subroutine subroutine

then .l+J

STIMON STIMOF IONON IONOF COMENT

~0 go go go go

- I:OUTAB(J)+5.10.“t

INTAB(I

- i)go

ro to IO IO w

suite suite suite suite suite

end end end end end

10 unfin.

i=O tlnftt: suite: deuxfin:

COMP t COMP K-1;COMP+O; I-1+1:

else OUTAB(./) + OUTAP (J); J -J

end II - 0; while II < J I do begin II-ll + 1; write OUTAB (II) in output end write “end file” in DATA FIFO file;

DATA

+ 5.10.“* INTABUJend + I go lo deuxfin:

if.

FIFO file;

end (i) Brief description of subroutine: DEBUT: finds the end of both data 6le and description tile. IONON, IONOF, SIMON, COMENT: each one puts questions tion file. (ii) Example of a data tile and its alphanumeric descriptive index: Data file: GS4302fDATA oo.an36oo OROtXMO oo.oocQ32 00.020@23

00.002804 oo.oo4o4o

to the experimenter

Description

element

to the descrip-

file: GS4302IDESCRI

1st

5kh\

. . . . OOC0.55 STIM.ON.FQT.GS4302

999999

END

Caption of Tabk I: for symbols, see Ref. (8). Gs4302: +tttificsbon code of the neuron. Fortmu in FORTRAN @asc: Data tik: (F6.9) Description

. .._.__........ 1NT:lOO MMA

._ ._,.....,... .,., IPI: 20MSEC DUR:.IMSEC

OF FILE...**.*...*.....*.*...*..*....**.*..*.......*.~.

End tile marks

(iii)

and adds.sn

file: (16, 9A6).

it: if if. if. if

NEURONAL

INTERSPIKE

INTERVAL

ANALYSIS

219

duration are stored by each subroutine in the form of an alphanumeric chain of characters, line by line, each line containing also the address of the occurrence of the manipulation in the data file. This data file is made of the successive IS1 converted into seconds. Since some IS1 were originally subdivided by the insertion of a code symbol, the inserted code symbols are as such ignored while the values located just before and after them are being concatenated (Table I). Data and description file are stored on the same digital magnetic tape. These data can be processed later on by means of a set of programs in order to get precise information about the spike train characteristics before, during, and after any manipulation applied to the brain. RESULTS Two subsets of programs were used, the first providing us with general information and the.second with answers to more specific questions. The general information consisted in classical spike rate versus time graphs and in an estimation of an IS1 probability density function (Pdf). Sequential spike rates were computed in the following way. The IS1 were read from the data file, counted and summed. Whenever the summed IS1 became a integer multiple of a given binwidth, usually 1 set, a line of “X” characters was printed, proportional in length to the counted spike number and the counter was reset. The mean and the standard deviation of these spike rates were computed and printed each time the data file address became equal to, or greater than, a file pointer in the descriptive file. The corresponding descriptive alphanumeric text was also printed. The other general information was obtained through the estimation of the Pdf. In fact, the Pdf of the IS1 could also be estimated from an histogram, but Sanderson and Kobler have proved that “a Parzen estimate using a Gaussian weighing function reduced the number of ZSZ required to achieve a given approximation by a factor of5 to 10” (9). We wrote a program to compute and plot IS1 Pdf on a Benson graphic plotter, using for that purpose medium-range parameters taken from examples given by Sanderson and Kobler. To this program was added an automatic estimation of the location of the peaks, i.e., IS1 values, and of their relative height (Fig. 4). The method used was similar to that described by Vibert and Caille who have also reviewed the literature on the estimation of normal subpopulations in an heterogeneous sample (10). The computer was not only used to describe the neuronal discharge patterns, but also to answer more specific questions. The latter may vary depending on the hypothetized neuronal information code. Let us assume that the observed neuron carries information coded by some time-dependent rythmicity. This rythmicity may be found neither by visual nor by general computational information analysis (Fig. 5). Perkel et al. have suggested the plotting together of consecutive IS1 and higher-order IS1 histograms; this was called an autocorrelogram (1 I ). By computing such an autocorrelogram from randomly shutfled data and from a raw data autocorrelogram, the difference between both

220

SANDNER

AND DREYFUS

Pd

W-

ORDER

FIG. 4. Sequential probability density function (Pdf) estimations from N successive ISI, with a 50% overlap: The method is schematically illustrated in the upper part of this figure. The medium level three graphs were Pdf, the x axis being the IS1 in milliseconds. The diagram at the bottom summarized the sequential Pdf by plotting the main ISI and their standard deviations computed by the UNIVAC along they axis (M = mean, Sd = Standard deviation) and the sequential order, before, during and after a stimulation along the x axis. The relative height of the peaks were represented by the height of the triangles whose base was made of 2 x standard deviation.

should be able to reveal the spike rythmicity. Furthermore, when trying to uncover a slow shift in ISI, or a bursting discharge pattern or an alternating of long and short ISI, a serial correlation computation was used. This was also studied by Perkel et al. (II). Sj = Cj/U*,

[II

where 0’ = . . I , - 1, 0, 1, . . .), cj = HtTi - P)(Ti+j - PII p = estimated mean, o = estimated variance, Ti = ith ISI, Sj = serial correlation coefficient ofjth order. We computed the serial correlogram from the first to the 50th order. The @I regularity was estimated by the standard

NEURONAL

INTERSPIKE

INTERVAL

ANALYSIS

222

deviation where

SANDNER

AND

DREYFUS

of Sj (j = 1, 2, . . . , 50). IS1 were considered to be stable if Sj = 0,

sj = (l/50)

30 X 2 .i=l

Sj.

L21

Bursting or alternating IS1 give oscillatory serial correlograms. This oscillatory behavior can be appraised by the computation of an oscillation coefficient K, from a Fourier transform of the serial correlogram (12)

where G,,(f) Example

is the power spectrum of Sj,f is the frequency.

of Result

By way of illustration, we may consider a spike train recorded from the mesencephalon. Its rate was not affected by a medial hypothalamic electrical stimulation (Fig. 5). Nevertheless, one group of IS1 was shortened, as shown on a Pdf graph (Fig. 4). Moreover, the sequential correlogram showed an increased tendency towards irregular bursting activity; this was attested to by an increased variability and a lowered oscillation coefficient. DISCUSSION

The above result strengthens our opinion that a more complete analysis of a spike train can provide useful information. The formatted storage allows occasional data exchanges between laboratories as well as easy access to our data at any time. In the future, easy access to the stored data will be of importance for two purposes: -First, the stored data will serve to answer questions that will emerge out of future discussions. This has already happened since a colleague asked information about the time constants of IS1 changes at the beginning of any manipulation. To answer that question, a program was written but not yet tested; it computes and prints these time constants from a Log transform of IS1 changes. -Second, the stored data will be used in transinformation evaluations and for the testing of models related to functional relationships within and between neuronal systems (13-15). The overriding goal of this work is to compare the theoretical models described in the literature with a set of collected experimental data as numerous as possible. The advantages resulting from the use of electronic computers are evident. But any technical innovation has to be discussed also against a background of financial considerations taking into account equivalent solutions to the same problem. A number of neurophysiologists use “one-line” data acquisition systems made of a minicomputer and a disk mass storage (16). Such a device

NEURONAL

INTERSPIKE

INTERVAL

ANALYSIS

223

has indeed several advantages: it reduces the time needed for the data acquistion as it does not comprise both the data register and replay session necessary in our method. Another advantage is the absence of any loss of reliability due to the data transmission. But there are also marked disadvantages, the first one being the cost of the minicomputer, which is at least 10 times that of our terminal. Obviously, we have to pay each use of the UNIVAC computer, but the use of this computer in the Centre de Calcul de Strasbourg-Cronenbourg is cheap for the C.N.R.S. laboratories, and we have free charge access to a large library of mathematical and statistical subroutines. Another disadvantage of the use of a minicomputer lies in the low speed and low computation abilities of such a computer. The Fourier transform of the serial correlation coefficients takes, for example, less than a 3-min computation when processing 4 x lo4 data in the UNIVAC computer. ACKNOWLEDGMENTS We are grateful to Professor P. Karli suggestions, to Mr. J. L. Michot and Mr. software assistance, and to Dr. M. Mirat technical and commercial collaboration. This work was supported by grants from (No. 78/175).

for a critical reading of the manuscript and helpful A. Ravet for valuable hardware and microcomputer and Dr. H. Sabourin, from the C.C.S.C., for their the D.G.R.S.T. (No. 79 70279) and from the D.R.E.T.

REFERENCES 1. ZADEH, L. A. From circuit theory to system theory. Proc. Z.R.E., May, 856 (1962). 2. DELGADO, J. M. R., ROBERTS, W. W., AND MILLER, N. E. Learning motivated by electrical stimulation of the brain. Amer. J. Physiol. 179, 587 (1954). 3. SANDNER, G., SCHMITT, P., AND KARLI, P. Central gray and medial hypothalamus: Correlation between escape behavior and unit activity. Brain Res. 170,459 (1979). 4. BULLOCK, T. H. Representation of information in neurons and sites for molecular participation. Proc. Nat. Acad. Sci. USA. 60, 1058 (1968). 5. ECKHORN, R., GROSSER, 0. I., PELLNITZ, K., AND PGPEL, B. Efficiency of diierent neuronal codes: information transfer calculations for three different neuronal systems. Biol. Cyberner. 22, 49 (1976). 6. RAPOPORT, A., AND NORVARTH, W. J. The theoretical channel capacity of a single neuron as determined by various coding systems. Inform. Contr. 3, 335 (1960). 7. SCHMITT, P., SANDNER, G., AND KARLI, P. Caractdristiques fonctionnelles des systbmes de renforcement: etude comportementale. Physiol. Behav. 16, 419 (1975). 8. AHO, A. V., HOPCROFT, J. E., AND ULLMAN, J. D. “The Design and Analysis of Computer Algorithms”. Addison-Wesley, Menlo Park, Calif., 1976. 9. SANDERSON, A. C., AND KOBLER, B. Sequential interval histogram analysis of non-stationary neuronal spike trains. Biol. Cybernet. 22, 61 (1976). 10. VIBERT, J. F., AND CAILLE, D. Estimation of a number of normal subpopulations in an heterogeneous sample and of their parameters; a recurrent method applied to neurobiology. Pfliigers Arch. 373, 283 (1978). Il. PERKEL, D. H., GER~TEIN, G. L., AND MOORE, G. P. Neuronal spike trains and stochastic point processes. I. The single spike train. Biophys. J. 7, 391 (1967). 12. ETEVENON, P. “Etude mdthodologique de l’Clectroenc6phaIographie quantitative, application a quelques exemples.” Docteur es Sciences thesis, Imprimerie Copidith, Paris, 1978.

224

SANDNER

AND DREYFUS

S. I. A mathematical approach to neural systems. In “Systems Neuroscience” (J. Metzler, Ed.), pp. 67-117. Academic Press, New York, 1977. 14. PERKEL, D. H. A computer program for simulating a network of interacting neurons. I. Organization and physiological assumptions. Cornput. Biomed. Res. 9, 3 1 (1976). IS. POZINE, N. “Simulation des structures neuroniques.” MIR, Moscow, 1974. 16. HARDING, G. W., AND TOWE, A. L. An automated on-line, real-time laboratory for single neuron studies. Comput. Biomcd. Res. 9, 471 (1976). 13. AMARI,