Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- library(pomp)
- # state function
- gompertz.proc.sim <- function(x, t, params, delta.t, ...) {
- eps <- exp(rnorm(n = 1, mean = 0, sd = params["sigma"]))
- S <- exp(-params["r"] * delta.t)
- setNames(params["K"]^(1 - S) * x["X"]^S * eps, "X")
- }
- # obs
- gompertz.meas.sim <- function(x, t, params, ...){
- setNames(rlnorm(n = 1, meanlog = log(x["X"]), sdlog = params["tau"]), "Y")
- }
- # density
- gompertz.meas.dens <- function(y, x, t, params, log, ...) {
- dlnorm(x=y["Y"], meanlog = log(x["X"]), sdlog = params["tau"], log = log )
- }
- # build our 'pomp' object
- gompertz <- pomp(data = data.frame(time = 1:100, Y = NA), times = "time",
- rprocess = discrete_time(step.fun = gompertz.proc.sim, delta.t = 1),
- rmeasure = gompertz.meas.sim, t0 = 0
- )
- # assign values to params
- theta <- c(r = 0.1, K = 1, sigma = 0.1, tau = 0.1, X.0 = 1)
- # simulate values for our obs
- gompertz <- simulate(gompertz, params = theta )
- ```
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement