Skip to content

bug: hsar_cpp_arma_rho_0 crashes on toy data set for Nsim - burnin >= 2048 #5

@pol-db4drd2

Description

@pol-db4drd2

hsar crashes when W is NULL given enough iterations after burnin on a toy problem.

  • The message given is Error: max(): object has no elements
  • Failure does not occur if I set Nsim less than burnin + 2048 regardless of burnin
  • I haven't checked on what happens if I vary thinning
  • Failure does not occur if W is not NULL regardless of M

I've walked through it in the R debugger and can confirm the error happens at the call to hsar_cpp_arma_rho_0. I don't know how to debug C++ code, so I can't investigate this any deeper.

Reproducible example

set.seed(1234567890)

L  <-  0.25
R  <- -0.3

WL <- spdep::nb2mat(spdep::cell2nb(25, 25, "queen"))
WU <- spdep::nb2mat(spdep::cell2nb( 5,  5, "queen"))

us <- dplyr::tibble(I = 0:24)
us <- dplyr::mutate(
  us,
  Mu = rnorm(dplyr::n(), 0, 1),
   u = Matrix::solve(Matrix::Diagonal(25) - L * WU, .data$Mu)
)

testspace <- dplyr::tibble(
  i  = 1:625 - 1,
  j  = .data$i  %% 25,
  k  = .data$i %/% 25, # I guess this is in row-major order
  bj = 6 * ((.data$j / 24)^3 - (.data$j / 24)^2),
  bk = (.data$k / 24)^3 + (.data$k / 24)^2 + (.data$k / 24),
  J  = .data$j %/% 5,
  K  = .data$k %/% 5,
  I  = .data$J + 5 * .data$K
)

testspace <- dplyr::left_join(testspace, us)

testspace <- dplyr::mutate(
  testspace,
   X = stats::rnorm(dplyr::n(), 0, 1),
  My = rnorm(dplyr::n(), .data$bj + .data$bk - .data$X + .data$u, sqrt(2)),
   y = Matrix::solve(Matrix::Diagonal(625) - R * WL, .data$My)
)


del <- model.matrix(~ factor(I) - 1, testspace)

form <- y ~ j + k + X

baz <- HSAR::hsar(form, testspace, NULL, Matrix::Matrix(WU), Matrix::Matrix(del))

Expected behavior

Create an object named baz of class mcmc_hsar_rho_0

Observed behavior

Error: max(): object has no elements

Session info

> sessionInfo()
R version 4.4.1 (2024-06-14)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.10

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.28.so;  LAPACK version 3.12.0

locale:
 [1] LC_CTYPE=C.UTF-8       LC_NUMERIC=C           LC_TIME=C.UTF-8        LC_COLLATE=C.UTF-8     LC_MONETARY=C.UTF-8   
 [6] LC_MESSAGES=C.UTF-8    LC_PAPER=C.UTF-8       LC_NAME=C              LC_ADDRESS=C           LC_TELEPHONE=C        
[11] LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C   

time zone: Etc/UTC
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] HSAR_0.6.2     testthat_3.2.3

loaded via a namespace (and not attached):
 [1] xfun_0.52          htmlwidgets_1.6.4  devtools_2.4.5     remotes_2.5.0      lattice_0.22-7     vctrs_0.6.5       
 [7] tools_4.4.1        LearnBayes_2.15.1  spdep_1.3-11       sandwich_3.1-1     proxy_0.4-27       Matrix_1.7-3      
[13] KernSmooth_2.23-26 desc_1.4.3         lifecycle_1.0.4    compiler_4.4.1     deldir_2.0-4       brio_1.1.5        
[19] codetools_0.2-20   httpuv_1.6.16      htmltools_0.5.8.1  usethis_3.1.0      class_7.3-23       yaml_2.3.10       
[25] later_1.4.2        urlchecker_1.0.1   MASS_7.3-65        ellipsis_0.3.2     classInt_0.4-11    cachem_1.1.0      
[31] wk_0.9.4           spatialreg_1.3-6   sessioninfo_1.2.3  boot_1.3-31        multcomp_1.4-28    nlme_3.1-168      
[37] mime_0.13          digest_0.6.37      mvtnorm_1.3-3      sf_1.0-20          purrr_1.0.4        splines_4.4.1     
[43] rprojroot_2.0.4    fastmap_1.2.0      grid_4.4.1         cli_3.6.5          magrittr_2.0.3     survival_3.8-3    
[49] pkgbuild_1.4.7     TH.data_1.1-3      e1071_1.7-16       withr_3.0.2        promises_1.3.2     sp_2.2-0          
[55] spData_2.3.4       rmarkdown_2.29     zoo_1.8-14         evaluate_1.0.3     memoise_2.0.1      coda_0.19-4.1     
[61] shiny_1.10.0       knitr_1.50         miniUI_0.1.2       s2_1.1.7           profvis_0.4.0      rlang_1.1.6       
[67] Rcpp_1.0.14        xtable_1.8-4       glue_1.8.0         DBI_1.2.3          pkgload_1.4.0      rstudioapi_0.17.1 
[73] R6_2.6.1           fs_1.6.6           units_0.8-7       

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions