Advertisement
Guest User

Untitled

a guest
Jul 19th, 2019
96
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 0.94 KB | None | 0 0
  1. library(pomp)
  2. # state function
  3. gompertz.proc.sim <- function(x, t, params, delta.t, ...) {
  4. eps <- exp(rnorm(n = 1, mean = 0, sd = params["sigma"]))
  5. S <- exp(-params["r"] * delta.t)
  6. setNames(params["K"]^(1 - S) * x["X"]^S * eps, "X")
  7. }
  8. # obs
  9. gompertz.meas.sim <- function(x, t, params, ...){
  10. setNames(rlnorm(n = 1, meanlog = log(x["X"]), sdlog = params["tau"]), "Y")
  11. }
  12. # density
  13. gompertz.meas.dens <- function(y, x, t, params, log, ...) {
  14. dlnorm(x=y["Y"], meanlog = log(x["X"]), sdlog = params["tau"], log = log )
  15. }
  16. # build our 'pomp' object
  17. gompertz <- pomp(data = data.frame(time = 1:100, Y = NA), times = "time",
  18. rprocess = discrete_time(step.fun = gompertz.proc.sim, delta.t = 1),
  19. rmeasure = gompertz.meas.sim, t0 = 0
  20. )
  21. # assign values to params
  22. theta <- c(r = 0.1, K = 1, sigma = 0.1, tau = 0.1, X.0 = 1)
  23. # simulate values for our obs
  24. gompertz <- simulate(gompertz, params = theta )
  25. ```
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement