Transform spike times from a `glm`

or `gam`

fitted model as
defined by Ogata (1988) and Brown et al (2002). If the model structure
is "correct" and if the model parameters are properly estimated the
result of the time transformation should be the realization of a
Poisson process with rate 1.

1 | ```
transformedTrain(obj, target = obj$data$event, select)
``` |

`obj` |
An object returned by |

`target` |
A binary (0,1) vector of integers with the same length
as |

`select` |
A character string defining a condition to be fulfilled
by the event in order to be selected, like: |

The `fitted.values`

component of `obj`

contains the
(estimated) probability to observe a spike in each time bin where the covariates
required by the fitted model were defined. It is then straightforward
to show using the concept of *product integral* (Kalbfleisch and
Prentice, 2002; Andersen et al, 1993),provided that the time bin width
is small enough to have a very small probability in each bin, that the
cumulated sum of these probabilities is the expected number of events
observed up to a given time. This expected number of events
which is returned by `transformedTrain`

. It is also the result of
the "time transformation" proposed by Ogata (1988) and brought to the
spike train analysis field under the name "time rescaling (theorem)"
by Brown et al (2002).

`transformedTrain`

can also be used to transform the times of the
spikes of neurons whose spike trains were simultaneously recorded and
discretized *in exactly the same way* as the neuron used to
generate `obj`

. This is useful to explore the possibility of
functional interactions between a putative pre-synaptic neuron (whose
spike train would correspond to argument `target`

) and a
post-synaptic one used to generate `obj`

.

`transformedTrain`

returns an object of class
`transformedTrain`

inheriting from class `spikeTrain`

. The
object is fundamentally a numeric vector with strictly increasing
elements containing the transformed times (or the expected number of events).

As mentioned only the spikes for which the covariates of the model are
available have their times transformed. That practically means that
the length of the `transformedTrain`

object returned by function
`transformedTrain`

*can be shorter* than the length of the
original `spikeTrain`

object (or more precisely than the number
of spikes defined in `target`

). If one works with a model
involving the elapsed times since the last three spikes then the
fourth spike of the train will be the first to be transformed. You
should therefore expect some left truncation of the data at the
beginning of each acquisition epoch.

Christophe Pouzat christophe.pouzat@gmail.com

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

Brown, E. N., Barbieri, R., Ventura, V., Kass, R. E. and Frank,
L. M. (2002) The time-rescaling theorem and its application to neural
spike train data analysis. *Neural Computation* **14**:
325-346.

Kalbfleisch, John D. and Prentice, Ross L. (2002) *The
Statistical Analysis of Failure Time Data.* Wiley Interscience.

Andersen, Per Kragh, Borgan, Ornulf, Gill, Richard D. and Keiding,
Niels (1993) *Statistical Models Based on Counting Processes.*
Springer-Verlag.

`plot.transformedTrain`

,
`summary.transformedTrain`

,
`mkGLMdf`

,
`data.frame`

,
`glm`

,
`mgcv`

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 | ```
## Not run:
## Let us consider neuron 1 of the CAL2S data set
data(CAL2S)
CAL2S <- lapply(CAL2S,as.spikeTrain)
CAL2S[["neuron 1"]]
renewalTestPlot(CAL2S[["neuron 1"]])
summary(CAL2S[["neuron 1"]])
## Make a data frame with a 4 ms time resolution
cal2Sdf <- mkGLMdf(CAL2S,0.004,0,60)
## keep the part relative to neuron 1, 2 and 3 separately
n1.cal2sDF <- cal2Sdf[cal2Sdf$neuron=="1",]
n2.cal2sDF <- cal2Sdf[cal2Sdf$neuron=="2",]
n3.cal2sDF <- cal2Sdf[cal2Sdf$neuron=="3",]
## remove unnecessary data
rm(cal2Sdf)
## Extract the elapsed time since the second to last and
## third to last for neuron 1. Normalise the result.
n1.cal2sDF[c("rlN.1","rsN.1","rtN.1")] <- brt4df(n1.cal2sDF,"lN.1",2,c("rlN.1","rsN.1","rtN.1"))
## load mgcv library
library(mgcv)
## fit a model with a tensorial product involving the last
## three spikes and using a cubic spline basis for the last two
## To gain time use a fixed df regression spline
n1S.fitA <- gam(event ~ te(rlN.1,rsN.1,bs="cr",fx=TRUE) + rtN.1,data=n1.cal2sDF,family=binomial(link="logit"))
## transform time
N1.Lambda <- transformedTrain(n1S.fitA)
## check out the resulting spike train using the fact
## that transformedTrain objects inherit from spikeTrain
## objects
N1.Lambda
## Use more formal checks
summary(N1.Lambda)
plot(N1.Lambda,which=c(1,2,4,5),ask=FALSE)
## Transform spike trains of neuron 2 and 3
N2.Lambda <- transformedTrain(n1S.fitA,n2.cal2sDF$event)
N3.Lambda <- transformedTrain(n1S.fitA,n3.cal2sDF$event)
## Check interactions
summary(N2.Lambda %frt% N1.Lambda)
summary(N3.Lambda %frt% N1.Lambda)
plot(N2.Lambda %frt% N1.Lambda,ask=FALSE)
plot(N3.Lambda %frt% N1.Lambda,ask=FALSE)
## End(Not run)
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.