Nothing
library(spatPomp)
## test error message
try(bm2(U=3,N=2,unit_specific_names="rho",shared_names="sigma"))
try(bm2_kalman_logLik(bm2(U=1,N=2)))
b2_U <- 4
set.seed(0)
b2 <- bm2(U=b2_U,N=2,unit_specific_names="rho")
set.seed(0)
b2alt <- bm2(U=b2_U,N=2,shared_names=c("tau","sigma","X_0"))
if(any(obs(b2)!=obs(b2alt))) stop("two different ways to specify the same bm2 model should match")
b2_bpfilter <- bpfilter(b2,Np=5,block_size=1)
paste("bpfilter logLik for bm2 model:",logLik(b2_bpfilter))
# here there are no transformations so use small rw.sd. to avoid negatives
b2_rw_list <- rep(list(0.001),times=b2_U)
names(b2_rw_list) <-paste0("rho",1:b2_U)
b2_rw.sd <- do.call(rw_sd,b2_rw_list)
b2_units = seq_len(b2_U)
b2_nblocks = b2_U/2
b2_block_list = split(b2_units, sort(b2_units %% b2_nblocks))
b2_block_list <- lapply(b2_block_list, as.integer)
b2_ibpf <- ibpf(b2,
params=coef(b2),
sharedParNames=NULL,
unitParNames="rho",
Nbpf=2,
spat_regression=0.1,
Np=5,
rw.sd=b2_rw.sd,
cooling.fraction.50=0.5,
block_list=b2_block_list
)
paste("ibpf logLik for b2 model:",logLik(b2_ibpf))
paste("kf logLik for b2:",bm2_kalman_logLik(b2))
#######################################################################
# test ibpf with argument re-use and replacement
# check that results match after re-use and replacement
#######################################################################
set.seed(5)
b2_ibpf <- ibpf(b2,Np=5,block_size=1,Nbpf=2,
rw.sd=b2_rw.sd,
cooling.frac=0.5, spat_regression=0.1,
unitParNames="rho",sharedParNames=NULL)
paste("bm2 ibpf loglik: ",round(logLik(b2_ibpf),10))
set.seed(5)
b2_ibpf_repeat <- ibpf(b2_ibpf,params=coef(b2), unitParNames="rho",
sharedParNames=NULL,spat_regression=0.1)
paste("check ibpf on ibpfd_spatPomp: ",
logLik(b2_ibpf_repeat)==logLik(b2_ibpf))
set.seed(5)
b2_ibpf_bpfilterd <- ibpf(b2_bpfilter,Np=5,block_size=1,Nbpf=2,
rw.sd=b2_rw.sd,
cooling.frac=0.5, spat_regression=0.1,
unitParNames="rho",sharedParNames=NULL)
paste("check ibpf on bpfilterd_spatPomp: ",
logLik(b2_ibpf_bpfilterd)==logLik(b2_ibpf))
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.