tests/measure.R

In kidusasfaw/pomp: Statistical Inference for Partially Observed Markov Processes

options(digits=3)

library(pomp2)
library(magrittr)

ou2() -> ou2

po <- window(ou2,end=10)

set.seed(3434388L)
simulate(po,nsim=5,format="arrays") -> y
y %>% extract2("states") -> x
y %>% extract2("obs") %>% extract(,1,) -> y
t <- time(po)
p <- coef(po)

dmeasure(po,x=x,y=y,times=t,params=p) -> L
dmeasure(po,x=x,y=y,times=t,params=p,log=T) -> ll
dmeasure(po,x=x,y=y,times=t,log=T,
  params=array(data=p,dim=length(p),dimnames=list(names(p)))) -> ll2
dmeasure(po,x=array(data=x,dim=c(dim(x),1,1),
  dimnames=list(rownames(x),NULL,NULL,NULL,NULL)),
  y=y,times=t,log=T,
  params=array(data=p,dim=length(p),dimnames=list(names(p)))) -> ll3
dmeasure(po,x=array(data=x[,1,1],dim=nrow(x),dimnames=list(rownames(x))),
  y=y[,1],times=t[1],params=p,log=TRUE) -> ll4
stopifnot(
  all.equal(ll[,1:3],log(L[,1:3])),
  identical(dim(ll),c(5L,10L)),
  all.equal(ll,ll2),
  all.equal(ll,ll3),
  all.equal(ll[1,1],ll4[1,1])
)

try(dmeasure("ou2",x=x,y=y,times=t,params=p))
try(dmeasure(x=x,y=y,times=t,params=p))
try(dmeasure(x,y=y,times=t,params=p))
try(dmeasure(po,x=x,y=y,times=t))
try(dmeasure(po,x=x,y=y,params=p))
try(dmeasure(po,x=x,times=t,params=p))
try(dmeasure(po,y=y,times=t,params=p))
try(dmeasure(po,x=as.numeric(x),y=y,times=t,params=p))
try(dmeasure(po,x=x,y=as.numeric(y),times=t,params=p))
try(dmeasure(po,x=x,y=y,times=NULL,params=p))
try(dmeasure(po,x=x[,,1],y=y[,1,drop=FALSE],times=t[1],params=p))
invisible(dmeasure(po,x=x[,,1,drop=FALSE],y=y[,1],times=t[1],params=p))
stopifnot(
  all.equal(dmeasure(po,x=x[,1,,drop=FALSE],y=y,times=t,params=p),
    dmeasure(po,x=x[,1,],y=y,times=t,params=p))
)
try(dmeasure(po,x=x,y=y[1,,drop=FALSE],times=t,params=p))
try(dmeasure(po,x=x[1,,,drop=FALSE],y=y,times=t,params=p))
k <- which(names(p)=="tau")
try(dmeasure(po,x=x,y=y,times=t,params=p[-k]))

pp <- parmat(p,5)
try(dmeasure(po,x=x,y=y,times=t,params=pp[,1:3]))
dmeasure(po,x=x,y=y,times=t,params=pp) -> d
stopifnot(dim(d)==c(5,10),names(dimnames(d))==c("rep","time"))

rmeasure(po,x=x,times=t,params=p) -> y
stopifnot(
  dim(y)==c(2,5,10),
  names(dimnames(y))==c("variable","rep","time")
)

try(rmeasure("ou2",x=x,times=t,params=p))
try(rmeasure(x=x,times=t,params=p))
try(rmeasure(x,times=t,params=p))
try(rmeasure(po,x=x,times=t))
try(rmeasure(po,x=x,params=p))
try(rmeasure(po,x=as.numeric(x),times=t,params=p))
try(rmeasure(po,x=x,times=NULL,params=p))
try(rmeasure(po,x=x[,,1],times=t[1],params=p))
invisible(rmeasure(po,x=x[,,1,drop=FALSE],times=t[1],params=p))
try(rmeasure(po,x=x[1,,,drop=FALSE],times=t,params=p))
k <- which(names(p)=="tau")
try(rmeasure(po,x=x,y=y,times=t,params=p[-k]))

pp <- parmat(p,5)
try(rmeasure(po,x=x,times=t,params=pp[,1:3]))
rmeasure(po,x=x,times=t,params=pp) -> y
stopifnot(dim(y)==c(2,5,10),names(dimnames(y))==c("variable","rep","time"))

po %>% pomp(
  rmeasure=function(...)c(1,2,3),
  dmeasure=function(...,log)c(3,2)
) -> po1
try(po1 %>% rmeasure(x=x,times=t,params=p))
try(po1 %>% dmeasure(x=x,y=y[,1,],times=t,params=p))

po %>% pomp(
  rmeasure=NULL,
  dmeasure=NULL
) -> po1
po1 %>% rmeasure(x=x,times=t,params=p) %>% is.na() %>% stopifnot()
po1 %>% dmeasure(x=x,y=y[,1,],times=t,params=p) %>% is.na() %>% stopifnot()

sir() -> sir
po <- window(sir,end=0.5)

set.seed(3434388L)
simulate(po,nsim=5,format="arrays") -> y
y %>% extract2("states") -> x
y %>% extract2("obs") %>% extract(,1,,drop=FALSE) -> y
t <- time(po)
p <- coef(po)

po %>% dmeasure(x=x,y=y,params=p,times=t,log=TRUE) -> d
po %>% rmeasure(x=x,params=p,times=t) -> yy

po %>% pomp(dmeasure=function(...,log)1) %>% dmeasure(x=x,y=y,params=p,times=t) -> d

try({
  pp <- p
  names(pp)[3] <- NA
  po %>% rmeasure(x=x,params=pp,times=t)
})

try({
  pp <- p
  names(pp)[3] <- ""
  po %>% rmeasure(x=x,params=pp,times=t)
})

try({
  pp <- p
  names(pp) <- NULL
  po %>% rmeasure(x=x,params=pp,times=t)
})