Nothing
#' Sampling from a Brownian bridge path give a grid time
#'
#' A Brownian bridge is a continuous-time stochastic process
#' B(t) whose probability distribution is the conditional
#' probability distribution of a standard Wiener process W(t)
#' subject to the condition (when standardized) that W(T) = 0,
#' so that the process is pinned to the same value at both t = 0 and t = T.
#' The implementation here is a generalized Brownian bridge
#' that allows start point and end point at different locations.
#'
#' @param time time points at which observations are to be simulated
#' @param start_pt the start point location of Brownian bridge
#' @param end_pt the end point location of Brownian brige
#' @param sigma volatility parameter of the Brownian motion
#'
#' @return
#' A \code{data.frame} whose first column is the time points
#' and second column is coordinate of the locations.
#'
#' @examples
#' ## Brownian bridge starting from location 0 and ending at location 1
#' ## with sigma 0.1 from time 0 to time 10
#' plot(rBB(seq(0, 10, length.out = 100), 0, 1, 0.1), type = "l")
#'
#' @export
rBB <- function(time, start_pt, end_pt, sigma) {
time <- time - time[1]
ntime <- length(time)
u <- time[1]
t <- time[ntime]
alpha <- start_pt
beta <- end_pt
sim_BB <- numeric(ntime)
for (i in seq_len(ntime)) {
s <- time[i]
cart_mean <- ((t-s)*alpha + (s-u)*beta)/(t-u)
cart_sd <- sqrt(sigma^2 * (s-u)*(t-s)/(t-u))
sim_BB[i] <- rnorm(1, cart_mean, cart_sd)
}
data.frame(time = time, BB = sim_BB)
}
sim1mrb.bbz <- function(s, sigma, time, brtimes, t0moving=TRUE,
start_pt, end_pt) {
#### time: time points in [0, s]
tt <- sort(unique(c(time, brtimes)))
if (time[length(time)] < brtimes[length(brtimes)]) {
tt <- tt[seq_len(length(tt)-1)]
}
nt <- length(tt)
nb <- length(brtimes)
status <- as.integer(t0moving)
tend <- brtimes[1]
j <- 1
mov_segments <- numeric(0)
for (i in 2:nt) {
if (tt[i] <= tend) { ## status unchanged
if (status == 1) ## moving
mov_segments <- c(mov_segments, tt[i] - tt[i - 1])
}
if (tt[i] == tend) {
status <- 1 - status ## switch status
j <- j + 1
tend <- brtimes[j]
}
}
sim_BB <- rBB(c(0, cumsum(mov_segments)), start_pt, end_pt, sigma)[, 2]
x <- rep(NA, nt)
x[1] <- start_pt
status <- as.integer(t0moving)
tend <- brtimes[1]
j <- 1
k <- 2
for (i in 2:nt) {
if (tt[i] <= tend) { ## status unchanged
if (status == 1) ## moving
{
x[i] <- sim_BB[k]
k <- k + 1
}
else ## resting
{x[i] <- x[i - 1]}
}
if (tt[i] == tend) {
status <- 1 - status ## switch status
j <- j + 1
tend <- brtimes[j]
}
}
x[tt %in% time]
}
#' Sampling from a Moving-Resting bridge
#'
#' A moving-resting process consists of two states: moving and resting.
#' The transition between the two states is modeled by an alternating
#' renewal process, with exponentially distributed duration. An animal
#' stays at the same location while resting, and moves according to a
#' Brownian motion while moving. A moving-resting bridge is an extension
#' of Brownian bridge that uses moving-resting process to bridge given
#' starting point and ending point.
#'
#' @param time time points at which observations are to be simulated
#' @param start_pt the start point location of Brownian bridge
#' (numeric for one dimension, vector for multiple dimension)
#' @param end_pt the end point location of Brownian brige
#' (numeric for one dimension, vector for multiple dimension)
#' @param lamM rate parameter of the exponential duration while moving
#' @param lamR rate parameter of the exponential duration while resting
#' @param sigma volatility parameter of the Brownian motion while moving
#' @param s0 the state at time 0, must be one of "m" or "r", for moving and
#' resting, respectively
#'
#' @return
#' A \code{data.frame} whose first column is the time points and whose
#' other columns are coordinates of the locations.
#'
#' @examples
#' tgrid <- seq(0, 10, length=1001)
#' ## make it irregularly spaced
#' tgrid <- sort(sample(tgrid, 800))
#'
#' ## a 2-dim moving-resting bridge starting from c(0, 0)
#' ## and ending at c(10, -10)
#' dat <- rMRB(tgrid, start_pt = c(0, 0), end_pt = c(10, -10),
#' lamM = 1, lamR = 1, sigma = 1, "m")
#' par(mfrow = c(2, 1))
#' plot(dat[,1], dat[,2], xlab="t", ylab="X(t)", type='l')
#' plot(dat[,1], dat[,3], xlab="t", ylab="X(t)", type='l')
#' par(mfrow = c(1, 1))
#'
#' @export
rMRB <- function(time, start_pt, end_pt,
lamM, lamR, sigma, s0) {
time <- time - time[1]
stopifnot(s0 %in% c("m", "r"))
t0moving <- (s0 == "m")
stopifnot(length(start_pt) == length(end_pt))
dim = length(start_pt)
lam1 <- if (t0moving) lamM else lamR
lam2 <- if (t0moving) lamR else lamM
tmax <- time[length(time)]
brtimes <- sim1mr.times.bbz(tmax, lam1, lam2)
coord <- matrix(NA, ncol = dim, nrow = length(time))
for (i in seq_len(dim)) {
coord[, i] <- sim1mrb.bbz(tmax, sigma, time, brtimes[, 1], t0moving,
start_pt[i], end_pt[i])
}
data.frame(time = time, coord)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.