likelihood: Use Traditional KDE or Silverman's KDE-FFT Method to...
In TasCL/cpda: Probability Density Approximation

Description Usage Arguments Details Value References See Also Examples

Our KDE-FFT method uses the following steps to calculate likelihood:

draws Monte Carlo simulations to construct a simulated histogram,
uses FFT to transform the histogram to spectral domain,
applies standard Gaussian kernels to smooth it,
uses inverse FFT to get simulated PDF, and
interpolates data linearly on the sPDF to obtain likelihood.

1	likelihood(xtilde, x, h = 0, m = 0.8, p = 10, N_s = 0, pw = FALSE)

`xtilde`	vector of simulations.
`x`	vector of data.
`h`	the bandwidth for KDE. If not given, the function uses Silverman's ‘rule of thumb’, Silverman (1986, page 48, eqn (3.31)).
`m`	a multiplier to adjust bandwidth proportationally. Default is 0.8.
`p`	a precision parameter defining the number of grid as power of 2. The default grid size is 1024. That p is set at 10 (i.e., 2^10).
`N_s`	number of simulations. This is to adapt the simulation when one uses PDA to calculate likelihood in a defective probability density function, such as the LBA model or the DDM. The function will count the number of simulations from `xtilde`. If `xtilde` stores simulations from a a defective distribution, one has to supply N_s.
`pw`	whether to calculate gaussian kernel directly or using FFT algorithm. pw = FALSE uses KDE-FFT method.

Although stats::density also uses KDE-FFT method, our method includes simulation and interpolation steps.

A vector.

Holmes, W. (2015). A practical guide to the Probability Density Approximation (PDA) with improved implementation and error characterization. Journal of Mathematical Psychology, 68-69, 13–24, https://dx.doi.org/10.1016/j.jmp.2015.08.006.

Silverman, B. W. (1982). Algorithm as 176: Kernel density estimation using the fast Fourier transform. Journal of the Royal Statistical Society. Series C (Applied Statistics), 31(1), 93-99. https://dx.doi.org/10.2307/2347084.

density, bandwidth.nrd, bw.nrd, bandwidth.nrd.

###################
## Example 1     ##
###################
x   <- seq(-3, 3, length.out=100) ## Data
sam <- rnorm(1e5)              ## Monte Carlo simulation
xlabel <- "Observations"
ylabel <- "Density"

sam  <- rnorm(n)
den1 <- dnorm(x)
den2 <- likelihood(x, sam, pw=TRUE)
den3 <- likelihood(x, sam, pw=FALSE)

par(mar=c(4, 5.3, 0.82, 1))
plot(x, den1, type="l", lty="dotted", xlab=xlabel, ylab=ylabel, cex.lab=3,
     cex.axis=1.5, lwd=3)

lines(x, den2, col="blue", lty="dashed", lwd=2)
lines(x, den3, col="red", lwd=2)

## Define bandwidth using bw.nrd0
bandwidth <- 0.8*bw.nrd0(sam)          
den1 <- likelihood(x, sam, h=bandwidth, pw=TRUE)
den2 <- likelihood(x, sam, h=bandwidth, pw=FALSE)

plot(x, den1, type="l", lty="dotted", xlab=xlabel, ylab=ylabel, cex.lab=3,
     cex.axis=1.5, lwd=3)
lines(x, den2, col="blue", lty="dashed", lwd=2)
lines(x, den3, col="red", lwd=2)


###################
## Example 2     ##
###################
rm(list=ls())
## Assuming that this is empirical data
y <- rtdists::rLBA(1e3, A=.5, b=1, t0=.25, mean_v=c(2.4, 1.2), sd_v=c(1,1))
rt1  <- y[y$response==1,"rt"]
rt2  <- y[y$response==2,"rt"]
srt1 <- sort(rt1)
srt2 <- sort(rt2)

xlabel <- "RT"; ylabel <- "Density"
par(mfrow=c(1,2))
hist(rt1, "fd", freq=F, xlim=c(0.1,1.5), main="Choice 1", 
    xlab=xlabel, ylab=ylabel)
hist(rt2, "fd", freq=F, xlim=c(0.1,1.5), main="Choice 2",
    xlab=xlabel, ylab=ylabel)
par(mfrow=c(1,1))

## Now draw simulations from another rLBA
n <- 1e5
samp <- rtdists::rLBA(n, A=.5, b=1, t0=.25, mean_v=c(2.4, 1.2), sd_v=c(1,1))
samp1 <- samp[samp[,2]==1, 1]
samp2 <- samp[samp[,2]==2, 1]
fft1  <- likelihood(srt1, samp1, n=n)  
fft2  <- likelihood(srt2, samp2, n=n)

## Calculate theoretical densities
den0 <- rtdists::dLBA(y$rt, y$response, A=.5, b=1, t0=.25, 
   mean_v=c(2.4, 1.2), sd_v=c(1, 1))
   
df0 <- cbind(y, den0)
df1 <- df0[df0[,2]==1,]
df2 <- df0[df0[,2]==2,]
den1 <- df1[order(df1[,1]),3]
den2 <- df2[order(df2[,1]),3]
plot(srt1,  den1, type="l")
lines(srt2, den1)
lines(srt1, fft1, col="red")
lines(srt2, fft2, col="red")