Re: [r-inla] GMRFLib_2order_approx: rescue NAN/INF values in logl for idx

[Helpdesk]

update to the most recent testing version

try adding this option

r <- inla( ...., safe=TRUE)



On Tue, 2022-01-25 at 02:36 -0800, Michel Valette wrote:
> Good morning,
> I am trying to run a spatio-temporal model for understanding the
> impact of several covariates on occurrence of point data. For doing
> this, I am using the inlabru package and so far, used only fixed
> effects for fitting the covariate in the model. I am running this
> model on a cluster since the dataset is quite large (I can join a
> subsample of data if required).
> My model seemed to run past until last week when I updated the R-INLA
> version and updated R to version 4.1. When I am including only the
> fixed effects or only the spatio-temporal interaction, my model seems
> to be running fine. However, when I include both type of effects, I
> have a lot of error message indicating:
> GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=2161
> GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=926
> GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=1063
> Eventually, the model run until returning this message:
>                GitID: file: error-handler.c 
> a42508e2a097af0d8b23b353dc7f0038f5ffc12c - Wed Jan 19 12:18:35 2022
> +0100
>             Error:12 Reason: The Newton-Raphson optimizer did not
> converge
>                Message: Condition `lambda < 1.0 / lambda_lim' is not
> TRUE
>                Line:2543 Function:
> GMRFLib_init_GMRF_approximation_store__intern
> sh: line 1: 3830997 Aborted                
> '~/R/library/INLA/bin/linux/64bit/inla.mkl.run' -s -v -t256:1 -B0 -P
> classic '/tmp/pbs.5024938.pbs/RtmpA5N48w/file3a65831f24537a/Model.ini'
> > '/tmp/pbs.5024938.pbs/RtmpA5N48w/file3a65831f24537a/Logfile.txt'
> Warning in iinla(model = info[["model"]], lhoods = info[["lhoods"]],
> options = info[["options"]]) :
>   iinla: Problem in inla: Error in inla.inlaprogram.has.crashed() :
>   The inla-program exited with an error. Unless you interupted it
> yourself, please rerun with verbose=TRUE and check the output
> carefully.
>   If this does not help, please contact the developers at
> <he...@r-inla.org>.
> iinla: Problem in inla: 1: lgcp(formula.1, modis2, ips, domain =
> list(coordinates = mesh,
>        year = modis2$year), options = list(control.fixed = list(mean =
> list(0),
>        prec = list(0.01))))
> 2: bru(components, lik, options = options)
> 3: iinla(model = info[["model"]], lhoods = info[["lhoods"]], options =
> info[[ [...]
> 4: try_callstack(...)
> 5: do.call(INLA::inla, inla.options.merged, envir =
> environment(model$effects))
> 6: (function (formula = NULL, family = "gaussian", contrasts = NULL,
>        data = NULL, quantiles = c(0.025, 0.5, 0.975), E = NULL,
>        offset = NULL, scale = NULL, weights = NULL, Ntrials = NULL,
>        strata = NULL, lp.scale = NULL, link.covariates = NULL, verbose
> = FALSE,
>     [... truncated]
> iinla: Giving up and returning last successfully obtained result for
> diagnostic purposes.
>
> I was wondering if you have any clue on where I might have an issue? I
> have tried to include different type of matern  model(both simple or
> with penalized complexity priors), but the model seem to keep
> failing. 
> Best regards,
>
> Michel  
>  
> sessionInfo()
> R version 4.1.2 (2021-11-01)
> Platform: x86_64-conda-linux-gnu (64-bit)
> Running under: CentOS Linux 8 (Core)
>  
> Matrix products: default
> BLAS/LAPACK:
> /rds/general/user/mlv120/home/anaconda3/envs/Renv2/lib/libopenblasp-
> r0.3.18.so
>  
> locale:
>  [1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C             
>  [3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_GB.UTF-8   
>  [5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_GB.UTF-8  
>  [7] LC_PAPER=en_GB.UTF-8       LC_NAME=C                
>  [9] LC_ADDRESS=C               LC_TELEPHONE=C           
> [11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C      
>  
> attached base packages:
> [1] parallel  stats     graphics grDevices utils     datasets 
> methods 
> [8] base    
>  
> other attached packages:
> [1] raster_3.5-11 inlabru_2.4.0 sf_1.0-5     INLA_22.01.19 sp_1.4-
> 6    
> [6] foreach_1.5.1 Matrix_1.4-0
>  
> loaded via a namespace (and not attached):
>  [1] Rcpp_1.0.8         magrittr_2.0.1     splines_4.1.2     
> units_0.7-2      
>  [5] tidyselect_1.1.1   lattice_0.20-45    R6_2.5.1          
> rlang_0.4.12     
>  [9] fansi_1.0.2        dplyr_1.0.7        tools_4.1.2       
> grid_4.1.2       
> [13] utf8_1.2.2         KernSmooth_2.23-20 terra_1.5-12      
> e1071_1.7-9      
> [17] DBI_1.1.2          ellipsis_0.3.2     iterators_1.0.13  
> class_7.3-20     
> [21] tibble_3.1.6       lifecycle_1.0.1    crayon_1.4.2      
> purrr_0.3.4      
> [25] vctrs_0.3.8        codetools_0.2-18   glue_1.6.0        
> proxy_0.4-26     
> [29] compiler_4.1.2     pillar_1.6.4       generics_0.1.1    
> classInt_0.4-3   
> [33] pkgconfig_2.0.3  
>  
> --
> You received this message because you are subscribed to the Google
> Groups "R-inla discussion group" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to r-inla-discussion...@googlegroups.com.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/r-inla-discussion-group/8e901be9-52ed-4dcd-9dd9-7638189e9ae0n%40googlegroups.com
> .

--
Håvard Rue
he...@r-inla.org

[Michel Valette]

Dear Pr Håvard,

Thank you for the quick answer. I have tried the suggested improvement, making sure it was using the last version of INLA (22.01.19) and oncluding the suggested improvement, using the inlabru call:

r <- lgcp(formula.1, modis2,ips,domain=list(coordinates=mesh,year=modis2$year),options = list(safe=TRUE))

and

r <- lgcp(formula.1, modis2,ips,domain=list(coordinates=mesh,year=modis2$year),options = list(safe=TRUE, inla.mode="experimental" ))

But the same message of error are appearing again.

GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=260
GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=261

if you want to give a look I can upload some data.

Yours faithfully,

Michel

[INLA help]

The warning itself is ok, the question is if the program is able to recover 

—

Haavard Rue

HelpDesk 

help@r-inla. org

To view this discussion on the web, visit https://groups.google.com/d/msgid/r-inla-discussion-group/a42880ba-249d-456a-861a-f02049759cdfn%40googlegroups.com.

[Michel Valette]

Alright, I'll let this run for a while and come back when I have the final output (previous versions took around 30 hours to run).

thank you so much for the answer!

Garanti sans virus. www.avast.com

[Helpdesk]

if the code is able to recover, then a few of these does not harm.

are the results ok, or does it all end in an error? if err, please send
data/code.

[Michel Valette]

Dear Pr. Håvard,

I have tried to let the code run for a while but this error message is multiplying (repeating on hundreds/thousands of lines), but I am not sure of what this is pointing out? Eventually, it crash and I have this message:

GitID: file: error-handler.c  a42508e2a097af0d8b23b353dc7f0038f5ffc12c - Wed Jan 19 12:18:35 2022 +0100
        Error:12 Reason: The Newton-Raphson optimizer did not converge
        Message: Condition `lambda < 1.0 / lambda_lim' is not TRUE
        Line:2543 Function: GMRFLib_init_GMRF_approximation_store__intern

 

