predict.dapp: Predict Method for Dynamic Admixture of Poisson Process

In neuromplex: Neural Multiplexing Analysis

Description

Summarizes predictive draws of weight curves from a fitted DAPP model

Usage

 
## S3 method for class 'dapp'
predict(object, tilt.prior = FALSE,
    mesh.tilt = 0.1, nprior = object$mcmc["nsamp"], ...)

Arguments

object	a fitted model of the class 'dapp'
tilt.prior	logical giving whether the prior should be tilted to mimic an analysis done with a uniform prior on the range(alpha)
mesh.tilt	a tuning parameter that controls how exactly tilting is done. Shorter mesh value gives tighter match but will require more Monte Carlo simulations
nprior	number of prior draws to be used for display
...	no addiitonal parameters used at this point

Details

This function is intended to be mostly used through predict.dapp.

Value

Gives prior and posterior summaries of the range and average predicted alpha curves. Also gives the same for the posterior draws of alpha for each recorded AB trial.

See Also

dapp, plot.dapp and summary.dapp.

Examples

## Not run: 
## generate 25 A and 30 B trials with rate functions
##    lambda.A(t) = 160*exp(-2*t/1000) + 40*exp(-0.2*t/1000)
##    lambda.B(t) = 40*exp(-2*t/1000)
## where time t is measured in ms. Then, generate 40 AB trials,
## roughly half with flat weight curves with a constant intensity
## either close to A, or close to B or close to the 50-50 mark,
## (equally likely). The remaining curves are sinusoidal
## that snake between 0.01 and 0.99 with a period randomly
## drawn between 400 and 1000

ntrials <- c(nA=25, nB=30, nAB=40)
flat.range <- list(A=c(0.85, 0.95),
                   B=c(0.05, 0.15),
                   mid=c(0.45,0.55))
flat.mix <- c(A=1/3, B=1/3, mid=1/3)
wavy.span <- c(0.01, 0.99)
wavy.period <- c(400, 1000)

T.horiz <- 1000
rateB <- 40 * exp(-2*(1:T.horiz)/T.horiz)
rateA <- 4*rateB + 40 * exp(-0.2*(1:T.horiz)/T.horiz)

synth.data <- synthesis.dapp(ntrials = ntrials, pr.flat = 0.5,
                             intervals = flat.range, wts = flat.mix,
                             span = wavy.span, period.range = wavy.period,
                             lambda.A=rateA, lambda.B=rateB)

## Visualize data and generated binned spike counts
spike.counts <- mplex.preprocess(synth.data$spiketimes, visualize=TRUE)

## Fit the DAPP model to data
fit.post <- dapp(spike.counts, verbose=FALSE)

## Prediction
pp <- predict(fit.post)

## Visualizing (range, ave) of alpha(t) for each recorded AB trial
te <- pp$trial.est
ggplot(te, aes(x=ave, y=range)) +
    stat_density_2d(aes(fill = ..level..), h=0.2, geom = "polygon") +
    scale_fill_viridis_c() +
    theme_bw() +
    facet_wrap(~as.factor(trial))

## Post process results to assign second order stochasticity labels
summary(fit.post)

## End(Not run)

neuromplex documentation built on April 22, 2021, 5:11 p.m.