psth computes and plot.psth plots a peri-stimulus time histogram (called PST,
post-stimulus time histogram by Gerstein and Kiang (1960)) from
repeated presentations of a stimulation. Confidence bands can be
obtained using the Poisson approximation.
Usage
psth(repeatedTrain, breaks = 20, include.lowest = TRUE,
right = TRUE, plot = TRUE, CI = 0.95, ...)
## S3 method for class 'psth'
plot(x, stimTimeCourse = NULL, colStim = "grey80",
colCI = NULL, xlab, ylab, main, xlim, ylim, lwd = 2,
col = 1, ...)
Arguments
repeatedTrain
a repeatedTrain object or a list which can be
coerced to such an object.
x
a psth object.
stimTimeCourse
NULL (default) or a two elements vector
specifying the time boundaries (in s) of a stimulus presentation.
colStim
the background color used for the stimulus.
breaks
a numeric. A single number is interpreted has the number
of bins; a vector of length 2 is interpreted as the bin width and
the step to use (see details); otherwise interpreted as the
position of the "breaks" between bins.
include.lowest
corresponding argument of hist.
right
corresponding argument of hist.
plot
corresponding argument of hist.
CI
The coverage probability of the confidence intervals.
colCI
if not NULL (default) a confidence band is
plotted with the specified color; two dashed lines are plotted otherwise.
xlim
a numeric (default value supplied). See
plot.
ylim
a numeric (default value supplied). See plot.
xlab
a character (default value supplied). See plot.
ylab
a character (default value supplied). See plot.
main
a character (default value supplied). See plot.
lwd
line width used to plot the estimated density. See plot.
col
color used to plot the estimated density. See plot.
...
see plot.
Details
When confidence bands are requested they are obtained from the
qunatiles of the Poisson distribution.
When a 2 elements vector is used as breaks argument it is
interpreted as specifying a bin width (first element if elements are
unnamed, "bw" element otherwise) and a step (second element if
elements are unnamed, "step" element otherwise). The idea is
then to obtain a smoother looking PSTH by counting spikes within
overlapping bins. That is if the center of the ith bin is xi the one
of the (i+1)th bin will be xi + step.
Value
When plot is set to FALSE in psth, a list of
class psth is returned and no plot
is generated. This list has the following components:
freq
a vector containing the instantaneous firing rate.
ciUp
a vector with the upper limit of the confidence band.
ciLow
a vector with the lower limit of the confidence band.
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.
counts
a matrix with as many rows as components in
repeatedTrain and as many columns as bins. Each element of
the matrix contains the number of spikes falling in a given trial in
a given bin.
nbTrials
the number of stimulations.
call
the matched call.
When plot is set to TRUE nothing is returned and a plot
is generated as a side effect. Of course the same occurs upon calling
plot.psth with a psth object argument.
## Load Vanillin responses data (first cockroach data set)
data(CAL1V)
## convert them into repeatedTrain objects
## The stimulus command is on between 4.49 s and 4.99s
CAL1V <- lapply(CAL1V,as.repeatedTrain)
## look at the individual raster plots
plot(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),main="N1")
## Create a simple black and white PSTH for neuron 1
psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=20)
## Rebuilt the same PSTH but with red confidence bands
psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=20,colCI=2)
## Make the PSTH smoother
psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=c(bw=0.5,step=0.05),colCI=2)
## Make a plot with PSTHs from 4 neurons superposed
## First get lists containing PSTHs from each neuron
psth1 <- psth(CAL1V[["neuron 1"]],breaks=c(bw=0.5,step=0.05),plot=FALSE)
psth2 <- psth(CAL1V[["neuron 2"]],breaks=c(bw=1,step=0.1),plot=FALSE)
psth3 <- psth(CAL1V[["neuron 3"]],breaks=c(bw=0.5,step=0.05),plot=FALSE)
psth4 <- psth(CAL1V[["neuron 4"]],breaks=c(bw=2,step=0.2),plot=FALSE)
## Get the maximal frequency to display
maxFreq <- max(max(psth1$ciUp),max(psth2$ciUp),max(psth3$ciUp),max(psth4$ciUp))
## Build plot
plot(c(0,10),c(0,75),type="n",
xaxs="i",yaxs="i",xlab="Time (s)",
ylab="Freq. (Hz)",
main="PSTHs from 4 simultaneously recorded neurons",
sub="20 stimulations with vanillin were used.")
## Add rectangle corresponding to stimulation command
rect(4.49,0,4.99,75,col="grey80",lty=0)
## Add the neurons PSTHs as confidence bands
polygon(c(psth1$mids,rev(psth1$mids)),c(psth1$ciLow,rev(psth1$ciUp)),col=1,border=NA)
polygon(c(psth2$mids,rev(psth2$mids)),c(psth2$ciLow,rev(psth2$ciUp)),col=2,border=NA)
polygon(c(psth3$mids,rev(psth3$mids)),c(psth3$ciLow,rev(psth3$ciUp)),col=3,border=NA)
polygon(c(psth4$mids,rev(psth4$mids)),c(psth4$ciLow,rev(psth4$ciUp)),col=4,border=NA)
legend(0.1,maxFreq,legend=paste("neuron",1:4),lty=1,col=1:4,bty="n")
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(STAR)
Loading required package: survival
Loading required package: mgcv
Loading required package: nlme
This is mgcv 1.8-12. For overview type 'help("mgcv-package")'.
Loading required package: R2HTML
Loading required package: gss
Loading required package: codetools
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/STAR/psth.Rd_%03d_medium.png", width=480, height=480)
> ### Name: psth
> ### Title: Compute and Plot Peri-Stimulus Time Histogram
> ### Aliases: psth plot.psth
> ### Keywords: ts survival
>
> ### ** Examples
>
> ## Load Vanillin responses data (first cockroach data set)
> data(CAL1V)
> ## convert them into repeatedTrain objects
> ## The stimulus command is on between 4.49 s and 4.99s
> CAL1V <- lapply(CAL1V,as.repeatedTrain)
> ## look at the individual raster plots
> plot(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),main="N1")
> ## Create a simple black and white PSTH for neuron 1
> psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=20)
> ## Rebuilt the same PSTH but with red confidence bands
> psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=20,colCI=2)
> ## Make the PSTH smoother
> psth(CAL1V[["neuron 1"]],stimTimeCourse=c(4.49,4.99),breaks=c(bw=0.5,step=0.05),colCI=2)
> ## Make a plot with PSTHs from 4 neurons superposed
> ## First get lists containing PSTHs from each neuron
> psth1 <- psth(CAL1V[["neuron 1"]],breaks=c(bw=0.5,step=0.05),plot=FALSE)
> psth2 <- psth(CAL1V[["neuron 2"]],breaks=c(bw=1,step=0.1),plot=FALSE)
> psth3 <- psth(CAL1V[["neuron 3"]],breaks=c(bw=0.5,step=0.05),plot=FALSE)
> psth4 <- psth(CAL1V[["neuron 4"]],breaks=c(bw=2,step=0.2),plot=FALSE)
> ## Get the maximal frequency to display
> maxFreq <- max(max(psth1$ciUp),max(psth2$ciUp),max(psth3$ciUp),max(psth4$ciUp))
> ## Build plot
> plot(c(0,10),c(0,75),type="n",
+ xaxs="i",yaxs="i",xlab="Time (s)",
+ ylab="Freq. (Hz)",
+ main="PSTHs from 4 simultaneously recorded neurons",
+ sub="20 stimulations with vanillin were used.")
> ## Add rectangle corresponding to stimulation command
> rect(4.49,0,4.99,75,col="grey80",lty=0)
> ## Add the neurons PSTHs as confidence bands
> polygon(c(psth1$mids,rev(psth1$mids)),c(psth1$ciLow,rev(psth1$ciUp)),col=1,border=NA)
> polygon(c(psth2$mids,rev(psth2$mids)),c(psth2$ciLow,rev(psth2$ciUp)),col=2,border=NA)
> polygon(c(psth3$mids,rev(psth3$mids)),c(psth3$ciLow,rev(psth3$ciUp)),col=3,border=NA)
> polygon(c(psth4$mids,rev(psth4$mids)),c(psth4$ciLow,rev(psth4$ciUp)),col=4,border=NA)
> legend(0.1,maxFreq,legend=paste("neuron",1:4),lty=1,col=1:4,bty="n")
>
>
>
>
>
> dev.off()
null device
1
>