Since this doesn't seem to occur when I include only the spde effect I guess it might be some issue of the interactions between the spde and other fixed effects.I also notice that the model work fine on a sub-sample but crash when I try to scale up.

I am joining the code and will send you the data by mail.

[Patricia Marques Moralejo Bermudi]

Dear Michele, 

I was wondering if you managed to figure out how to solve it, please? I have the same problem.

Thank you in advance,
Patricia

[Patricia Marques Moralejo Bermudi]

Michel*  my apologies for the typo

[Helpdesk (Haavard Rue)]

the new version of R-INLA should be better to avoid this. Its usually is
because the data and the model does not fit well...

On Thu, 2023-01-19 at 09:42 -0800, 'Patricia Marques Moralejo Bermudi'

[aynaz lotfata]

Can you please let me know how you resolved this error "GMRFLib_2order_approx: rescue NAN/INF values in logl for idx=323838"?

thanks

[Helpdesk (Haavard Rue)]

its just a warning. usually it will be able to fix it, itself.

unless this appears all the time, its ok

[Dmytro Perepolkin]

Are the indices mentioned by this warning message referring to the observation in my data? I see the same three indices mentioned several times in the warning message. 

If so, I suspect it could be observations on the edge of the mesh. Is there a way to diagnose what happens at those index locations in the data?

[Helpdesk (Haavard Rue)]

its an overflow in the likelihood contribution from that linear
predictor (this index is in is 0-based indexing wheras the R-indexing is
1-based...).

usually it will be able to correct itself

[Roel van Klink]

I have approximately 4000 of these warnings, after which Inla crashed and dumped the 16 gb core. Interestingly, a few minutes earlier, it also dumped the core 2 times (each 2gb)

with the following error:

GitId: 882c43d68fb015d72b699bcc458cf277dfe110c2 - Mon Apr 24 16:19:16 2023 +0100
Error:3 Reason: Matrix is (numerical) singular
Line:264 Function: GMRFLib_comp_posdef_inverse

Which matrix is this referring to? Could there be a problem with the way I specified my random effects?

thanks

[Helpdesk (Haavard Rue)]

it seems to be here:

inverse_hessian = Calloc(ISQR(nhyper), double);
Memcpy(inverse_hessian, hessian, ISQR(nhyper) * sizeof(double));
GMRFLib_comp_posdef_inverse(inverse_hessian, nhyper);



so something with hyperparameters not going well.

can you share code/data so I can run it here? this is, I guess, a user input
error

[Roel van Klink]

I can send data and code, but could the problem be that I coded my 3 nested random slopes and intercepts like this:

level 1: 1378

level 2: 1378:281

level 3: 1378:281:24

I found this somewhere online.

dataframe of 12000 lines.

thanks

[Helpdesk (Haavard Rue)]

The singularity is introduced in the model spesification... so it is something
you spesify that I guess you do not want

[Roel van Klink]

I have attached my dataframe and the r-code I use to run it on the HPC (4x50gb cores)

I'll need to change the threads to 4, but any other tips on how to make this run more efficiently would also be welcome.

[Helpdesk (Haavard Rue)]

I attached a revised code.

comments

- 'sd.res' was to high

- when using replicate=r then its replcated max(r) times. if r=c(1,4,100), then
this is 100 times while only 3 is strictly needed. the code fix this

- when using iYear as weights in the iid, its easier if its iYear is scaled.

- changed the prior for the 'ar1' so its contract to a constant, you may want to
change it if you like

- added proper priors for the fixed effects, in the debugging phase and it might
not be needed anymore

- I added argument 'E' to act as an overall offset, as the mean is so off zero.
you can also use argument E=E for the poisson or offset=log(E)


let me know how this goes

btw, now its faster because of the 'replicate' thing

best
H