R Under development (unstable) (2024-02-17 r85935 ucrt) -- "Unsuffered Consequences" Copyright (C) 2024 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 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(spatPomp) Loading required package: pomp > set.seed(0) > > b2_U <- 4 > b2 <- bm2(U=b2_U,N=2,unit_specific_names="rho") > > b2_bpfilter <- bpfilter(b2,Np=5,block_size=1) > paste("bpfilter logLik for bm2 model:",logLik(b2_bpfilter)) [1] "bpfilter logLik for bm2 model: -12.9121635525516" > > # 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)) [1] "ibpf logLik for b2 model: -14.1817035202796" > > paste("kf logLik for b2:",bm2_kalman_logLik(b2)) [1] "kf logLik for b2: -13.2725611269272" > > ####################################################################### > # 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)) [1] "bm2 ibpf loglik: -13.2406696135" > > 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)) [1] "check ibpf on ibpfd_spatPomp: TRUE" > > 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)) [1] "check ibpf on bpfilterd_spatPomp: TRUE" > > > > > proc.time() user system elapsed 1.39 0.17 4.95