diff --git a/DESCRIPTION b/DESCRIPTION index 6437e060..d6f921fe 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,8 +1,8 @@ Package: pomp Type: Package Title: Statistical Inference for Partially Observed Markov Processes -Version: 5.9.0.0 -Date: 2024-06-02 +Version: 5.9.1.0 +Date: 2024-06-07 Authors@R: c(person(given=c("Aaron","A."),family="King",role=c("aut","cre"),email="kingaa@umich.edu",comment=c(ORCID="0000-0001-6159-3207")), person(given=c("Edward","L."),family="Ionides",role="aut",comment=c(ORCID="0000-0002-4190-0174")) , person(given="Carles",family="Bretó",role="aut",comment=c(ORCID="0000-0003-4695-4902")), diff --git a/R/rprocess_spec.R b/R/rprocess_spec.R index c3ff9568..8203d740 100644 --- a/R/rprocess_spec.R +++ b/R/rprocess_spec.R @@ -52,7 +52,7 @@ ##' They differ as to how this is done. ##' Specifically, ##' \enumerate{ -##' \item \code{onestep} takes a single step to go from any given time \code{t1} to any later time \code{t2} (\code{t1 < t2}). +##' \item \code{onestep} takes a single step to go from any given time \code{t1} to any later time \code{t2} (\code{t1 <= t2}). ##' Thus, this plug-in is designed for use in situations where a closed-form solution to the process exists. ##' \item To go from \code{t1} to \code{t2}, \code{euler} takes \code{n} steps of equal size, where ##' \preformatted{ diff --git a/inst/NEWS b/inst/NEWS index 1c6b092a..dd909969 100644 --- a/inst/NEWS +++ b/inst/NEWS @@ -1,5 +1,16 @@ _N_e_w_s _f_o_r _p_a_c_k_a_g_e '_p_o_m_p' +_C_h_a_n_g_e_s _i_n '_p_o_m_p' _v_e_r_s_i_o_n _5._9._1: + + • The rprocess plugin ‘onestep’ now always takes exactly one + step to get from one observation time to the next, even when + the interval between them is zero. That is, if ‘P’ is a + ‘pomp’ object with an rprocess built using + ‘rprocess=onestep(f)’, ‘t0=timezero(P)’, and ‘t=time(P)’, + then ‘f’ will be called exactly once for the interval + (‘t0’,‘t[1]’) and once for each of the intervals + (‘t[i-1]’,‘t[i]’), ‘i=2,...,length(t)’. + _C_h_a_n_g_e_s _i_n '_p_o_m_p' _v_e_r_s_i_o_n _5._8._4: • To add additional elements for use by the basic model diff --git a/inst/NEWS.Rd b/inst/NEWS.Rd index f74466e3..3943ddb0 100644 --- a/inst/NEWS.Rd +++ b/inst/NEWS.Rd @@ -1,5 +1,11 @@ \name{NEWS} \title{News for package `pomp'} +\section{Changes in \pkg{pomp} version 5.9.1}{ + \itemize{ + \item The rprocess plugin \code{onestep} now always takes exactly one step to get from one observation time to the next, even when the interval between them is zero. + That is, if \code{P} is a \sQuote{pomp} object with an rprocess built using \code{rprocess=onestep(f)}, \code{t0=timezero(P)}, and \code{t=time(P)}, then \code{f} will be called exactly once for the interval (\code{t0},\code{t[1]}) and once for each of the intervals (\code{t[i-1]},\code{t[i]}), \code{i=2,...,length(t)}. + } +} \section{Changes in \pkg{pomp} version 5.8.4}{ \itemize{ \item To add additional elements for use by the basic model components (i.e., \dQuote{userdata}), the \code{userdata} argument should now be used. diff --git a/man/rprocess_spec.Rd b/man/rprocess_spec.Rd index 9e4414f8..1bd8a547 100644 --- a/man/rprocess_spec.Rd +++ b/man/rprocess_spec.Rd @@ -100,7 +100,7 @@ The simulator plug-ins \code{discrete_time}, \code{euler}, and \code{onestep} al They differ as to how this is done. Specifically, \enumerate{ -\item \code{onestep} takes a single step to go from any given time \code{t1} to any later time \code{t2} (\code{t1 < t2}). +\item \code{onestep} takes a single step to go from any given time \code{t1} to any later time \code{t2} (\code{t1 <= t2}). Thus, this plug-in is designed for use in situations where a closed-form solution to the process exists. \item To go from \code{t1} to \code{t2}, \code{euler} takes \code{n} steps of equal size, where \preformatted{ diff --git a/src/euler.c b/src/euler.c index 85344e53..de3c4ee5 100644 --- a/src/euler.c +++ b/src/euler.c @@ -208,7 +208,7 @@ SEXP euler_simulator switch (method) { case onestep: default: // one step dt = time[step]-t; - nstep = (dt > 0) ? 1 : 0; + nstep = 1; break; case discrete: // fixed step dt = deltat; diff --git a/src/rprocess.c b/src/rprocess.c index 5255cb29..8b233680 100644 --- a/src/rprocess.c +++ b/src/rprocess.c @@ -111,7 +111,7 @@ SEXP do_rprocess (SEXP object, SEXP xstart, SEXP tstart, SEXP times, SEXP params double deltat = 1.0; PROTECT(fn = GET_SLOT(rproc,install("step.fn"))); PROTECT(X = euler_simulator(fn,xstart,tstart,times,params,deltat,type, - accumvars,covar,args,gnsi)); + accumvars,covar,args,gnsi)); nprotect += 2; } break; @@ -122,7 +122,7 @@ SEXP do_rprocess (SEXP object, SEXP xstart, SEXP tstart, SEXP times, SEXP params PROTECT(fn = GET_SLOT(rproc,install("step.fn"))); deltat = *(REAL(AS_NUMERIC(GET_SLOT(rproc,install("delta.t"))))); PROTECT(X = euler_simulator(fn,xstart,tstart,times,params,deltat,type, - accumvars,covar,args,gnsi)); + accumvars,covar,args,gnsi)); nprotect += 2; } break; diff --git a/tests/steps-01.png b/tests/steps-01.png deleted file mode 100644 index 44cc8d8f..00000000 Binary files a/tests/steps-01.png and /dev/null differ diff --git a/tests/steps-02.png b/tests/steps-02.png deleted file mode 100644 index 7717ab29..00000000 Binary files a/tests/steps-02.png and /dev/null differ diff --git a/tests/steps.R b/tests/steps.R index a5d954f4..37227bad 100644 --- a/tests/steps.R +++ b/tests/steps.R @@ -1,9 +1,10 @@ options(digits=3) -library(pomp) suppressPackageStartupMessages({ + library(tidyr) library(dplyr) }) +library(pomp) tm <- freeze(sort(runif(n=20,max=50)),seed=1930502785L) simulate( @@ -98,3 +99,45 @@ stopifnot( as.numeric(states(sm,"steps")) ) ) + +simulate( + t0=0, + times=c(0,0,0,1,3,3,7), + rinit=\(t0,...)c(x=t0,n=0), + rprocess=onestep(\(t,x,n,delta.t,...) c(x=t+delta.t,n=n+1)), + rmeasure=\(x,n,...) c(y1=x,y2=n), + dmeasure=\(...,log) if (log) 0 else 1 +) -> s1 + +simulate( + t0=0, + times=c(0,0,0,1,3,3,7), + rinit=Csnippet("x=t; n=0;"), + rprocess=onestep(Csnippet("x=t+dt; n=n+1;")), + rmeasure=Csnippet("y1 = x; y2 = n;"), + dmeasure=Csnippet("lik = (give_log) ? 0 : 1;"), + statenames=c("x","n"), + obsnames=c("y1","y2") +) -> s2 + +s1 |> + as.data.frame() |> + with( + stopifnot( + time==x, + n==seq_along(time) + ) + ) + +s2 |> + pfilter(Np=1,save.states=TRUE) |> + saved_states(format="list") |> + melt() |> + pivot_wider() |> + mutate(time=time(s2)[.L1]) |> + with( + stopifnot( + time==x, + n==seq_along(time) + ) + ) diff --git a/tests/steps.Rout.save b/tests/steps.Rout.save index ddcb7adc..146928f3 100644 --- a/tests/steps.Rout.save +++ b/tests/steps.Rout.save @@ -19,10 +19,11 @@ Type 'q()' to quit R. > options(digits=3) > -> library(pomp) > suppressPackageStartupMessages({ ++ library(tidyr) + library(dplyr) + }) +> library(pomp) > > tm <- freeze(sort(runif(n=20,max=50)),seed=1930502785L) > simulate( @@ -120,3 +121,45 @@ Warning message: + ) + ) > +> simulate( ++ t0=0, ++ times=c(0,0,0,1,3,3,7), ++ rinit=\(t0,...)c(x=t0,n=0), ++ rprocess=onestep(\(t,x,n,delta.t,...) c(x=t+delta.t,n=n+1)), ++ rmeasure=\(x,n,...) c(y1=x,y2=n), ++ dmeasure=\(...,log) if (log) 0 else 1 ++ ) -> s1 +> +> simulate( ++ t0=0, ++ times=c(0,0,0,1,3,3,7), ++ rinit=Csnippet("x=t; n=0;"), ++ rprocess=onestep(Csnippet("x=t+dt; n=n+1;")), ++ rmeasure=Csnippet("y1 = x; y2 = n;"), ++ dmeasure=Csnippet("lik = (give_log) ? 0 : 1;"), ++ statenames=c("x","n"), ++ obsnames=c("y1","y2") ++ ) -> s2 +> +> s1 |> ++ as.data.frame() |> ++ with( ++ stopifnot( ++ time==x, ++ n==seq_along(time) ++ ) ++ ) +> +> s2 |> ++ pfilter(Np=1,save.states=TRUE) |> ++ saved_states(format="list") |> ++ melt() |> ++ pivot_wider() |> ++ mutate(time=time(s2)[.L1]) |> ++ with( ++ stopifnot( ++ time==x, ++ n==seq_along(time) ++ ) ++ ) +>