Description Usage Arguments References Examples
Bayesian estimation of a model dY_t = b(φ_j,t,Y_t)dt + γ \widetilde{s}(t,Y_t)dW_t, φ_j\sim N(μ, Ω), Y_{t_0}=y_0(φ, t_0).
1 2 3 |
model.class |
class of the hierarchical diffusion model including all required information, see |
t |
list or vector of time points |
data |
list or matrix of observation variables |
nMCMC |
length of Markov chain |
propSd |
vector of proposal variances for φ |
adapt |
if TRUE (default), proposal variance is adapted |
proposal |
proposal density: "normal" (default) or "lognormal" (for positive parameters) |
Hermann, S., K. Ickstadt and C. H. Mueller (2016). Bayesian Prediction of Crack Growth Based on a Hierarchical Diffusion Model. Applied Stochastic Models in Business and Industry, DOI: 10.1002/asmb.2175.
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 | mu <- 2; Omega <- 0.4; phi <- matrix(rnorm(21, mu, sqrt(Omega)))
model <- set.to.class("mixedDiffusion",
parameter = list(phi = phi, mu = mu, Omega = Omega, gamma2 = 0.1),
b.fun = function(phi, t, x) phi*x, sT.fun = function(t, x) x)
t <- seq(0, 1, by = 0.01)
data <- simulate(model, t = t, plot.series = TRUE)
est <- estimate(model, t, data[1:20,], 100) # nMCMC should be much larger
plot(est)
# OU
b.fun <- function(phi, t, y) phi[1]-phi[2]*y; y0.fun <- function(phi, t) phi[3]
mu <- c(10, 5, 0.5); Omega <- c(0.9, 0.01, 0.01)
phi <- sapply(1:3, function(i) rnorm(21, mu[i], sqrt(Omega[i])))
model <- set.to.class("mixedDiffusion",
parameter = list(phi = phi, mu = mu, Omega = Omega, gamma2 = 0.1),
y0.fun = y0.fun, b.fun = b.fun, sT.fun = function(t, x) 1)
t <- seq(0, 1, by = 0.01)
data <- simulate(model, t = t, plot.series = TRUE)
est <- estimate(model, t, data[1:20,], 100) # nMCMC should be much larger
plot(est)
##
t.list <- list()
for(i in 1:20) t.list[[i]] <- t
t.list[[21]] <- t[1:50]
data.list <- list()
for(i in 1:20) data.list[[i]] <- data[i,]
data.list[[21]] <- data[21, 1:50]
est <- estimate(model, t.list, data.list, 100)
pred <- predict(est, t = t[50:101], which.series = "current", ind.pred = 21,
b.fun.mat = function(phi, t, y) phi[,1]-phi[,2]*y)
|
10 of 51 predictions are calculated
20 of 51 predictions are calculated
30 of 51 predictions are calculated
40 of 51 predictions are calculated
50 of 51 predictions are calculated
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.