[Pomp-commits] r1204 - in pkg/mif2: . R inst man tests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Jun 15 15:13:10 CEST 2015
Author: kingaa
Date: 2015-06-15 15:13:10 +0200 (Mon, 15 Jun 2015)
New Revision: 1204
Removed:
pkg/mif2/.Rbuildignore
pkg/mif2/.Rinstignore
pkg/mif2/DESCRIPTION
pkg/mif2/NAMESPACE
pkg/mif2/R/generics.R
pkg/mif2/R/mif2-methods.R
pkg/mif2/R/mif2.R
pkg/mif2/inst/NEWS
pkg/mif2/inst/NEWS.Rd
pkg/mif2/man/mif2-methods.Rd
pkg/mif2/man/mif2.Rd
pkg/mif2/tests/ou2-mif2.R
pkg/mif2/tests/ou2-mif2.Rout.save
Log:
- remove mif2 as separate package
Deleted: pkg/mif2/.Rbuildignore
===================================================================
--- pkg/mif2/.Rbuildignore 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/.Rbuildignore 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,4 +0,0 @@
-inst/doc/Makefile
-inst/doc/(.+?)\.bst$
-inst/doc/(.+?)\.R$
-inst/doc/(.+?)\.png$
Deleted: pkg/mif2/.Rinstignore
===================================================================
--- pkg/mif2/.Rinstignore 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/.Rinstignore 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,3 +0,0 @@
-inst/doc/Makefile
-inst/doc/fullnat.bst
-inst/doc/(.+?)\.rda$
Deleted: pkg/mif2/DESCRIPTION
===================================================================
--- pkg/mif2/DESCRIPTION 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/DESCRIPTION 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,18 +0,0 @@
-Package: mif2
-Type: Package
-Title: MIF2 algorithm development package
-Version: 1.0-4
-Date: 2014-07-15
-Authors at R: c(person(given=c("Aaron","A."),family="King",
- role=c("aut","cre"),email="kingaa at umich.edu"),
- person(given=c("Edward","L."),family="Ionides",role=c("aut")),
- person(given="Dao",family="Nguyen",role=c("ctb"))
- )
-URL: http://pomp.r-forge.r-project.org
-Description: MIF2 algorithm: for eventual integration into pomp
-Depends: R(>= 3.0.0), pomp(>= 0.53-1)
-Imports: methods
-License: GPL(>= 2)
-LazyData: true
-BuildVignettes: true
-Collate: generics.R mif2.R mif2-methods.R
Deleted: pkg/mif2/NAMESPACE
===================================================================
--- pkg/mif2/NAMESPACE 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/NAMESPACE 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,6 +0,0 @@
-import(pomp)
-import(methods)
-
-exportClasses(mif2d.pomp)
-exportMethods(mif2,continue,plot,conv.rec,logLik)
-export()
Deleted: pkg/mif2/R/generics.R
===================================================================
--- pkg/mif2/R/generics.R 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/R/generics.R 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1 +0,0 @@
-setGeneric('mif2',function(object,...)standardGeneric("mif2"))
Deleted: pkg/mif2/R/mif2-methods.R
===================================================================
--- pkg/mif2/R/mif2-methods.R 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/R/mif2-methods.R 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,17 +0,0 @@
-## this file contains short definitions of methods for the 'mif2d.pomp' class
-
-setMethod('logLik','mif2d.pomp',function(object,...)object at loglik)
-
-setMethod('conv.rec','mif2d.pomp',
- function (object, pars, transform = FALSE, ...) {
- conv.rec.internal(object=object,pars=pars,transform=transform,...)
- }
- )
-
-setMethod(
- "plot",
- "mif2d.pomp",
- function (x, y = NULL, ...) {
- compare.mif(x)
- }
- )
Deleted: pkg/mif2/R/mif2.R
===================================================================
--- pkg/mif2/R/mif2.R 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/R/mif2.R 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,466 +0,0 @@
-## MIF algorithm functions
-
-## define the mif class
-setClass(
- 'mif2d.pomp',
- contains='pfilterd.pomp',
- representation=representation(
- transform = "logical",
- Nmif = 'integer',
- perturb.fn = 'function',
- conv.rec = 'matrix'
- )
- )
-
-
-mif2.pfilter <- function (object, params, Np,
- tol, max.fail,
- pred.mean, pred.var, filter.mean,
- mifiter, perturb.fn,
- transform, verbose,
- .getnativesymbolinfo = TRUE) {
-
- ptsi.inv <- ptsi.for <- gnsi.rproc <- gnsi.dmeas <- as.logical(.getnativesymbolinfo)
- transform <- as.logical(transform)
- mifiter <- as.integer(mifiter)
- Np <- as.integer(Np)
- Np <- c(Np,Np[1])
-
- times <- time(object,t0=TRUE)
- ntimes <- length(times)-1
-
- paramnames <- rownames(params)
- npars <- nrow(params)
-
- ## perturb parameters
- params <- apply(params,2L,perturb.fn,mifiter=mifiter,timeno=1L)
-
- ## transform parameters if necessary
- if (transform) {
- tparams <- partrans(object,params,dir="forward",
- .getnativesymbolinfo=ptsi.for)
- ptsi.for <- FALSE
- }
-
- ## get initial states
- x <- init.state(
- object,
- params=if (transform) tparams else params,
- )
- statenames <- rownames(x)
- nvars <- nrow(x)
-
- loglik <- rep(NA,ntimes)
- eff.sample.size <- numeric(ntimes)
- nfail <- 0
-
- ## set up storage for prediction means, variances, etc.
- if (pred.mean)
- pred.m <- matrix(
- data=0,
- nrow=nvars+npars,
- ncol=ntimes,
- dimnames=list(c(statenames,paramnames),NULL)
- )
- else
- pred.m <- array(data=numeric(0),dim=c(0,0))
-
- if (pred.var)
- pred.v <- matrix(
- data=0,
- nrow=nvars+npars,
- ncol=ntimes,
- dimnames=list(c(statenames,paramnames),NULL)
- )
- else
- pred.v <- array(data=numeric(0),dim=c(0,0))
-
- if (filter.mean)
- filt.m <- matrix(
- data=0,
- nrow=nvars+length(paramnames),
- ncol=ntimes,
- dimnames=list(c(statenames,paramnames),NULL)
- )
- else
- filt.m <- array(data=numeric(0),dim=c(0,0))
-
- for (nt in seq_len(ntimes)) {
-
- if (nt > 1) {
- ## perturb parameters
- params <- apply(params,2L,perturb.fn,mifiter=mifiter,timeno=nt)
- ## transform parameters if necessary
- if (transform) {
- tparams <- partrans(object,params,dir="forward",
- .getnativesymbolinfo=ptsi.for)
- }
- }
-
- ## advance the state variables according to the process model
- X <- try(
- rprocess(
- object,
- xstart=x,
- times=times[c(nt,nt+1)],
- params=if (transform) tparams else params,
- offset=1,
- .getnativesymbolinfo=gnsi.rproc
- ),
- silent=FALSE
- )
- if (inherits(X,'try-error'))
- stop(sQuote("mif2.pfilter")," error: process simulation error")
- gnsi.rproc <- FALSE
-
- ## determine the weights
- weights <- try(
- dmeasure(
- object,
- y=object at data[,nt,drop=FALSE],
- x=X,
- times=times[nt+1],
- params=if (transform) tparams else params,
- log=FALSE,
- .getnativesymbolinfo=gnsi.dmeas
- ),
- silent=FALSE
- )
- if (inherits(weights,'try-error'))
- stop(sQuote("mif2.pfilter")," error: error in calculation of weights")
- if (any(!is.finite(weights))) {
- stop(sQuote("mif2.pfilter")," error: ",sQuote("dmeasure"),
- " returns non-finite value")
- }
- gnsi.dmeas <- FALSE
-
- ## compute prediction mean, prediction variance, filtering mean,
- ## effective sample size, log-likelihood
- ## also do resampling if filtering has not failed
- xx <- try(
- .Call(
- pomp:::pfilter_computations,
- X,params,Np[nt+1],
- FALSE,numeric(0),
- pred.mean,pred.var,filter.mean,
- FALSE,weights,tol
- ),
- silent=FALSE
- )
- if (inherits(xx,'try-error')) {
- stop(sQuote("mif2.pfilter")," error",call.=FALSE)
- }
- all.fail <- xx$fail
- loglik[nt] <- xx$loglik
- eff.sample.size[nt] <- xx$ess
-
- x <- xx$states
- params <- xx$params
-
- if (pred.mean)
- pred.m[,nt] <- xx$pm
- if (pred.var)
- pred.v[,nt] <- xx$pv
- if (filter.mean)
- filt.m[,nt] <- xx$fm
-
- if (all.fail) { ## all particles are lost
- nfail <- nfail+1
- if (verbose)
- message("filtering failure at time t = ",times[nt+1])
- if (nfail>max.fail)
- stop(sQuote("mif2.pfilter")," error: too many filtering failures",call.=FALSE)
- }
-
- if (verbose && (nt%%5==0))
- cat("mif2 pfilter timestep",nt,"of",ntimes,"finished\n")
-
- }
-
- if (nfail>0)
- warning(sprintf(ngettext(nfail,msg1="%d filtering failure occurred in ",
- msg2="%d filtering failures occurred in "),nfail),
- sQuote("mif2.pfilter"),call.=FALSE)
-
- new(
- "pfilterd.pomp",
- object,
- pred.mean=pred.m,
- pred.var=pred.v,
- filter.mean=filt.m,
- paramMatrix=params,
- eff.sample.size=eff.sample.size,
- cond.loglik=loglik,
- saved.states=list(),
- saved.params=list(),
- seed=as.integer(NULL),
- Np=as.integer(head(Np,-1)),
- tol=tol,
- nfail=as.integer(nfail),
- loglik=sum(loglik)
- )
-}
-
-mif2.internal <- function (object, Nmif,
- start, Np, perturb.fn,
- tol, max.fail,
- verbose, transform, .ndone = 0L,
- .getnativesymbolinfo = TRUE,
- ...) {
-
- gnsi <- as.logical(.getnativesymbolinfo)
- transform <- as.logical(transform)
-
- Nmif <- as.integer(Nmif)
- if (Nmif<0)
- stop("mif error: ",sQuote("Nmif")," must be a positive integer",call.=FALSE)
-
- if (transform)
- start <- partrans(object,start,dir="inverse")
-
- ntimes <- length(time(object))
-
- start.names <- names(start)
- if (is.null(start.names))
- stop("mif2 error: ",sQuote("start")," must be a named vector")
-
- if (!is.function(perturb.fn))
- stop("mif2 error: ",sQuote("perturb.fn")," must be a function")
-
- if (is.function(Np)) {
- Np <- try(
- vapply(seq.int(from=1,to=ntimes,by=1),Np,numeric(1)),
- silent=FALSE
- )
- if (inherits(Np,"try-error"))
- stop("if ",sQuote("Np")," is a function, it must return a single positive integer")
- }
- if (length(Np)==1)
- Np <- rep(Np,times=ntimes)
- else if (length(Np)!=ntimes)
- stop(sQuote("Np")," must have length 1 or length ",ntimes)
- if (any(Np<=0))
- stop("number of particles, ",sQuote("Np"),", must always be positive")
- if (!is.numeric(Np))
- stop(sQuote("Np")," must be a number, a vector of numbers, or a function")
- Np <- as.integer(Np)
-
- conv.rec <- matrix(
- data=NA,
- nrow=Nmif+1,
- ncol=length(start)+2,
- dimnames=list(
- seq(.ndone,.ndone+Nmif),
- c('loglik','nfail',names(start))
- )
- )
- conv.rec[1L,] <- c(NA,NA,start)
-
- if (.ndone > 0) { # call is from 'continue'
- paramMatrix <- object at paramMatrix
- } else if (Nmif > 0) { # initial call
- paramMatrix <- array(data=start,dim=c(length(start),Np[1L]),
- dimnames=list(names(start),NULL))
- } else { # no work to do
- paramMatrix <- array(data=numeric(0),dim=c(0,0))
- }
-
- object <- as(object,"pomp")
-
- for (n in seq_len(Nmif)) { ## iterate the filtering
-
- pfp <- try(
- mif2.pfilter(
- object=object,
- params=paramMatrix,
- Np=Np,
- tol=tol,
- max.fail=max.fail,
- pred.mean=(n==Nmif),
- pred.var=(n==Nmif),
- filter.mean=(n==Nmif),
- mifiter=.ndone+n,
- perturb.fn=perturb.fn,
- transform=transform,
- verbose=verbose,
- .getnativesymbolinfo=gnsi
- ),
- silent=FALSE
- )
- if (inherits(pfp,"try-error"))
- stop("mif2 particle-filter error")
-
- gnsi <- FALSE
-
- theta <- rowMeans(pfp at paramMatrix[,,drop=FALSE])
-
- conv.rec[n+1,-c(1,2)] <- theta
- conv.rec[n,c(1,2)] <- c(pfp at loglik,pfp at nfail)
-
- if (verbose) cat("MIF2 iteration ",n," of ",Nmif," completed\n")
-
- } ### end of main loop
-
- ## back transform the parameter estimate if necessary
- if (transform) theta <- partrans(pfp,theta,dir="forward")
-
- new(
- "mif2d.pomp",
- pfp,
- params=theta,
- paramMatrix=pfp at paramMatrix,
- transform=transform,
- Nmif=Nmif,
- perturb.fn=perturb.fn,
- conv.rec=conv.rec,
- tol=tol
- )
-}
-
-setMethod(
- "mif2",
- signature=signature(object="pomp"),
- function (object, Nmif = 1,
- start, Np, perturb.fn,
- tol = 1e-17, max.fail = Inf,
- verbose = getOption("verbose"),
- transform = FALSE,
- ...) {
-
- if (missing(start)) start <- coef(object)
- if (length(start)==0)
- stop(
- "mif2 error: ",sQuote("start")," must be specified if ",
- sQuote("coef(object)")," is NULL",
- call.=FALSE
- )
-
- if (missing(Np))
- stop("mif2 error: ",sQuote("Np")," must be specified",call.=FALSE)
-
- if (missing(perturb.fn))
- stop("mif2 error: ",sQuote("perturb.fn")," must be specified",call.=FALSE)
- perturb.fn <- match.fun(perturb.fn)
- if (!all(c('params','mifiter','timeno','...')%in%names(formals(perturb.fn)))) {
- stop(
- "mif2 error: ",
- sQuote("perturb.fn"),
- " must be a function of prototype ",
- sQuote("perturb.fn(params,mifiter,timeno,...)"),
- call.=FALSE
- )
- }
-
- mif2.internal(
- object=object,
- Nmif=Nmif,
- start=start,
- Np=Np,
- perturb.fn=perturb.fn,
- tol=tol,
- max.fail=max.fail,
- transform=transform,
- verbose=verbose,
- ...
- )
-
- }
- )
-
-
-setMethod(
- "mif2",
- signature=signature(object="pfilterd.pomp"),
- function (object, Nmif = 1, Np, tol, ...) {
-
- if (missing(Np)) Np <- object at Np
- if (missing(tol)) tol <- object at tol
-
- mif2(
- object=as(object,"pomp"),
- Nmif=Nmif,
- Np=Np,
- tol=tol,
- ...
- )
- }
- )
-
-setMethod(
- "mif2",
- signature=signature(object="mif2d.pomp"),
- function (object, Nmif, start,
- Np, perturb.fn,
- tol, transform,
- ...) {
-
- if (missing(Nmif)) Nmif <- object at Nmif
- if (missing(start)) start <- coef(object)
- if (missing(perturb.fn)) perturb.fn <- object at perturb.fn
- if (missing(transform)) transform <- object at transform
-
- if (missing(Np)) Np <- object at Np
- if (missing(tol)) tol <- object at tol
-
- mif2(
- object=as(object,"pomp"),
- Nmif=Nmif,
- start=start,
- Np=Np,
- perturb.fn=perturb.fn,
- tol=tol,
- transform=transform,
- ...
- )
- }
- )
-
-setMethod(
- 'continue',
- signature=signature(object='mif2d.pomp'),
- function (object, Nmif = 1,
- ...) {
-
- ndone <- object at Nmif
-
- obj <- mif2(
- object=object,
- Nmif=Nmif,
- .ndone=ndone,
- ...
- )
-
- object at conv.rec[ndone+1,c('loglik','nfail')] <- obj at conv.rec[1L,c('loglik','nfail')]
- obj at conv.rec <- rbind(
- object at conv.rec,
- obj at conv.rec[-1L,colnames(object at conv.rec)]
- )
- obj at Nmif <- as.integer(ndone+Nmif)
-
- obj
- }
- )
-
-default.cooling <- function (object, fraction, par.sd, ic.sd) {
-
- nT <- length(time(object))
- theta <- (1-fraction)/fraction/(50*nT-1)
-
- function (params, mifiter, timeno, ...) {
- pert <- params
- sigma <- 1/(1+theta*((mifiter-1)*nT+timeno-1))
- if (timeno==1) {
- pert[names(ic.sd)] <- rnorm(
- n=length(ic.sd),
- mean=pert[names(ic.sd)],
- sd=ic.sd*sigma
- )
- }
- pert[names(par.sd)] <- rnorm(
- n=length(par.sd),
- mean=pert[names(par.sd)],
- sd=par.sd*sigma
- )
- pert
- }
-}
Deleted: pkg/mif2/inst/NEWS
===================================================================
--- pkg/mif2/inst/NEWS 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/inst/NEWS 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,6 +0,0 @@
-_N_e_w_s _f_o_r _p_a_c_k_a_g_e '_m_i_f_2'
-
-_C_h_a_n_g_e_s _i_n '_m_i_f_2' _v_e_r_s_i_o_n _0._5_0-_1:
-
- • ‘mif2’ is now an add-on to ‘pomp’.
-
Deleted: pkg/mif2/inst/NEWS.Rd
===================================================================
--- pkg/mif2/inst/NEWS.Rd 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/inst/NEWS.Rd 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,7 +0,0 @@
-\name{NEWS}
-\title{News for package `mif2'}
-\section{Changes in \pkg{mif2} version 0.50-1}{
- \itemize{
- \item \pkg{mif2} is now an add-on to \pkg{pomp}.
- }
-}
Deleted: pkg/mif2/man/mif2-methods.Rd
===================================================================
--- pkg/mif2/man/mif2-methods.Rd 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/man/mif2-methods.Rd 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,56 +0,0 @@
-\name{mif2-methods}
-\docType{methods}
-\alias{mif2-methods}
-\alias{logLik,mif2d.pomp-method}
-\alias{logLik-mif2d.pomp}
-\alias{conv.rec,mif2d.pomp-method}
-\alias{conv.rec-mif2d.pomp}
-\alias{plot-mif2d.pomp}
-\alias{plot,mif2d.pomp-method}
-\title{Methods of the "mif2d.class" class}
-\description{Methods of the \code{mif2d.class} class.}
-\usage{
-\S4method{logLik}{mif2d.pomp}(object, \dots)
-\S4method{conv.rec}{mif2d.pomp}(object, pars, transform = FALSE, \dots)
-\S4method{plot}{mif2d.pomp}(x, y = NULL, \dots)
-}
-\arguments{
- \item{object, x}{The \code{mif2d.pomp} object.}
- \item{pars}{Names of parameters.}
- \item{y}{Ignored.}
- \item{transform}{
- optional logical;
- should the parameter transformations be applied?
- See \code{\link[=coef-pomp]{coef}} for details.
- }
- \item{\dots}{
- Further arguments (either ignored or passed to underlying functions).
- }
-}
-\section{Methods}{
- \describe{
- \item{conv.rec}{
- \code{conv.rec(object, pars = NULL)} returns the columns of the convergence-record matrix corresponding to the names in \code{pars}.
- By default, all rows are returned.
- }
- \item{logLik}{
- Returns the value in the \code{loglik} slot.
- }
- \item{plot}{
- Plots a series of diagnostic plots.
- }
- }
-}
-\references{
- E. L. Ionides, C. Bret\\'o, & A. A. King,
- Inference for nonlinear dynamical systems,
- Proc. Natl. Acad. Sci. U.S.A., 103:18438--18443, 2006.
-
- A. A. King, E. L. Ionides, M. Pascual, and M. J. Bouma,
- Inapparent infections and cholera dynamics,
- Nature, 454:877--880, 2008.
-}
-\author{Aaron A. King \email{kingaa at umich dot edu}}
-\seealso{\code{\link{mif}}, \code{\link{pomp}}, \code{\link{pomp-class}}, \code{\link{pfilter}}}
-\keyword{models}
-\keyword{ts}
Deleted: pkg/mif2/man/mif2.Rd
===================================================================
--- pkg/mif2/man/mif2.Rd 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/man/mif2.Rd 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,76 +0,0 @@
-\name{mif2}
-\docType{methods}
-\alias{mif2}
-\alias{mif2d.pomp-class}
-\alias{mif2,mif2d.pomp-method}
-\alias{mif2-mif2d.pomp}
-\alias{mif2,pfilterd.pomp-method}
-\alias{mif2-pfilterd.pomp}
-\alias{mif2,pomp-method}
-\alias{mif2-pomp}
-\alias{continue,mif2d.pomp-method}
-\alias{continue-mif2}
-\title{The MIF2 algorithm}
-\description{
- The MIF2 algorithm for estimating the parameters of a partially-observed Markov process.
-}
-\usage{
-\S4method{mif2}{pomp}(object, Nmif = 1, start, Np, perturb.fn,
- tol = 1e-17, max.fail = Inf, verbose = getOption("verbose"),
- transform = FALSE, \dots)
-\S4method{mif2}{pfilterd.pomp}(object, Nmif = 1, Np, tol, \dots)
-\S4method{mif2}{mif2d.pomp}(object, Nmif, start, Np, perturb.fn,
- tol, transform, \dots)
-\S4method{continue}{mif2d.pomp}(object, Nmif = 1, \dots)
-}
-\arguments{
- \item{object}{
- An object of class \code{pomp}.
- }
- \item{Nmif}{
- The number of MIF iterations to perform.
- }
- \item{start}{
- named numerical vector;
- the starting guess of the parameters.
- }
- \item{Np}{
- the number of particles to use in filtering.
- This may be specified as a single positive integer, in which case the same number of particles will be used at each timestep.
- Alternatively, if one wishes the number of particles to vary across timestep, one may specify \code{Np} either as a vector of positive integers (of length \code{length(time(object,t0=TRUE))}) or as a function taking a positive integer argument.
- In the latter case, \code{Np(k)} must be a single positive integer, representing the number of particles to be used at the \code{k}-th timestep:
- \code{Np(0)} is the number of particles to use going from \code{timezero(object)} to \code{time(object)[1]},
- \code{Np(1)}, from \code{timezero(object)} to \code{time(object)[1]},
- and so on, while when \code{T=length(time(object,t0=TRUE))},
- \code{Np(T)} is the number of particles to sample at the end of the time-series.
- }
- \item{perturb.fn}{
- \code{perturb.fn(params,mifiter,timeno,\dots)}
- }
- \item{tol}{
- See the description under \code{\link{pfilter}}.
- }
- \item{max.fail}{
- See the description under \code{\link{pfilter}}.
- }
- \item{verbose}{
- logical; if TRUE, print progress reports.
- }
- \item{transform}{
- logical;
- if \code{TRUE}, optimization is performed on the transformed scale.
- }
- \item{\dots}{
- additional arguments that override the defaults.
- }
-}
-\references{
- E. L. Ionides, A. Bhadra, Y. Atchad{\\'e}, & A. A. King,
- Iterated filtering,
- Annals of Statistics, 39:1776--1802, 2011.
-}
-\author{Aaron A. King \email{kingaa at umich dot edu}}
-\seealso{
- \code{\link{pomp}}, \code{\link{pomp-class}}, \code{\link{pfilter}}.
-}
-\keyword{ts}
Deleted: pkg/mif2/tests/ou2-mif2.R
===================================================================
--- pkg/mif2/tests/ou2-mif2.R 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/tests/ou2-mif2.R 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,144 +0,0 @@
-library(mif2)
-
-pompExample(ou2)
-
-pdf(file="ou2-mif2.pdf")
-
-p.truth <- coef(ou2)
-guess2 <- guess1 <- p.truth
-guess1[c('x1.0','x2.0','alpha.2','alpha.3')] <- 0.9*guess1[c('x1.0','x2.0','alpha.2','alpha.3')]
-guess2[c('x1.0','x2.0','alpha.2','alpha.3')] <- 1.2*guess1[c('x1.0','x2.0','alpha.2','alpha.3')]
-
-set.seed(64857673L)
-mif1a <- mif(
- ou2,
- Nmif=100,
- start=guess1,
- pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
- rw.sd=c(
- x1.0=0.5,x2.0=0.5,
- alpha.2=0.1,alpha.3=0.1),
- transform=F,
- Np=1000,
- var.factor=1,
- ic.lag=10,
- cooling.type="hyperbolic",
- cooling.fraction=0.05,
- method="mif2",
- tol=1e-8
- )
-
-mif2a <- mif(ou2,Nmif=100,start=guess1,
- pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
- rw.sd=c(
- x1.0=0.5,x2.0=0.5,
- alpha.2=0.1,alpha.3=0.1),
- transform=F,
- Np=1000,
- var.factor=1,
- ic.lag=10,
- cooling.type="geometric",
- cooling.fraction=0.95^50,
- max.fail=100,
- method="mif",
- tol=1e-8
- )
-
-plot(c(mif1a,mif2a))
-
-set.seed(64857673L)
-mif1b <- mif(
- ou2,
- Nmif=50,
- start=guess1,
- pars=c('alpha.2','alpha.3'),
- ivps=c('x1.0','x2.0'),
- rw.sd=c(
- x1.0=0.5,x2.0=0.5,
- alpha.2=0.1,alpha.3=0.1),
- transform=F,
- Np=1000,
- var.factor=1,
- ic.lag=10,
- cooling.type="hyperbolic",
- cooling.fraction=0.05,
- method="mif2"
- )
-mif1b <- continue(mif1b,Nmif=50)
-
-mif2b <- mif(
- ou2,
- Nmif=50,
- start=guess1,
- pars=c('alpha.2','alpha.3'),
- ivps=c('x1.0','x2.0'),
- rw.sd=c(
- x1.0=0.5,x2.0=0.5,
- alpha.2=0.1,alpha.3=0.1),
- transform=F,
- Np=1000,
- var.factor=1,
- ic.lag=10,
- cooling.whatsit=200,
- cooling.type="geometric",
- cooling.factor=0.95,
- max.fail=100,
- method="mif"
- )
-mif2b <- continue(mif2b,Nmif=50)
-
-mif2c <- mif(
- ou2,
- Nmif=50,
- start=guess1,
- pars=c('alpha.2','alpha.3'),
- ivps=c('x1.0','x2.0'),
- rw.sd=c(
- x1.0=0.5,x2.0=0.5,
- alpha.2=0.1,alpha.3=0.1),
- transform=F,
- Np=1000,
- var.factor=1,
- cooling.type="hyperbolic",
- cooling.fraction=0.05,
- max.fail=100,
- method="mif2"
- )
-mif2c <- continue(mif2c,Nmif=50)
-
-plot(c(mif1b,mif2b))
-plot(c(mif1a,mif1b))
-plot(c(mif2a,mif2b))
-plot(c(mif1b,mif2c))
-
-mif3a <- mif2(
- ou2,
- Nmif=50,
- start=guess1,
- perturb.fn=function(params,mifiter,timeno,...){
- pert <- params
- ic.sd <- c(x1.0=0.5,x2.0=0.5)
- par.sd <- c(alpha.2=0.1,alpha.3=0.1)
- frac <- 0.05
- nT <- length(time(ou2))
- theta <- (1-frac)/frac/(50*nT-1)
- sigma <- 1/(1+theta*((mifiter-1)*nT+timeno-1))
- if (timeno==1) {
- pert[names(ic.sd)] <- rnorm(
- n=length(ic.sd),
- mean=pert[names(ic.sd)],
- sd=ic.sd*sigma
- )
- }
- pert[names(par.sd)] <- rnorm(
- n=length(par.sd),
- mean=pert[names(par.sd)],
- sd=par.sd*sigma
- )
- pert
- },
- transform=FALSE,
- Np=1000
- )
-
-dev.off()
Deleted: pkg/mif2/tests/ou2-mif2.Rout.save
===================================================================
--- pkg/mif2/tests/ou2-mif2.Rout.save 2015-06-12 12:27:56 UTC (rev 1203)
+++ pkg/mif2/tests/ou2-mif2.Rout.save 2015-06-15 13:13:10 UTC (rev 1204)
@@ -1,179 +0,0 @@
-
-R version 3.1.1 (2014-07-10) -- "Sock it to Me"
-Copyright (C) 2014 The R Foundation for Statistical Computing
-Platform: x86_64-unknown-linux-gnu (64-bit)
-
-R is free software and comes with ABSOLUTELY NO WARRANTY.
-You are welcome to redistribute it under certain conditions.
-Type 'license()' or 'licence()' for distribution details.
-
-R is a collaborative project with many contributors.
-Type 'contributors()' for more information and
-'citation()' on how to cite R or R packages in publications.
-
-Type 'demo()' for some demos, 'help()' for on-line help, or
-'help.start()' for an HTML browser interface to help.
-Type 'q()' to quit R.
-
-> library(mif2)
-Loading required package: pomp
-Loading required package: mvtnorm
-Loading required package: subplex
-Loading required package: nloptr
-Loading required package: deSolve
-Loading required package: coda
-Loading required package: lattice
->
-> pompExample(ou2)
-newly created pomp object(s):
- ou2
->
-> pdf(file="ou2-mif2.pdf")
->
-> p.truth <- coef(ou2)
-> guess2 <- guess1 <- p.truth
-> guess1[c('x1.0','x2.0','alpha.2','alpha.3')] <- 0.9*guess1[c('x1.0','x2.0','alpha.2','alpha.3')]
-> guess2[c('x1.0','x2.0','alpha.2','alpha.3')] <- 1.2*guess1[c('x1.0','x2.0','alpha.2','alpha.3')]
->
-> set.seed(64857673L)
-> mif1a <- mif(
-+ ou2,
-+ Nmif=100,
-+ start=guess1,
-+ pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
-+ rw.sd=c(
-+ x1.0=0.5,x2.0=0.5,
-+ alpha.2=0.1,alpha.3=0.1),
-+ transform=F,
-+ Np=1000,
-+ var.factor=1,
-+ ic.lag=10,
-+ cooling.type="hyperbolic",
-+ cooling.fraction=0.05,
-+ method="mif2",
-+ tol=1e-8
-+ )
->
-> mif2a <- mif(ou2,Nmif=100,start=guess1,
-+ pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
-+ rw.sd=c(
-+ x1.0=0.5,x2.0=0.5,
-+ alpha.2=0.1,alpha.3=0.1),
-+ transform=F,
-+ Np=1000,
-+ var.factor=1,
-+ ic.lag=10,
-+ cooling.type="geometric",
-+ cooling.fraction=0.95^50,
-+ max.fail=100,
-+ method="mif",
-+ tol=1e-8
-+ )
->
-> plot(c(mif1a,mif2a))
->
-> set.seed(64857673L)
-> mif1b <- mif(
-+ ou2,
-+ Nmif=50,
-+ start=guess1,
-+ pars=c('alpha.2','alpha.3'),
-+ ivps=c('x1.0','x2.0'),
-+ rw.sd=c(
-+ x1.0=0.5,x2.0=0.5,
-+ alpha.2=0.1,alpha.3=0.1),
-+ transform=F,
-+ Np=1000,
-+ var.factor=1,
-+ ic.lag=10,
-+ cooling.type="hyperbolic",
-+ cooling.fraction=0.05,
-+ method="mif2"
-+ )
-> mif1b <- continue(mif1b,Nmif=50)
->
-> mif2b <- mif(
-+ ou2,
-+ Nmif=50,
-+ start=guess1,
-+ pars=c('alpha.2','alpha.3'),
-+ ivps=c('x1.0','x2.0'),
-+ rw.sd=c(
-+ x1.0=0.5,x2.0=0.5,
-+ alpha.2=0.1,alpha.3=0.1),
-+ transform=F,
-+ Np=1000,
-+ var.factor=1,
-+ ic.lag=10,
-+ cooling.whatsit=200,
-+ cooling.type="geometric",
-+ cooling.factor=0.95,
-+ max.fail=100,
-+ method="mif"
-+ )
-Warning message:
-'cooling.factor' is deprecated.
-See '?mif' for instructions on specifying the cooling schedule.
-> mif2b <- continue(mif2b,Nmif=50)
->
-> mif2c <- mif(
-+ ou2,
-+ Nmif=50,
-+ start=guess1,
-+ pars=c('alpha.2','alpha.3'),
-+ ivps=c('x1.0','x2.0'),
-+ rw.sd=c(
-+ x1.0=0.5,x2.0=0.5,
-+ alpha.2=0.1,alpha.3=0.1),
-+ transform=F,
-+ Np=1000,
-+ var.factor=1,
-+ cooling.type="hyperbolic",
-+ cooling.fraction=0.05,
-+ max.fail=100,
-+ method="mif2"
-+ )
-> mif2c <- continue(mif2c,Nmif=50)
->
-> plot(c(mif1b,mif2b))
-> plot(c(mif1a,mif1b))
-> plot(c(mif2a,mif2b))
-> plot(c(mif1b,mif2c))
->
-> mif3a <- mif2(
-+ ou2,
-+ Nmif=50,
-+ start=guess1,
-+ perturb.fn=function(params,mifiter,timeno,...){
-+ pert <- params
-+ ic.sd <- c(x1.0=0.5,x2.0=0.5)
-+ par.sd <- c(alpha.2=0.1,alpha.3=0.1)
-+ frac <- 0.05
-+ nT <- length(time(ou2))
-+ theta <- (1-frac)/frac/(50*nT-1)
-+ sigma <- 1/(1+theta*((mifiter-1)*nT+timeno-1))
-+ if (timeno==1) {
-+ pert[names(ic.sd)] <- rnorm(
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/pomp -r 1204
More information about the pomp-commits
mailing list