[Pomp-commits] r864 - in branches/mif2: . R inst inst/include man src tests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Nov 5 15:05:35 CET 2013
Author: kingaa
Date: 2013-11-05 15:05:32 +0100 (Tue, 05 Nov 2013)
New Revision: 864
Added:
branches/mif2/src/SSA.f90
Removed:
branches/mif2/inst/LICENSE
branches/mif2/src/SSA.f
Modified:
branches/mif2/.Rbuildignore
branches/mif2/DESCRIPTION
branches/mif2/NAMESPACE
branches/mif2/R/builder.R
branches/mif2/R/mif.R
branches/mif2/R/pfilter-methods.R
branches/mif2/R/pfilter.R
branches/mif2/R/profile-design.R
branches/mif2/inst/NEWS
branches/mif2/inst/include/pomp.h
branches/mif2/man/pfilter-methods.Rd
branches/mif2/man/profile-design.Rd
branches/mif2/src/SSA_wrapper.c
branches/mif2/src/eulermultinom.c
branches/mif2/src/lookup_table.c
branches/mif2/src/pomp.h
branches/mif2/tests/ou2-mif.Rout.save
branches/mif2/tests/pfilter.R
branches/mif2/tests/pfilter.Rout.save
Log:
- bring mif2 branch up to date with changes to trunk
Modified: branches/mif2/.Rbuildignore
===================================================================
--- branches/mif2/.Rbuildignore 2013-06-11 16:28:47 UTC (rev 863)
+++ branches/mif2/.Rbuildignore 2013-11-05 14:05:32 UTC (rev 864)
@@ -2,5 +2,3 @@
inst/doc/(.+?)\.bst$
inst/doc/(.+?)\.R$
inst/doc/(.+?)\.png$
-^.*\.Rproj$
-^\.Rproj\.user$
Modified: branches/mif2/DESCRIPTION
===================================================================
--- branches/mif2/DESCRIPTION 2013-06-11 16:28:47 UTC (rev 863)
+++ branches/mif2/DESCRIPTION 2013-11-05 14:05:32 UTC (rev 864)
@@ -1,27 +1,28 @@
Package: pomp
Type: Package
Title: Statistical inference for partially observed Markov processes
-Version: 0.44-1
-Date: 2013-03-26
-Author: Aaron A. King, Edward L. Ionides, Carles Breto, Steve Ellner, Bruce Kendall, Dao Nguyen, Helen Wearing, Matthew J. Ferrari, Michael Lavine, Daniel C. Reuman
-Maintainer: Aaron A. King <kingaa at umich.edu>
-Author at R: c(person(given=c("Aaron","A."),family="King",role=c("aut","cre"),email="kingaa at umich.edu"),
+Version: 0.46-1
+Date: 2013-10-30
+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=c("Carles"),family="Breto",role=c("aut")),
person(given=c("Stephen","P."),family="Ellner",role=c("ctb")),
person(given=c("Matthew","J."),family="Ferrari",role=c("ctb")),
person(given=c("Bruce","E."),family="Kendall",role=c("ctb")),
person(given=c("Michael"),family="Lavine",role=c("ctb")),
- person(given=c("Dao"),family="Nguyen",role=c("ctb")),
person(given=c("Daniel","C."),family="Reuman",role=c("ctb")),
person(given=c("Helen"),family="Wearing",role=c("ctb")),
- person(given=c("Simon","N."),family="Wood",role=c("ctb")))
+ person(given=c("Simon","N."),family="Wood",role=c("ctb")),
+ person(given="Dao",family="Nguyen",role=c("ctb")) )
+Maintainer: Aaron A. King <kingaa at umich.edu>
+Author: Aaron A. King, Edward L. Ionides, Carles Breto, Steve Ellner, Bruce Kendall, Helen Wearing, Matthew J. Ferrari, Michael Lavine, Daniel C. Reuman, Simon N. Wood, Dao Nguyen
URL: http://pomp.r-forge.r-project.org
Description: Inference methods for partially-observed Markov processes
-Depends: R(>= 2.14.1), stats, methods, graphics, mvtnorm, subplex, deSolve
+Depends: R(>= 2.14.1), stats, graphics, mvtnorm, subplex, deSolve
+Imports: methods
License: GPL(>= 2)
-LazyLoad: true
-LazyData: false
+LazyLoad: yes
+LazyData: no
BuildVignettes: no
Collate: aaa.R authors.R version.R eulermultinom.R plugins.R
parmat.R slice-design.R profile-design.R sobol.R bsplines.R sannbox.R
Modified: branches/mif2/NAMESPACE
===================================================================
--- branches/mif2/NAMESPACE 2013-06-11 16:28:47 UTC (rev 863)
+++ branches/mif2/NAMESPACE 2013-11-05 14:05:32 UTC (rev 864)
@@ -62,6 +62,7 @@
export(
as.data.frame.pomp,
+ as.data.frame.pfilterd.pomp,
reulermultinom,
deulermultinom,
rgammawn,
Modified: branches/mif2/R/builder.R
===================================================================
--- branches/mif2/R/builder.R 2013-06-11 16:28:47 UTC (rev 863)
+++ branches/mif2/R/builder.R 2013-11-05 14:05:32 UTC (rev 864)
@@ -107,7 +107,6 @@
decl <- list(
periodic_bspline_basis_eval="\tvoid (*periodic_bspline_basis_eval)(double,double,int,int,double*);\nperiodic_bspline_basis_eval = (void (*)(double,double,int,int,double*)) R_GetCCallable(\"pomp\",\"periodic_bspline_basis_eval\");\n",
- reulermultinom="\tvoid (*reulermultinom)(int,double,double*,double,double*);\nreulermultinom = (void (*)(int,double,double*,double,double*)) R_GetCCallable(\"pomp\",\"reulermultinom\");\n",
get_pomp_userdata="\tconst SEXP (*get_pomp_userdata)(const char *);\npomp_get_userdata = (const SEXP (*)(const char*)) R_GetCCallable(\"pomp\",\"get_pomp_userdata\");\n",
get_pomp_userdata_int="\tconst int * (*get_pomp_userdata_int)(const char *);\npomp_get_userdata_int = (const int *(*)(const char*)) R_GetCCallable(\"pomp\",\"get_pomp_userdata_int\");\n",
get_pomp_userdata_double="\tconst double * (*get_pomp_userdata_double)(const char *);\npomp_get_userdata_double = (const double *(*)(const char*)) R_GetCCallable(\"pomp\",\"get_pomp_userdata_double\");\n"
Modified: branches/mif2/R/mif.R
===================================================================
--- branches/mif2/R/mif.R 2013-06-11 16:28:47 UTC (rev 863)
+++ branches/mif2/R/mif.R 2013-11-05 14:05:32 UTC (rev 864)
@@ -2,55 +2,55 @@
default.pomp.particles.fun <- function (Np, center, sd, ...) {
matrix(
- data=rnorm(
- n=Np*length(center),
- mean=center,
- sd=sd
- ),
- nrow=length(center),
- ncol=Np,
- dimnames=list(
- names(center),
- NULL
- )
- )
+ data=rnorm(
+ n=Np*length(center),
+ mean=center,
+ sd=sd
+ ),
+ nrow=length(center),
+ ncol=Np,
+ dimnames=list(
+ names(center),
+ NULL
+ )
+ )
}
cooling.function <- function (type, perobs, fraction, ntimes) {
switch(
- type,
- geometric={
- factor <- fraction^(1/50)
- if (perobs) {
- function (nt, m) {
- alpha <- factor^(nt/ntimes+m-1)
- list(alpha=alpha,gamma=alpha^2)
- }
- } else {
- function (nt, m) {
- alpha <- factor^(m-1)
- list(alpha=alpha,gamma=alpha^2)
- }
- }
- },
- hyperbolic={
- if (perobs) {
- scal <- (50*ntimes*fraction-1)/(1-fraction)
- function (nt, m) {
- alpha <- (1+scal)/(scal+nt+ntimes*(m-1))
- list(alpha=alpha,gamma=alpha^2)
- }
- } else {
- scal <- (50*fraction-1)/(1-fraction)
- function (nt, m) {
- alpha <- (1+scal)/(scal+m-1)
- list(alpha=alpha,gamma=alpha^2)
- }
-
- }
- },
- stop("unrecognized cooling schedule type ",sQuote(type))
- )
+ type,
+ geometric={
+ factor <- fraction^(1/50)
+ if (perobs) {
+ function (nt, m) {
+ alpha <- factor^(nt/ntimes+m-1)
+ list(alpha=alpha,gamma=alpha^2)
+ }
+ } else {
+ function (nt, m) {
+ alpha <- factor^(m-1)
+ list(alpha=alpha,gamma=alpha^2)
+ }
+ }
+ },
+ hyperbolic={
+ if (perobs) {
+ scal <- (50*ntimes*fraction-1)/(1-fraction)
+ function (nt, m) {
+ alpha <- (1+scal)/(scal+nt+ntimes*(m-1))
+ list(alpha=alpha,gamma=alpha^2)
+ }
+ } else {
+ scal <- (50*fraction-1)/(1-fraction)
+ function (nt, m) {
+ alpha <- (1+scal)/(scal+m-1)
+ list(alpha=alpha,gamma=alpha^2)
+ }
+
+ }
+ },
+ stop("unrecognized cooling schedule type ",sQuote(type))
+ )
}
mif.cooling <- function (factor, n) { # default geometric cooling schedule
@@ -95,10 +95,10 @@
if (length(start)==0)
stop(
- "mif error: ",sQuote("start")," must be specified if ",
- sQuote("coef(object)")," is NULL",
- call.=FALSE
- )
+ "mif error: ",sQuote("start")," must be specified if ",
+ sQuote("coef(object)")," is NULL",
+ call.=FALSE
+ )
if (transform)
start <- partrans(object,start,dir="inverse")
@@ -112,7 +112,7 @@
stop("mif error: ",sQuote("rw.sd")," must be a named non-negative numerical vector",call.=FALSE)
if (!all(rw.names%in%start.names))
stop("mif error: all the names of ",sQuote("rw.sd")," must be names of ",sQuote("start"),call.=FALSE)
- #rw.names <- names(rw.sd[rw.sd>0])
+ #rw.names <- names(rw.sd[rw.sd>0])
if (length(rw.names) == 0)
stop("mif error: ",sQuote("rw.sd")," must have one positive entry for each parameter to be estimated",call.=FALSE)
@@ -120,7 +120,7 @@
if (missing(ivps)) stop("mif error: ",sQuote("ivps")," must be specified",call.=FALSE)
if (
- !is.character(pars) ||
+ !is.character(pars) ||
!is.character(ivps) ||
!all(pars%in%start.names) ||
!all(ivps%in%start.names) ||
@@ -128,32 +128,32 @@
any(ivps%in%pars) ||
!all(pars%in%rw.names) ||
!all(ivps%in%rw.names)
- )
+ )
stop(
- "mif error: ",
- sQuote("pars")," and ",sQuote("ivps"),
- " must be mutually disjoint subsets of ",
- sQuote("names(start)"),
- " and must have a positive random-walk SDs specified in ",
- sQuote("rw.sd"),
- call.=FALSE
- )
+ "mif error: ",
+ sQuote("pars")," and ",sQuote("ivps"),
+ " must be mutually disjoint subsets of ",
+ sQuote("names(start)"),
+ " and must have a positive random-walk SDs specified in ",
+ sQuote("rw.sd"),
+ call.=FALSE
+ )
if (!all(rw.names%in%c(pars,ivps))) {
extra.rws <- rw.names[!(rw.names%in%c(pars,ivps))]
warning(
- ngettext(length(extra.rws),"mif warning: the variable ",
- "mif warning: the variables "),
- paste(sQuote(extra.rws),collapse=", "),
- ngettext(length(extra.rws)," has positive random-walk SD specified, but is included in neither ",
- " have positive random-walk SDs specified, but are included in neither "),
- sQuote("pars")," nor ",sQuote("ivps"),
- ngettext(length(extra.rws),". This random walk SD will be ignored.",
- ". These random walk SDs will be ignored."),
- call.=FALSE
- )
+ ngettext(length(extra.rws),"mif warning: the variable ",
+ "mif warning: the variables "),
+ paste(sQuote(extra.rws),collapse=", "),
+ ngettext(length(extra.rws)," has positive random-walk SD specified, but is included in neither ",
+ " have positive random-walk SDs specified, but are included in neither "),
+ sQuote("pars")," nor ",sQuote("ivps"),
+ ngettext(length(extra.rws),". This random walk SD will be ignored.",
+ ". These random walk SDs will be ignored."),
+ call.=FALSE
+ )
}
- #rw.sd <- rw.sd[c(pars,ivps)]
+ #rw.sd <- rw.sd[c(pars,ivps)]
rw.names <- colnames(rw.sd)
rwsdMat <-rw.sd
if (missing(particles))
@@ -163,9 +163,9 @@
if (missing(Np)) stop("mif error: ",sQuote("Np")," must be specified",call.=FALSE)
if (is.function(Np)) {
Np <- try(
- vapply(seq.int(from=0,to=ntimes,by=1),Np,numeric(1)),
- silent=FALSE
- )
+ vapply(seq.int(from=0,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")
}
@@ -185,20 +185,20 @@
stop("mif error: ",sQuote("ic.lag")," must be a positive integer",call.=FALSE)
if (ic.lag>ntimes) {
warning(
- "mif warning: ",sQuote("ic.lag")," = ",ic.lag," > ",ntimes,
- " = length(time(",sQuote("object"),"))",
- " is nonsensical. Setting ",sQuote("ic.lag")," = ",ntimes,".",
- call.=FALSE
- )
+ "mif warning: ",sQuote("ic.lag")," = ",ic.lag," > ",ntimes,
+ " = length(time(",sQuote("object"),"))",
+ " is nonsensical. Setting ",sQuote("ic.lag")," = ",ntimes,".",
+ call.=FALSE
+ )
ic.lag <- length(time(object))
}
if ((length(pars)==0)&&(ic.lag<length(time(object)))) {
warning(
- "mif warning: only IVPs are to be estimated, yet ",sQuote("ic.lag")," = ",ic.lag,
- " < ",ntimes," = length(time(",sQuote("object"),")),",
- " so unnecessary work is to be done.",
- call.=FALSE
- )
+ "mif warning: only IVPs are to be estimated, yet ",sQuote("ic.lag")," = ",ic.lag,
+ " < ",ntimes," = length(time(",sQuote("object"),")),",
+ " so unnecessary work is to be done.",
+ call.=FALSE
+ )
}
## the following deals with the deprecated option 'cooling.factor'
@@ -225,11 +225,11 @@
stop("mif error: ",sQuote("cooling.fraction")," must be a number between 0 and 1",call.=FALSE)
cooling <- cooling.function(
- type=cooling.type,
- perobs=(method=="mif2"),
- fraction=cooling.fraction,
- ntimes=ntimes
- )
+ type=cooling.type,
+ perobs=(method=="mif2"),
+ fraction=cooling.fraction,
+ ntimes=ntimes
+ )
if ((method=="mif2")&&(Np[1L]!=Np[ntimes+1]))
stop("the first and last values of ",sQuote("Np")," must agree when method = ",sQuote("mif2"))
@@ -245,26 +245,26 @@
conv.rec <- matrix(
- data=NA,
- nrow=Nmif+1,
- ncol=length(theta)+2,
- dimnames=list(
- seq(.ndone,.ndone+Nmif),
- c('loglik','nfail',names(theta))
- )
- )
+ data=NA,
+ nrow=Nmif+1,
+ ncol=length(theta)+2,
+ dimnames=list(
+ seq(.ndone,.ndone+Nmif),
+ c('loglik','nfail',names(theta))
+ )
+ )
conv.rec[1L,] <- c(NA,NA,theta)
if (!all(is.finite(theta[c(pars,ivps)]))) {
stop(
- sQuote("mif")," cannot estimate non-finite parameters.\n",
- "The following ",if (transform) "transformed ", "parameters are non-finite: ",
- paste(
- sQuote(c(pars,ivps)[!is.finite(theta[c(pars,ivps)])]),
- collapse=","
- ),
- call.=FALSE
- )
+ sQuote("mif")," cannot estimate non-finite parameters.\n",
+ "The following ",if (transform) "transformed ", "parameters are non-finite: ",
+ paste(
+ sQuote(c(pars,ivps)[!is.finite(theta[c(pars,ivps)])]),
+ collapse=","
+ ),
+ call.=FALSE
+ )
}
obj <- as(object,"pomp")
@@ -286,14 +286,14 @@
## initialize the parameter portions of the particles
P <- try(
- particles(
- tmp.mif,
- Np=Np[1L],
- center=theta,
- sd=sigma.n*var.factor
- ),
- silent = FALSE
- )
+ particles(
+ tmp.mif,
+ Np=Np[1L],
+ center=theta,
+ sd=sigma.n*var.factor
+ ),
+ silent = FALSE
+ )
if (inherits(P,"try-error"))
stop("mif error: error in ",sQuote("particles"),call.=FALSE)
@@ -303,27 +303,27 @@
}
colnames(rwsdMat) <- names(start)
pfp <- try(
- pfilter.internal(
- object=obj,
- params=P,
- Np=Np,
- tol=tol,
- max.fail=max.fail,
- pred.mean=(n==Nmif),
- pred.var=((method=="mif")||(n==Nmif)),
- filter.mean=TRUE,
- cooling=cooling,
- cooling.m=.ndone+n,
- .mif2=(method=="mif2"),
- .rw.sd=rwsdMat*cool.sched$alpha,
- .transform=transform,
- save.states=FALSE,
- save.params=FALSE,
- verbose=verbose,
- .getnativesymbolinfo=gnsi
- ),
- silent=FALSE
- )
+ pfilter.internal(
+ object=obj,
+ params=P,
+ Np=Np,
+ tol=tol,
+ max.fail=max.fail,
+ pred.mean=(n==Nmif),
+ pred.var=((method=="mif")||(n==Nmif)),
+ filter.mean=TRUE,
+ cooling=cooling,
+ cooling.m=.ndone+n,
+ .mif2=(method=="mif2"),
+ .rw.sd=rwsdMat*cool.sched$alpha,
+ .transform=transform,
+ save.states=FALSE,
+ save.params=FALSE,
+ verbose=verbose,
+ .getnativesymbolinfo=gnsi
+ ),
+ silent=FALSE
+ )
if (inherits(pfp,"try-error"))
stop("mif error: error in ",sQuote("pfilter"),call.=FALSE)
@@ -331,26 +331,26 @@
## update parameters
switch(
- method,
- mif={ # original Ionides et al. (2006) average
- theta <- .Call(mif_update,pfp,theta,cool.sched$gamma,var.factor,rwsdMat[1,],pars)
- theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
- },
- unweighted={ # unweighted average
- theta[pars] <- rowMeans(pfp at filter.mean[pars,,drop=FALSE])
- theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
- },
- fp={ # fixed-point iteration
- theta[pars] <- pfp at filter.mean[pars,ntimes,drop=FALSE]
- theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
- },
- mif2={ # "efficient" iterated filtering
- paramMatrix <- pfp at paramMatrix
- theta[pars] <- rowMeans(paramMatrix[pars,,drop=FALSE])
- theta[ivps] <- pfp at filter.mean[ivps,ntimes]
- },
- stop("unrecognized method ",sQuote(method))
- )
+ method,
+ mif={ # original Ionides et al. (2006) average
+ theta <- .Call(mif_update,pfp,theta,cool.sched$gamma,var.factor,rwsdMat[1,],pars)
+ theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
+ },
+ unweighted={ # unweighted average
+ theta[pars] <- rowMeans(pfp at filter.mean[pars,,drop=FALSE])
+ theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
+ },
+ fp={ # fixed-point iteration
+ theta[pars] <- pfp at filter.mean[pars,ntimes,drop=FALSE]
+ theta[ivps] <- pfp at filter.mean[ivps,ic.lag]
+ },
+ mif2={ # "efficient" iterated filtering
+ paramMatrix <- pfp at paramMatrix
+ theta[pars] <- rowMeans(paramMatrix[pars,,drop=FALSE])
+ theta[ivps] <- pfp at filter.mean[ivps,ntimes]
+ },
+ stop("unrecognized method ",sQuote(method))
+ )
conv.rec[n+1,-c(1,2)] <- theta
conv.rec[n,c(1,2)] <- c(pfp at loglik,pfp at nfail)
@@ -363,266 +363,266 @@
if (transform) theta <- partrans(pfp,theta,dir="forward")
new(
- "mif",
- pfp,
- transform=transform,
- params=theta,
- ivps=ivps,
- pars=pars,
- Nmif=Nmif,
- particles=particles,
- var.factor=var.factor,
- ic.lag=ic.lag,
- random.walk.sd=rwsdMat,
- tol=tol,
- conv.rec=conv.rec,
- method=method,
- cooling.type=cooling.type,
- cooling.fraction=cooling.fraction,
- paramMatrix=if (method=="mif2") paramMatrix else array(data=numeric(0),dim=c(0,0))
- )
+ "mif",
+ pfp,
+ transform=transform,
+ params=theta,
+ ivps=ivps,
+ pars=pars,
+ Nmif=Nmif,
+ particles=particles,
+ var.factor=var.factor,
+ ic.lag=ic.lag,
+ random.walk.sd=rwsdMat,
+ tol=tol,
+ conv.rec=conv.rec,
+ method=method,
+ cooling.type=cooling.type,
+ cooling.fraction=cooling.fraction,
+ paramMatrix=if (method=="mif2") paramMatrix else array(data=numeric(0),dim=c(0,0))
+ )
}
setGeneric('mif',function(object,...)standardGeneric("mif"))
setMethod(
- "mif",
- signature=signature(object="pomp"),
- function (object, Nmif = 1,
- start,
- pars, ivps = character(0),
- particles, rw.sd,
- Np, ic.lag, var.factor,
- cooling.type = c("geometric","hyperbolic"),
- cooling.fraction, cooling.factor,
- method = c("mif","unweighted","fp","mif2"),
- tol = 1e-17, max.fail = Inf,
- verbose = getOption("verbose"),
- transform = FALSE,
- ...) {
-
- transform <- as.logical(transform)
- method <- match.arg(method)
-
-
- ntimes <- length(time(object))
- if (missing(start)) start <- coef(object)
- if (missing(rw.sd))
- stop("mif error: ",sQuote("rw.sd")," must be specified",call.=FALSE)
- if (missing(pars)) {
- stop("mif error: ",sQuote("par")," must be specified",call.=FALSE)
- }
- rw.names <- c(pars,ivps)
- dtheta <- length(start)
- rwsdMat <- matrix(rep(0, dtheta*(ntimes+1)), ncol=dtheta, nrow = (ntimes+1))
- colnames(rwsdMat) <- names(start)
-
- if (is.matrix(rw.sd)) rwsdMat<-rw.sd
- else if (is.list(rw.sd))
- {
- for (i in 1:length(rw.names))
- { if (rw.names[i] %in% ivps)
- {
- if (length((rw.sd[[rw.names[i]]])==1))
- {
- rwsdMat[1,rw.names[i]] <- rw.sd[[rw.names[i]]][1]
- }
- else
- {
- stop(sQuote("rw.sd")," must have length 1 or length ",ntimes+1)
-
- }
- }
- else if (rw.names[i] %in% pars)
- {
- if (length(rw.sd[[rw.names[i]]])==1)
- {
- rwsdMat[,rw.names[i]] <- rep(rw.sd[[rw.names[i]]],(ntimes+1))
- }
- else if (length(rw.sd[[rw.names[i]]])==(ntimes+1))
- {
- rwsdMat[,rw.names[i]] <- rw.sd[[rw.names[i]]]
- }
- else
- {
- stop(sQuote("rw.sd")," must have length 1 or length ",ntimes+1)
+ "mif",
+ signature=signature(object="pomp"),
+ function (object, Nmif = 1,
+ start,
+ pars, ivps = character(0),
+ particles, rw.sd,
+ Np, ic.lag, var.factor,
+ cooling.type = c("geometric","hyperbolic"),
+ cooling.fraction, cooling.factor,
+ method = c("mif","unweighted","fp","mif2"),
+ tol = 1e-17, max.fail = Inf,
+ verbose = getOption("verbose"),
+ transform = FALSE,
+ ...) {
+ transform <- as.logical(transform)
+ method <- match.arg(method)
+
+
+ ntimes <- length(time(object))
+ if (missing(start)) start <- coef(object)
+ if (missing(rw.sd))
+ stop("mif error: ",sQuote("rw.sd")," must be specified",call.=FALSE)
+ if (missing(pars)) {
+ stop("mif error: ",sQuote("par")," must be specified",call.=FALSE)
+ }
+ rw.names <- c(pars,ivps)
+ dtheta <- length(start)
+ rwsdMat <- matrix(rep(0, dtheta*(ntimes+1)), ncol=dtheta, nrow = (ntimes+1))
+ colnames(rwsdMat) <- names(start)
+
+ if (is.matrix(rw.sd)) rwsdMat<-rw.sd
+ else if (is.list(rw.sd))
+ {
+ for (i in 1:length(rw.names))
+ { if (rw.names[i] %in% ivps)
+ {
+ if (length((rw.sd[[rw.names[i]]])==1))
+ {
+ rwsdMat[1,rw.names[i]] <- rw.sd[[rw.names[i]]][1]
+ }
+ else
+ {
+ stop(sQuote("rw.sd")," must have length 1 or length ",ntimes+1)
+
+ }
+ }
+ else if (rw.names[i] %in% pars)
+ {
+ if (length(rw.sd[[rw.names[i]]])==1)
+ {
+ rwsdMat[,rw.names[i]] <- rep(rw.sd[[rw.names[i]]],(ntimes+1))
+ }
+ else if (length(rw.sd[[rw.names[i]]])==(ntimes+1))
+ {
+ rwsdMat[,rw.names[i]] <- rw.sd[[rw.names[i]]]
+ }
+ else
+ {
+ stop(sQuote("rw.sd")," must have length 1 or length ",ntimes+1)
+
+ }
+ }
+
+ }
+
+ }
+ else if (is.vector(rw.sd))
+ {
+ if (missing(pars)) {
+ rw.names <- names(rw.sd)[rw.sd>0]
+ pars <- rw.names[!(rw.names%in%ivps)]
+ }
+ for (i in 1:length(rw.names))
+ { if (rw.names[i] %in% ivps)
+ {
+ rwsdMat[1,rw.names[i]] <- rw.sd[rw.names[i]]
+
+ }
+ else if (rw.names[i] %in% pars)
+ {
+ rwsdMat[,rw.names[i]] <- rep(rw.sd[rw.names[i]],(ntimes+1))
+ }
+
+ }
+
+
+
+ }
+
+
+
+ if (missing(Np))
+ stop("mif error: ",sQuote("Np")," must be specified",call.=FALSE)
+ if (missing(ic.lag) && length(ivps)>0)
+ stop("mif error: ",sQuote("ic.lag")," must be specified if ",sQuote("ivps")," are",call.=FALSE)
+ if (missing(var.factor))
+ stop("mif error: ",sQuote("var.factor")," must be specified",call.=FALSE)
+
+ cooling.type <- match.arg(cooling.type)
+
+ if (missing(particles)) { # use default: normal distribution
+ particles <- default.pomp.particles.fun
+ } else {
+ particles <- match.fun(particles)
+ if (!all(c('Np','center','sd','...')%in%names(formals(particles))))
+ stop(
+ "mif error: ",
+ sQuote("particles"),
+ " must be a function of prototype ",
+ sQuote("particles(Np,center,sd,...)"),
+ call.=FALSE
+ )
+ }
+
+ mif.internal(
+ object=object,
+ Nmif=Nmif,
+ start=start,
+ pars=pars,
+ ivps=ivps,
+ particles=particles,
+ rw.sd=rwsdMat,
+ Np=Np,
+ cooling.type=cooling.type,
+ cooling.factor=cooling.factor,
+ cooling.fraction=cooling.fraction,
+ var.factor=var.factor,
+ ic.lag=ic.lag,
+ method=method,
+ tol=tol,
+ max.fail=max.fail,
+ verbose=verbose,
+ transform=transform
+ )
+
}
- }
-
- }
-
- }
- else if (is.vector(rw.sd))
- {
- if (missing(pars)) {
- rw.names <- names(rw.sd)[rw.sd>0]
- pars <- rw.names[!(rw.names%in%ivps)]
- }
- for (i in 1:length(rw.names))
- { if (rw.names[i] %in% ivps)
- {
- rwsdMat[1,rw.names[i]] <- rw.sd[rw.names[i]]
-
- }
- else if (rw.names[i] %in% pars)
- {
- rwsdMat[,rw.names[i]] <- rep(rw.sd[rw.names[i]],(ntimes+1))
- }
-
- }
-
-
-
- }
-
-
-
- if (missing(Np))
- stop("mif error: ",sQuote("Np")," must be specified",call.=FALSE)
- if (missing(ic.lag) && length(ivps)>0)
- stop("mif error: ",sQuote("ic.lag")," must be specified if ",sQuote("ivps")," are",call.=FALSE)
- if (missing(var.factor))
- stop("mif error: ",sQuote("var.factor")," must be specified",call.=FALSE)
-
- cooling.type <- match.arg(cooling.type)
-
- if (missing(particles)) { # use default: normal distribution
- particles <- default.pomp.particles.fun
- } else {
- particles <- match.fun(particles)
- if (!all(c('Np','center','sd','...')%in%names(formals(particles))))
- stop(
- "mif error: ",
- sQuote("particles"),
- " must be a function of prototype ",
- sQuote("particles(Np,center,sd,...)"),
- call.=FALSE
- )
- }
-
- mif.internal(
- object=object,
- Nmif=Nmif,
- start=start,
- pars=pars,
- ivps=ivps,
- particles=particles,
- rw.sd=rwsdMat,
- Np=Np,
- cooling.type=cooling.type,
- cooling.factor=cooling.factor,
- cooling.fraction=cooling.fraction,
- var.factor=var.factor,
- ic.lag=ic.lag,
- method=method,
- tol=tol,
- max.fail=max.fail,
- verbose=verbose,
- transform=transform
- )
-
- }
-)
+ )
setMethod(
- "mif",
- 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
-
- mif(
- object=as(object,"pomp"),
- Nmif=Nmif,
- Np=Np,
- tol=tol,
- ...
- )
- }
-)
+ "mif",
+ 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
+
+ mif(
+ object=as(object,"pomp"),
+ Nmif=Nmif,
+ Np=Np,
+ tol=tol,
+ ...
+ )
+ }
+ )
setMethod(
- "mif",
- signature=signature(object="mif"),
- function (object, Nmif,
- start,
- pars, ivps,
- particles, rw.sd,
- Np, ic.lag, var.factor,
- cooling.type, cooling.fraction,
- method,
- tol,
- transform,
- ...) {
-
- if (missing(Nmif)) Nmif <- object at Nmif
- if (missing(start)) start <- coef(object)
- if (missing(pars)) pars <- object at pars
- if (missing(ivps)) ivps <- object at ivps
- if (missing(particles)) particles <- object at particles
- if (missing(rw.sd)) rw.sd <- object at random.walk.sd
- if (missing(ic.lag)) ic.lag <- object at ic.lag
- if (missing(var.factor)) var.factor <- object at var.factor
- if (missing(cooling.type)) cooling.type <- object at cooling.type
- if (missing(cooling.fraction)) cooling.fraction <- object at cooling.fraction
- if (missing(method)) method <- object at method
- if (missing(transform)) transform <- object at transform
- transform <- as.logical(transform)
-
- if (missing(Np)) Np <- object at Np
- if (missing(tol)) tol <- object at tol
-
- mif(
- object=as(object,"pomp"),
- Nmif=Nmif,
- start=start,
- pars=pars,
- ivps=ivps,
- particles=particles,
- rw.sd=rw.sd,
- Np=Np,
- cooling.type=cooling.type,
- cooling.fraction=cooling.fraction,
- var.factor=var.factor,
- ic.lag=ic.lag,
- method=method,
- tol=tol,
- transform=transform,
- ...
- )
- }
-)
+ "mif",
+ signature=signature(object="mif"),
+ function (object, Nmif,
+ start,
+ pars, ivps,
+ particles, rw.sd,
+ Np, ic.lag, var.factor,
+ cooling.type, cooling.fraction,
+ method,
+ tol,
+ transform,
+ ...) {
+
+ if (missing(Nmif)) Nmif <- object at Nmif
+ if (missing(start)) start <- coef(object)
+ if (missing(pars)) pars <- object at pars
+ if (missing(ivps)) ivps <- object at ivps
+ if (missing(particles)) particles <- object at particles
+ if (missing(rw.sd)) rw.sd <- object at random.walk.sd
+ if (missing(ic.lag)) ic.lag <- object at ic.lag
+ if (missing(var.factor)) var.factor <- object at var.factor
+ if (missing(cooling.type)) cooling.type <- object at cooling.type
+ if (missing(cooling.fraction)) cooling.fraction <- object at cooling.fraction
+ if (missing(method)) method <- object at method
+ if (missing(transform)) transform <- object at transform
+ transform <- as.logical(transform)
+
+ if (missing(Np)) Np <- object at Np
+ if (missing(tol)) tol <- object at tol
+
+ mif(
+ object=as(object,"pomp"),
+ Nmif=Nmif,
+ start=start,
+ pars=pars,
+ ivps=ivps,
+ particles=particles,
+ rw.sd=rw.sd,
+ Np=Np,
+ cooling.type=cooling.type,
+ cooling.fraction=cooling.fraction,
+ var.factor=var.factor,
+ ic.lag=ic.lag,
+ method=method,
+ tol=tol,
+ transform=transform,
+ ...
+ )
+ }
+ )
setMethod(
- 'continue',
- signature=signature(object='mif'),
- function (object, Nmif = 1,
- ...) {
-
- ndone <- object at Nmif
-
- obj <- mif(
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/pomp -r 864
More information about the pomp-commits
mailing list