hist.lockedTrain: Auto- and Cross-Intensity Function Estimate for Spike Trains

Description

hist.lockedTrain constructs and plot.hist.lockedTrain plots estimates of what Cox and Lewis (1966) call the auto- or cross-intensity functions. The auto-intensity function is also called the renewal density by Cox and Lewis and by Perkel et al (1967) while it is called the intensity of the Palm distribution by Ogata (1988). The (estimate of) the cross-intensity function is called cross-correlation function by Perkel et al (1967b) and cross-correlation histogram by Brillinger et al (1976).

Usage

"hist"(x, bw, breaks = NULL, plot = TRUE, ...)
"plot"(x, style = c("Ogata", "Brillinger"), CI, unit = "s", xlab, ylab, xlim, ylim, type, pch, ...)

Arguments

a lockedTrain object returned by the lockedTrain function.

a non-negative numeric, the bin width.

breaks

a vector giving the breakpoints between cells. If NULL (default) breaks are built using argument bw and component laglim of obj.

plot

a logical. If TRUE a plot is generated as a side effect and nothing is returned, if FALSE a list of class hist.lockedTime is returned.

style

a character. The style of the plot, "Ogata" or "Brillinger".

a numeric vector with at most two elements. The coverage probability of the confidence intervals.

unit

a character. The unit in which the spike times are expressed.

xlab

a character. The x label. Default supplied.

ylab

a character. The y label. Default supplied.

xlim

a numeric. See plot. Default supplied.

ylim

a numeric. See plot. Default supplied.

type

see lines. Default supplied.

pch

see plot. Default supplied.

...

see plot.

Value

density: the density estimate. Equivalent of the component density returned by hist.
breaks: a numeric vector with the breaks in between which spikes were counted. Similar to the component of the same name returned by hist.
mids: a numeric vector with the mid points of breaks. . Similar to the component of the same name returned by hist.
bw: the bin width used.
nRef: the total number of reference spikes used.
refFreq: the mean frequency of the reference neuron.
testFreq: the mean frequency of the test neuron.
obsTime: the total observation time used (in s).
CCH: a logical which is TRUE if a cross-intensity was estimated and FALSE in the case of an auto-intensity.
call: the matched call.

Details

The intensity of the Palm distribution (Ogata, 1988, p 13) is estimated by:

m (s) = \frac{Prob (event in (t + s, t + s + Δ s) ∣ event at t)}{Δ s}

It is called renewal density by Perkel et al (1967) and defined by their Eq. 10, p 404. Under the null hypothesis of a stationary Poisson process it is a constant whose value is the mean discharge rate.

The cross-intensity function of two spike trains A and B is estimated by (Perkel et al, 1967b, p424, Eq. 4 and 5): $m_{A B} (s) = \frac{Prob (B event in (t + s, t + s + Δ s) ∣ A event at t)}{Δ s}$ The style argument of plot.hist.lockedTrain generates a plot looking like Fig. 6, p 18 of Ogata (1988) if set to "Ogata". Using style "Brillinger" plots the square root of the estimate.

References

Ogata, Yosihiko (1988) Statistical Models for Earthquake Occurrences and Residual Analysis for Point Processes. Journal of the American Statistical Association 83: 9-27.

D. R. Cox and P. A. W. Lewis (1966) The Statistical Analysis of Series of Events. John Wiley and Sons.

J. A. McFadden (1962) On the Lengths of Intervals in a Stationary Point Process. Journal of the Royal Statistical Society. Series B, 24: 364-382 Perkel D. H., Gerstein, G. L. and Moore G. P. (1967) Neural Spike Trains and Stochastic Point Processes. I. The Single Spike Train. Biophys. J., 7: 391-418. http://www.pubmedcentral.nih.gov/articlerender.fcgi?tool=pubmed&pubmedid=4292791

Perkel D. H., Gerstein, G. L. and Moore G. P. (1967b) Neural Spike Trains and Stochastic Point Processes. I. Simultaneous Spike Trains. Biophys. J., 7: 419-440. http://www.pubmedcentral.nih.gov/articlerender.fcgi?tool=pubmed&pubmedid=4292792

David R. Brillinger, Hugh L. Bryant and Jose P. Segundo (1976) Identification of synaptic interactions. Biol Cybern, 22: 213-228.

David R. Brillinger (1976) Estimation of the Second-Order Intensities of a Bivariate Stationary Point Process. Journal of the Royal Statistical Society. Series B (Methodological), 38, 60-66.

Examples

Run this code

## Reproduce Fig. 6 of Ogata 1988
data(ShallowShocks)
shalShocks <- lockedTrain(as.spikeTrain(ShallowShocks$Date),,c(0,500))
shalShocksH <- hist(shalShocks,5,plot=FALSE)
plot(shalShocksH,"Ogata",c(0.95,0.99),xlab="TIME LAG (DAYS)",ylab="NUMBER OF EVENTS PER DAY")
## Reproduce Fig. 7 of Ogata 1988
myBinNb <- 101
myMidPoints <- seq(from = 1, to = 6, length.out=myBinNb)
myMidPoints <- 10^myMidPoints/200
myBreaks <- c(0,myMidPoints[-myBinNb] + diff(myMidPoints) / 2)
shalShocksH2 <- hist(shalShocks,breaks=myBreaks,plot=FALSE)
yy <- abs(shalShocksH2$density - shalShocksH2$refFreq)
plot(shalShocksH2$mids[shalShocksH2$density>0],
     yy[shalShocksH2$density>0],
     pch = 1,
     xlim = c(0.001,10000),
     log = "xy",
     xlab = "TIME LAG (DAYS)",
     ylab = "NUMBER OF EVENTS PER DAY"
     )

Run the code above in your browser using DataLab

Last chance! 50% off unlimited learning