[Pomp-commits] r463 - in pkg: . R data inst inst/data-R inst/doc inst/examples inst/include man src tests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon May 9 21:58:29 CEST 2011
Author: kingaa
Date: 2011-05-09 21:58:28 +0200 (Mon, 09 May 2011)
New Revision: 463
Added:
pkg/inst/examples/gompertz.R
pkg/inst/examples/gompertz.c
pkg/tests/pomppomp.R
pkg/tests/pomppomp.Rout.save
Modified:
pkg/NAMESPACE
pkg/R/pomp.R
pkg/data/blowflies.rda
pkg/data/dacca.rda
pkg/data/euler.sir.rda
pkg/data/gillespie.sir.rda
pkg/data/gompertz.rda
pkg/data/ou2.rda
pkg/data/ricker.rda
pkg/data/rw2.rda
pkg/data/verhulst.rda
pkg/inst/ChangeLog
pkg/inst/NEWS
pkg/inst/data-R/euler.sir.R
pkg/inst/data-R/gompertz.R
pkg/inst/data-R/ou2.R
pkg/inst/data-R/ricker.R
pkg/inst/data-R/verhulst.R
pkg/inst/doc/advanced_topics_in_pomp.Rnw
pkg/inst/doc/advanced_topics_in_pomp.pdf
pkg/inst/doc/gompertz-trajmatch.rda
pkg/inst/doc/intro_to_pomp.Rnw
pkg/inst/doc/intro_to_pomp.pdf
pkg/inst/examples/logistic.R
pkg/inst/examples/sir.R
pkg/inst/include/pomp.h
pkg/man/pomp-class.Rd
pkg/man/pomp.Rd
pkg/src/pomp.h
pkg/src/skeleton.c
pkg/tests/fhn.R
pkg/tests/fhn.Rout.save
pkg/tests/logistic.R
pkg/tests/logistic.Rout.save
pkg/tests/ou2-simulate.Rout.save
pkg/tests/ricker.Rout.save
pkg/tests/rw2.R
pkg/tests/rw2.Rout.save
pkg/tests/sir.R
pkg/tests/sir.Rout.save
pkg/tests/skeleton.Rout.save
Log:
- the 'pomp' arguments 'skeleton.map' and 'skeleton.vectorfield' are now deprecated in favor of 'skeleton' and 'skeleton.type'
- 'pomp' is now a generic function with methods for data of type 'numeric', 'matrix', 'data.frame', and 'pomp'. the latter method allows for copying of 'pomp' objects with modification.
Modified: pkg/NAMESPACE
===================================================================
--- pkg/NAMESPACE 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/NAMESPACE 2011-05-09 19:58:28 UTC (rev 463)
@@ -42,6 +42,7 @@
)
exportMethods(
+ pomp,
plot,show,print,coerce,summary,logLik,window,"$",
dprocess,rprocess,rmeasure,dmeasure,init.state,skeleton,
data.array,obs,coef,"coef<-",time,"time<-",timezero,"timezero<-",
@@ -55,7 +56,6 @@
)
export(
- pomp,
as.data.frame.pomp,
reulermultinom,
deulermultinom,
Modified: pkg/R/pomp.R
===================================================================
--- pkg/R/pomp.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/R/pomp.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -21,11 +21,13 @@
paramnames = 'character',
covarnames = 'character',
PACKAGE = 'character',
- userdata = 'list',
- call = "call"
+ userdata = 'list'
)
)
+## this is the initial-condition setting function that is used by default
+## it simply finds all parameters in the vector 'params' that have a name ending in '.0'
+## and returns a vector with their values with names stripped of '.0'
default.initializer <- function (params, t0, ...) {
ivpnames <- grep("\\.0$",names(params),val=TRUE)
if (length(ivpnames)<1)
@@ -35,15 +37,17 @@
x
}
-## constructor of the pomp class
-pomp <- function (data, times, t0, ..., rprocess, dprocess,
- rmeasure, dmeasure, measurement.model,
- skeleton.map, skeleton.vectorfield, initializer, covar, tcovar,
- obsnames, statenames, paramnames, covarnames,
- PACKAGE) {
- ## save the call
- this.call <- match.call()
+## as of version 0.37-1 'pomp' is a generic function
+setGeneric("pomp",function(data,...)standardGeneric("pomp"))
+## basic constructor of the pomp class
+pomp.constructor <- function (data, times, t0, ..., rprocess, dprocess,
+ rmeasure, dmeasure, measurement.model,
+ skeleton = NULL, skeleton.type = c("map","vectorfield"),
+ initializer, covar, tcovar,
+ obsnames, statenames, paramnames, covarnames,
+ PACKAGE) {
+
## check the data
if (is.data.frame(data)) {
if (!is.character(times) || length(times)!=1 || !(times%in%names(data)))
@@ -96,21 +100,12 @@
if (missing(dmeasure))
dmeasure <- function(y,x,t,params,log,covars,...)stop(sQuote("dmeasure")," not specified")
- if (missing(skeleton.map)) {
- if (missing(skeleton.vectorfield)) {# skeleton is unspecified
- skeleton.type <- as.character(NA)
- skeleton <- pomp.fun(f=function(x,t,params,covars,...)stop(sQuote("skeleton")," not specified"))
- } else { # skeleton is a vectorfield (ordinary differential equation)
- skeleton.type <- "vectorfield"
- skeleton <- pomp.fun(f=skeleton.vectorfield,PACKAGE=PACKAGE,proto=quote(skeleton.vectorfield(x,t,params,...)))
- }
+ skeleton.type <- match.arg(skeleton.type)
+
+ if (is.null(skeleton)) {
+ skeleton <- pomp.fun(f=function(x,t,params,covars,...)stop(sQuote("skeleton")," not specified"))
} else {
- if (missing(skeleton.vectorfield)) { # skeleton is a map (discrete-time system)
- skeleton.type <- "map"
- skeleton <- pomp.fun(f=skeleton.map,PACKAGE=PACKAGE,proto=quote(skeleton.map(x,t,params,...)))
- } else { # a dynamical system cannot be both a map and a vectorfield
- stop("pomp error: it is not permitted to specify both ",sQuote("skeleton.map")," and ",sQuote("skeleton.vectorfield"))
- }
+ skeleton <- pomp.fun(f=skeleton,PACKAGE=PACKAGE,proto=quote(skeleton(x,t,params,...)))
}
if (missing(initializer)) {
@@ -243,8 +238,7 @@
paramnames = paramnames,
covarnames = covarnames,
PACKAGE = PACKAGE,
- userdata = list(...),
- call=this.call
+ userdata = list(...)
)
}
@@ -327,3 +321,223 @@
)
}
+## deal with the change-over from 'skeleton.map'/'skeleton.vectorfield' to
+## 'skeleton'/'skeleton.type'
+skeleton.jigger <- function (skeleton = NULL, skeleton.type,
+ skeleton.map = NULL, skeleton.vectorfield = NULL) {
+ if (!is.null(skeleton.map)) {
+ warning(
+ "the arguments ",sQuote("skeleton.map")," and ",sQuote("skeleton.vectorfield"),
+ " are now deprecated.",
+ " Use ",sQuote("skeleton")," and ",sQuote("skeleton.type")," instead.",
+ call.=FALSE
+ )
+ if (!is.null(skeleton.vectorfield))
+ stop("pomp error: it is not permitted to specify both ",sQuote("skeleton.vectorfield")," and ",sQuote("skeleton.map"))
+ if (!is.null(skeleton))
+ stop("pomp error: it is not permitted to specify both ",sQuote("skeleton.map")," and ",sQuote("skeleton"))
+ skeleton.type <- "map"
+ skeleton <- skeleton.map
+ }
+ if (!is.null(skeleton.vectorfield)) {
+ warning(
+ "the arguments ",sQuote("skeleton.map")," and ",sQuote("skeleton.vectorfield"),
+ " are now deprecated.",
+ " Use ",sQuote("skeleton")," and ",sQuote("skeleton.type")," instead.",
+ call.=FALSE
+ )
+ if (!is.null(skeleton.map))
+ stop("pomp error: it is not permitted to specify both ",sQuote("skeleton.vectorfield")," and ",sQuote("skeleton.map"))
+ if (!is.null(skeleton))
+ stop("pomp error: it is not permitted to specify both ",sQuote("skeleton.vectorfield")," and ",sQuote("skeleton"))
+ skeleton.type <- "vectorfield"
+ skeleton <- skeleton.vectorfield
+ }
+ list(fn=skeleton,type=skeleton.type)
+}
+
+
+setMethod(
+ "pomp",
+ signature(data="data.frame"),
+ function (data, times, t0, ..., rprocess, dprocess,
+ rmeasure, dmeasure, measurement.model,
+ skeleton = NULL, skeleton.type = c("map","vectorfield"),
+ skeleton.map = NULL, skeleton.vectorfield = NULL,
+ initializer, covar, tcovar,
+ obsnames, statenames, paramnames, covarnames,
+ PACKAGE) {
+ skel <- skeleton.jigger(
+ skeleton=skeleton,
+ skeleton.type=skeleton.type,
+ skeleton.map=skeleton.map,
+ skeleton.vectorfield=skeleton.vectorfield
+ )
+ pomp.constructor(
+ data=data,
+ times=times,
+ t0=t0,
+ rprocess=rprocess,
+ dprocess=dprocess,
+ rmeasure=rmeasure,
+ dmeasure=dmeasure,
+ measurement.model=measurement.model,
+ skeleton=skel$fn,
+ skeleton.type=skel$type,
+ initializer=initializer,
+ covar=covar,
+ tcovar=tcovar,
+ obsnames=obsnames,
+ statenames=statenames,
+ paramnames=paramnames,
+ covarnames=covarnames,
+ PACKAGE=PACKAGE,
+ ...
+ )
+ }
+ )
+
+setMethod(
+ "pomp",
+ signature(data="matrix"),
+ function (data, times, t0, ..., rprocess, dprocess,
+ rmeasure, dmeasure, measurement.model,
+ skeleton = NULL, skeleton.type = c("map","vectorfield"),
+ skeleton.map = NULL, skeleton.vectorfield = NULL,
+ initializer, covar, tcovar,
+ obsnames, statenames, paramnames, covarnames,
+ PACKAGE) {
+ skel <- skeleton.jigger(
+ skeleton=skeleton,
+ skeleton.type=skeleton.type,
+ skeleton.map=skeleton.map,
+ skeleton.vectorfield=skeleton.vectorfield
+ )
+ pomp.constructor(
+ data=data,
+ times=times,
+ t0=t0,
+ rprocess=rprocess,
+ dprocess=dprocess,
+ rmeasure=rmeasure,
+ dmeasure=dmeasure,
+ measurement.model=measurement.model,
+ skeleton=skel$fn,
+ skeleton.type=skel$type,
+ initializer=initializer,
+ covar=covar,
+ tcovar=tcovar,
+ obsnames=obsnames,
+ statenames=statenames,
+ paramnames=paramnames,
+ covarnames=covarnames,
+ PACKAGE=PACKAGE,
+ ...
+ )
+ }
+ )
+
+
+setMethod(
+ "pomp",
+ signature(data="numeric"),
+ function (data, times, t0, ..., rprocess, dprocess,
+ rmeasure, dmeasure, measurement.model,
+ skeleton = NULL, skeleton.type = c("map","vectorfield"),
+ skeleton.map = NULL, skeleton.vectorfield = NULL,
+ initializer, covar, tcovar,
+ obsnames, statenames, paramnames, covarnames,
+ PACKAGE) {
+ skel <- skeleton.jigger(
+ skeleton=skeleton,
+ skeleton.type=skeleton.type,
+ skeleton.map=skeleton.map,
+ skeleton.vectorfield=skeleton.vectorfield
+ )
+ pomp.constructor(
+ data=matrix(data,nrow=1,ncol=length(data)),
+ times=times,
+ t0=t0,
+ rprocess=rprocess,
+ dprocess=dprocess,
+ rmeasure=rmeasure,
+ dmeasure=dmeasure,
+ measurement.model=measurement.model,
+ skeleton=skel$fn,
+ skeleton.type=skel$type,
+ initializer=initializer,
+ covar=covar,
+ tcovar=tcovar,
+ obsnames=obsnames,
+ statenames=statenames,
+ paramnames=paramnames,
+ covarnames=covarnames,
+ PACKAGE=PACKAGE,
+ ...
+ )
+ }
+ )
+
+setMethod(
+ "pomp",
+ signature(data="pomp"),
+ function (data, times, t0, ..., rprocess, dprocess,
+ rmeasure, dmeasure, measurement.model,
+ skeleton, skeleton.type,
+ initializer, covar, tcovar,
+ obsnames, statenames, paramnames, covarnames,
+ PACKAGE) {
+ mmg <- !missing(measurement.model)
+ dmg <- !missing(dmeasure)
+ rmg <- !missing(rmeasure)
+ if (missing(times)) times <- data at times
+ if (missing(t0)) t0 <- data at t0
+ if (mmg) {
+ if (dmg||rmg)
+ warning(
+ "specifying ",sQuote("measurement.model"),
+ " overrides specification of ",sQuote("rmeasure")," and ",sQuote("dmeasure")
+ )
+ mm <- measform2pomp(measurement.model)
+ rmeasure <- mm$rmeasure
+ dmeasure <- mm$dmeasure
+ } else {
+ if (!rmg) rmeasure <- data at rmeasure
+ if (!dmg) dmeasure <- data at dmeasure
+ }
+ if (missing(rprocess)) rprocess <- data at rprocess
+ if (missing(dprocess)) dprocess <- data at dprocess
+ if (missing(initializer)) initializer <- data at initializer
+ if (missing(covar)) covar <- data at covar
+ if (missing(tcovar)) tcovar <- data at tcovar
+ if (missing(obsnames)) obsnames <- data at obsnames
+ if (missing(statenames)) statenames <- data at statenames
+ if (missing(paramnames)) paramnames <- data at paramnames
+ if (missing(covarnames)) covarnames <- data at covarnames
+ if (missing(PACKAGE)) PACKAGE <- data at PACKAGE
+ if (missing(skeleton.type)) skeleton.type <- data at skeleton.type
+ if (missing(skeleton)) skeleton <- data at skeleton
+
+ pomp.constructor(
+ data=data at data,
+ times=times,
+ t0=t0,
+ rprocess=rprocess,
+ dprocess=dprocess,
+ rmeasure=rmeasure,
+ dmeasure=dmeasure,
+ skeleton=skeleton,
+ skeleton.type=skeleton.type,
+ initializer=initializer,
+ covar=covar,
+ tcovar=tcovar,
+ obsnames=obsnames,
+ statenames=statenames,
+ paramnames=paramnames,
+ covarnames=covarnames,
+ PACKAGE=PACKAGE,
+ ...
+ )
+ }
+ )
+
Modified: pkg/data/blowflies.rda
===================================================================
(Binary files differ)
Modified: pkg/data/dacca.rda
===================================================================
(Binary files differ)
Modified: pkg/data/euler.sir.rda
===================================================================
(Binary files differ)
Modified: pkg/data/gillespie.sir.rda
===================================================================
(Binary files differ)
Modified: pkg/data/gompertz.rda
===================================================================
(Binary files differ)
Modified: pkg/data/ou2.rda
===================================================================
(Binary files differ)
Modified: pkg/data/ricker.rda
===================================================================
(Binary files differ)
Modified: pkg/data/rw2.rda
===================================================================
(Binary files differ)
Modified: pkg/data/verhulst.rda
===================================================================
(Binary files differ)
Modified: pkg/inst/ChangeLog
===================================================================
--- pkg/inst/ChangeLog 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/ChangeLog 2011-05-09 19:58:28 UTC (rev 463)
@@ -1,5 +1,17 @@
+2011-05-09 kingaa
+
+ * [r462] R/pomp-methods.R: - remove printing of call slot
+ * [r461] src/gompertz.c: - better commenting and some cleanup
+ * [r460] R/pomp-fun.R: - fix bug in error-checking
+ - allow pass-through of 'pomp.fun' objects
+ * [r459] DESCRIPTION, R/simulate-pomp.R, R/trajectory-pomp.R: -
+ remove the warnings that were added in 0.34-1 regarding changed
+ default behavior of 'simulate' and 'trajectory'
+
2011-05-05 kingaa
+ * [r458] inst/ChangeLog, inst/doc/advanced_topics_in_pomp.pdf,
+ inst/doc/intro_to_pomp.pdf:
* [r457] DESCRIPTION, data/blowflies.rda, data/dacca.rda,
inst/ChangeLog, inst/data-R/blowflies.R,
inst/data-R/blowflies.csv, inst/data-R/dacca.R,
Modified: pkg/inst/NEWS
===================================================================
--- pkg/inst/NEWS 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/NEWS 2011-05-09 19:58:28 UTC (rev 463)
@@ -1,16 +1,29 @@
-NEW IN VERSION 0.36-7:
+NEWS
+0.37-1
+ o The arguments 'skeleton.map' and 'skeleton.vectorfield' are now deprecated.
+ Specify a discrete-time skeleton using 'skeleton.type="map"' and a continuous-time skeleton via 'skeleton.type="vectorfield"'.
+
+ o 'pomp' is now a method rather than a function.
+ Calling it on a 'pomp' object copies the object; one can optionally modify some of the slots.
+
+ o The warning messages to do with the default behaviors of 'simulate' and 'trajectory' with respect to the arguments 'times' and 't0' have been removed.
+ These warnings were introduced in version 0.34-1 to alert users to changes in the default behavior of 'simulate' and 'trajectory' introduced with that version.
+ See the documentation ("pomp?simulate", "pomp?trajectory") for a description of the new default behaviors.
+
+0.36-7:
+
o 'trajectory' now gives a more informative error when no skeleton is present
o 'pred.mean', 'pred.var', and 'filter.mean' methods are now defined for 'pfilterd.pomp' objects (and not just 'mif' objects)
o A new 'data()'-loadable example is included: do '?blowflies' to read about it.
-NEW IN VERSION 0.36-6:
+0.36-6:
o The B-spline basis construction functions 'bspline.basis' and 'periodic.bspline.basis' now take a 'names' argument that allows one to give names to the basis functions.
-NEW IN VERSION 0.36-5:
+0.36-5:
o The "intro_to_pomp" vignette has been updated. The first part of the vignette centers on a new, more transparent first example.
@@ -25,7 +38,7 @@
o The rprocess plugin 'discrete.time.sim' now takes an optional 'delta.t' argument to allow for arbitrary (but still necessarily regular) spacing of observations. By default, delta.t=1.
-NEW IN VERSION 0.36-1:
+0.36-1:
o There are two new classes: 'traj.matched.pomp' and 'pfilterd.pomp' to hold the results of 'traj.match' and 'pfilter' computations, respectively. The old classes 'mif' and 'pmcmc' now inherit from 'pfilterd.pomp' in a natural way.
@@ -37,7 +50,7 @@
o There is a change to the interface to 'pfilter' in that 'save.states' results in the filling of the 'last.states' slot. Before version 0.36-1, the 'pfilter' returned a list. The element 'states' of that list corresponds to the slot 'last.states' in version 0.36-1. It was necessary to make this name-change in order to avoid a conflict with the 'states' slot inherited from the 'pomp' class.
-NEW IN VERSION 0.35-1:
+0.35-1:
o Added a facility for computation of "synthetic likelihood" (Wood 2010) using a robust estimate of the covariance. A new slot in objects of class 'probed.pomp' holds this value.
@@ -53,7 +66,7 @@
o The 'seed' slot in 'probed.pomp' and 'spect.pomp' objects has been removed, along with unwanted behavior to do with it.
-NEW IN VERSION 0.34-2:
+0.34-2:
o Added a method to coerce 'probed.pomp' objects to 'data.frame' objects.
@@ -62,7 +75,7 @@
o 'simulate' and 'rprocess' are accelerated somewhat through judicious use of 'memcpy'.
-NEW IN VERSION 0.34-1:
+0.34-1:
o The default behaviors of 'simulate' and 'trajectory' with respect to the argument 'times' has changed.
There is a new argument, 't0', which specifies the time at which the initial conditions obtain.
@@ -92,7 +105,7 @@
o The functions 'euler.simulate', 'onestep.simulate', and 'onestep.density', deprecated since version 0.29-1, have been removed.
-NEW IN VERSION 0.33-1:
+0.33-1:
o Major improvements in speed in the probe-matching routines have been realized by coding the computationally-intensive portions of these calculations in C.
@@ -100,31 +113,31 @@
In particular, the new function 'probe.marginal' regresses the marginal distributions of actual and simulated data against a reference distribution, returning as probes the regression coefficients.
'probe.marginal' and 'probe.acf' have been coded in C for speed.
-NEW IN VERSION 0.32-1:
+0.32-1:
o pomp now includes parameter estimation methods based on probe-matching and power-spectrum-matching as explained in papers by Kendall et al. (B. E. Kendall, C. J. Briggs, W. M. Murdoch, P. Turchin, S. P. Ellner, E. McCauley, R. M. Nisbet, S. N. Wood (1999) Why do populations cycle? A synthesis of statistical and mechanistic modeling approaches. Ecology, 80:1789--1805) and Reuman et al. (D.C. Reuman, R.A. Desharnais, R.F. Costantino, O. Ahmad, J.E. Cohen (2006) Power spectra reveal the influence of stochasticity on nonlinear population dynamics. Proc. Nat'l. Acad. Sci., U.S.A.} 103:18860--18865).
-NEW IN VERSION 0.31-2:
+0.31-2:
o new "window", "timezero", and "timezero<-" methods are available for manipulating 'pomp' objects.
-NEW IN VERSION 0.31-1:
+0.31-1:
o pomp now includes an implementation of the particle Markov chain Monte Carlo algorithm of Andrieu et al. (Andrieu, C., Doucet, A., & Holenstein, R. (2010) Particle Markov chain Monte Carlo methods. J. R. Stat. Soc. B 72:1--33).
This algorithm is implemented in the method 'pmcmc': see 'help(pmcmc)' for details.
Implementation primarily due to Ed Ionides.
-NEW IN VERSION 0.30-1:
+0.30-1:
o pomp now includes an implementation of the approximate Bayesian Sequential Monte Carlo algorithm of Liu & West (Liu, J. & West, M. (2001) Combining Parameter and State Estimation in Simulation-Based Filtering. In Doucet, A.; de Freitas, N. & Gordon, N. J. (eds.) Sequential Monte Carlo Methods in Practice, Springer, New York, pp. 197--224).
This algorithm is implemented in the method 'bsmc': see 'help(bsmc)' for details.
Thanks to Matt Ferrari and Michael Lavine for the implementation.
-NEW IN VERSION 0.29-5:
+0.29-5:
o The cholera transmission model from King, A. A., Ionides, E. L., Pascual, M., and Bouma, M. J. (2008) Inapparent infections and cholera dynamics, Nature 454:877--880 together with data from the Dacca district of historic Bengal are now included.
See 'help(dacca)' for details.
-NEW IN VERSION 0.29-1:
+0.29-1:
o This is a major revision.
@@ -145,41 +158,41 @@
o The changes are documented in detail in the file CHANGES_0.29-1.txt.
-NEW IN VERSION 0.22-1:
+0.22-1:
o The nonlinear forecasting method of Kendall, Ellner, et al. is now implemented in the package. This is an "indirect inference" method using quasi-likelihood as an objective function.
-NEW IN VERSION 0.21-4:
+0.21-4:
o Improved documentation in "euler.Rd".
-NEW IN VERSION 0.21-3:
+0.21-3:
o Trajectories of the deterministic skeleton can now be computed via the method "trajectory".
o Bug fix in stochastic Euler algorithm.
-NEW IN VERSION 0.20-8:
+0.20-8:
o Better error handling in mif
o Fix inaccuracy in documentation of 'pomp' arguments 'rmeasure' and 'dmeasure'
-NEW IN VERSION 0.20-7:
+0.20-7:
o Bug fix in convergence record storage
-NEW IN VERSION 0.20-6:
+0.20-6:
o Bug fix in "intro_to_pomp" vignette
-NEW IN VERSION 0.20-5:
+0.20-5:
o Better error handling and improved documentation for 'mif'
o Testing of examples in the 'examples' directory, including compilation, linking, and running of dynamically-loadable library on unix machines
-NEW IN VERSION 0.20-4:
+0.20-4:
o A bug in the Euler simulator was repaired.
@@ -188,7 +201,7 @@
o A potential bug in the Euler-multinomial simulator was repaired.
-NEW IN VERSION 0.20-2:
+0.20-2:
o The lookup_table (see file 'include/pomp.h') interface for C codes has been streamlined.
It is now slightly more straightforward to set up a 'struct lookup_table'.
@@ -202,7 +215,7 @@
o Several bugs fixed.
-NEW IN VERSION 0.20-1
+0.20-1
o This is a major revision.
Modified: pkg/inst/data-R/euler.sir.R
===================================================================
--- pkg/inst/data-R/euler.sir.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/data-R/euler.sir.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -13,7 +13,8 @@
delta.t=1/52/20,
PACKAGE="pomp"
),
- skeleton.vectorfield="_sir_ODE",
+ skeleton.type="vectorfield",
+ skeleton="_sir_ODE",
rmeasure="_sir_binom_rmeasure",
dmeasure="_sir_binom_dmeasure",
PACKAGE="pomp",
Modified: pkg/inst/data-R/gompertz.R
===================================================================
--- pkg/inst/data-R/gompertz.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/data-R/gompertz.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -10,7 +10,8 @@
),
rmeasure="_gompertz_normal_rmeasure",
dmeasure="_gompertz_normal_dmeasure",
- skeleton.map="_gompertz_skeleton",
+ skeleton.type="map",
+ skeleton="_gompertz_skeleton",
paramnames=c("log.r","log.K","log.sigma","log.tau"),
statenames=c("X"),
obsnames=c("Y")
Modified: pkg/inst/data-R/ou2.R
===================================================================
--- pkg/inst/data-R/ou2.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/data-R/ou2.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -56,7 +56,8 @@
},
dmeasure = "_ou2_normal_dmeasure",
rmeasure = "_ou2_normal_rmeasure",
- skeleton.map = function (x, t, params, ...) {
+ skeleton.type="map",
+ skeleton = function (x, t, params, ...) {
with(
as.list(c(x,params)),
{
Modified: pkg/inst/data-R/ricker.R
===================================================================
--- pkg/inst/data-R/ricker.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/data-R/ricker.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -10,10 +10,12 @@
),
rmeasure="_ricker_poisson_rmeasure",
dmeasure="_ricker_poisson_dmeasure",
- skeleton.map="_ricker_skeleton",
+ skeleton.type="map",
+ skeleton="_ricker_skeleton",
paramnames=c("log.r","log.sigma","log.phi"),
statenames=c("N","e"),
- obsnames=c("y")
+ obsnames=c("y"),
+ PACKAGE="pomp"
),
params=c(
log.r=3.8,
Modified: pkg/inst/data-R/verhulst.R
===================================================================
--- pkg/inst/data-R/verhulst.R 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/data-R/verhulst.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -33,7 +33,8 @@
}
),
measurement.model=obs~lnorm(meanlog=log(n),sdlog=log(1+tau)),
- skeleton.vectorfield=function(x,t,params,...){
+ skeleton.type="vectorfield",
+ skeleton=function(x,t,params,...){
with(
as.list(c(x,params)),
r*n*(1-n/K)
Modified: pkg/inst/doc/advanced_topics_in_pomp.Rnw
===================================================================
--- pkg/inst/doc/advanced_topics_in_pomp.Rnw 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/doc/advanced_topics_in_pomp.Rnw 2011-05-09 19:58:28 UTC (rev 463)
@@ -87,7 +87,8 @@
delta.t=1/52/20
),
## native routine for the skeleton:
- skeleton.vectorfield="_sir_ODE",
+ skeleton.type="vectorfield",
+ skeleton="_sir_ODE",
## binomial measurement model:
rmeasure="_sir_binom_rmeasure",
dmeasure="_sir_binom_dmeasure",
Modified: pkg/inst/doc/advanced_topics_in_pomp.pdf
===================================================================
(Binary files differ)
Modified: pkg/inst/doc/gompertz-trajmatch.rda
===================================================================
(Binary files differ)
Modified: pkg/inst/doc/intro_to_pomp.Rnw
===================================================================
--- pkg/inst/doc/intro_to_pomp.Rnw 2011-05-09 19:52:51 UTC (rev 462)
+++ pkg/inst/doc/intro_to_pomp.Rnw 2011-05-09 19:58:28 UTC (rev 463)
@@ -207,7 +207,7 @@
X <- x["X"]
## generate a log-normal random variable:
eps <- exp(rnorm(n=1,mean=0,sd=sigma))
- ## compute the state at time t+1:
+ ## compute the state at time t+delta.t:
S <- exp(-r*delta.t)
xnew <- c(X=unname(K^(1-S)*X^S*eps))
return(xnew)
@@ -352,20 +352,13 @@
To do so, we'll need to incorporate it into the \pomp\ object.
We can do this by specifying the \code{dmeasure} argument in another call to \pomp:
<<second-pomp-construction>>=
-dat <- as(gompertz,"data.frame")
-theta <- coef(gompertz)
gompertz <- pomp(
- data=dat[c("time","Y")],
- times="time",
- rprocess=discrete.time.sim(gompertz.proc.sim),
- rmeasure=gompertz.meas.sim,
- dmeasure=gompertz.meas.dens,
- t0=0
+ gompertz,
+ dmeasure=gompertz.meas.dens
)
coef(gompertz) <- theta
@
-Note that we've first extracted the data from our old \code{gompertz} and set up the new one with the same data and parameters.
-The calls to \code{coef} and \code{coef<-} in the lines above make sure the parameters have the same values they had before.
+The result of the above is a new \pomp\ object \code{gompertz} in every way identical to the one we had before, but with the measurement-model density function \code{dmeasure} now specified.
To compute the likelihood of the data, we can use the function \code{pfilter}.
This runs a plain vanilla particle filter (AKA sequential Monte Carlo) algorithm and results in an unbiased estimate of the likelihood.
@@ -527,7 +520,7 @@
X <- x["X"]
## generate a log-normal random variable:
eps <- exp(rnorm(n=1,mean=0,sd=sigma))
- ## compute the state at time t+1:
+ ## compute the state at time t+delta.t:
S <- exp(-r*delta.t)
xnew <- c(X=unname(K^(1-S)*X^S*eps))
return(xnew)
@@ -581,6 +574,12 @@
@
A \code{pomp} object corresponding to the one just created (but with the \code{rprocess}, \code{rmeasure}, and \code{dmeasure} bits coded in C for speed) can be loaded by executing \verb+data(gompertz)+.
+For your convenience, codes creating this \pomp\ object are included with the package.
+To view them, do:
+<<eval=F>>=
+file.show(system.file("examples/gompertz.R",package="pomp"))
+file.show(system.file("examples/gompertz.c",package="pomp"))
+@
\clearpage
\section{Estimating parameters using iterated filtering: \code{mif}}
@@ -718,7 +717,7 @@
The \code{skeleton} slot should be filled with the deterministic skeleton of the process model.
In the discrete-time case, this is the map
\begin{equation*}
- x\,\mapsto\,\expect{X_{t+1}\;\vert\;X_{t}=x,\theta}.
+ x\,\mapsto\,\expect{X_{t+{\Delta}t}\;\vert\;X_{t}=x,\theta}.
\end{equation*}
In the continuous-time case, this is the vectorfield
\begin{equation*}
@@ -750,7 +749,8 @@
rprocess=discrete.time.sim(gompertz.proc.sim),
rmeasure=gompertz.meas.sim,
dmeasure=gompertz.meas.dens,
- skeleton.map=gompertz.skel,
+ skeleton.type="map",
+ skeleton=gompertz.skel,
t0=0
)
coef(new.gompertz) <- theta
@@ -760,7 +760,7 @@
coef(new.gompertz,"log.sigma") <- log(0.05)
new.gompertz <- simulate(new.gompertz,seed=88737400L)
@
-We use the \code{skeleton.map} argument for discrete-time processes and \code{skeleton.vectorfield} for continuous-time processes.
+We use \code{skeleton.type="map"} for discrete-time processes (such as the Gompertz model) and \code{skeleton.type="vectorfield"} for continuous-time processes.
%% Note that we have turned off the process noise in \code{new.gompertz} (next to last line) so that trajectory matching is actually formally appropriate for this model.
The \pomp\ function \code{traj.match} calls the optimizer \code{optim} to minimize the discrepancy between the trajectory and the data.
Modified: pkg/inst/doc/intro_to_pomp.pdf
===================================================================
(Binary files differ)
Added: pkg/inst/examples/gompertz.R
===================================================================
--- pkg/inst/examples/gompertz.R (rev 0)
+++ pkg/inst/examples/gompertz.R 2011-05-09 19:58:28 UTC (rev 463)
@@ -0,0 +1,111 @@
+require(pomp)
+
+## First code up the Gompertz example in R:
+
+pomp(
+ data=data.frame(time=1:100,Y=NA),
+ times="time",
+ t0=0,
+ rprocess=discrete.time.sim( # a discrete-time process (see ?plugins)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/pomp -r 463
More information about the pomp-commits
mailing list