Setting high precision for a nonlinear latent variable: a comment

[guidofioravanti]

Dear Finn,

I have another question about model convergence.

I am testing a datafusion model similar to the one described in the paper of Villejo et al. (2024) "A Data Fusion Model for Meteorological Data using the INLA-SPDE method", available through arxiv.org.

The objective is to estimate a common latent field between ground (sparse) observations and model-based continuous gridded surfaces in Europe. As in the Villejo's paper, a multiplicative bias parameter alpha is introduced to take into account a bias between the two data sources. I understand that inlabru is able to manage this non-linearity through a Taylor expansion.

For this multiplicative parameter, I set in the component definition :alpha(1, mean.linear=1 and prec.linear=950).

My model converges only if the precision of the alpha parameter (prec.linear) is set to a high value (above 900). The higher the density of the mesh the higher the precision of the parameter must be to get model convergence.

Despite the high precision, the posterior marginal distribution of alpha has a reasonable value mean value of around 1.4, So, the high precision of alpha, on one hand, is necessary for model convergence, on the other hand does not impede the posterior distribution estimation of alpha.

I have examined the bru_convergence_plot and the posterior distributions of the fixed and random effects. All the diagnostics seem to be good. The model convergence met. And the alpha parameter has a reasonable posterior distribution.

Do you have any comment about the high precision for alpha for getting model convergence?

Best,
Guido

[finnlindgren]

Hi Guido,
I missed this when you posted it, sorry!

Join modelling with non-linear predictors (through the product of the "bias parameter" and the random field) can be sensitive to the specifics of the problem. E.g. if one data source dominates the other, or the model is mis-specified. For some problems, the joint latent mode is in a "bad" place, pushing the bias parameter away from zero, which perhaps is what may have happened in your case. While this happens, the combined effect my still be estimated appropriately; the problem can be a partial non-identifiability of the individual components. We are investing this further, and I would be interested in having a look at your specific case, to see if it's something that can be addressed, and/or properly explained.

I notice you set mean.linear = 1; your INLA version used for the experiment would be useful to know, as some versions had a bug in the internal variational Bayes code that didn't handle non-zero latent prior means properly (this has since been fixed).

The bru_convergence_plot() output would also be useful (if possible, using the latest inlabru development version, as the plotting code has been improved/updated recently), and the bru_log() output.

[guidofioravanti]

Hi Finn,

this time I was the one who mised the post. Thanks for your attention to may posts. Here at my reserach centre we are eager to test INLA for operational drought monitoring. This is, more ore less, the topic of my research.

The INLA version is INLA_24.02.09 built 2024-02-09 03:35:28 UTC.
The "inlabru" version is "2.10.1".

I think that what you suggest could be the case for my dataset: the data are not sufficient informative to identify the scaling parameter.

An interesting point is this one.

The following model (the one I had been tested when I got in touch with you) does not converge (unless to put a very high prior on the precision of the alpha parameter, as explained in my previous post):

      lik_data <- like(family = "gaussian", formula = spi ~ (Intercept+dem+gradient)+spde_data,data =subData,include=c("spde_data","spde_era5","dem","gradient","Intercept"))

      lik_era5 <- like(family = "gaussian", formula = spi ~ alpha* (Intercept+spde_era5+spde_data+dem+gradient),data=df_monthly_reanalysis,include=c("spde_era5","dem","gradient","Intercept"))    

This model instead converges..and quite quickly (but the validation scores are worse than the ones from an identical model with alpha=1):

      lik_data <- like(family = "gaussian", formula = spi ~alpha* (Intercept+dem+gradient+spde_era5)+spde_data,data =subData,include=c("spde_data","spde_era5","dem","gradient","Intercept"))

      lik_era5 <- like(family = "gaussian", formula = spi ~ (Intercept+spde_era5+dem+gradient),data=df_monthly_reanalysis,include=c("spde_era5","dem","gradient","Intercept"))    

      cmp1<- ~  -1+
        Intercept(1)+
        alpha(1,model="linear",mean.linear=1,prec.linear=1)+ #950
        dem(elev,model="linear")+
        spde_data(geometry,model=matern_data)+ 
        spde_era5(geometry,model=matern_era5)+ 
        gradient(eval_spatial(monthly_gradient,where=.data.,layer=yymmdd),model="linear") 
    

In the above code, subData is a set of stations (119) covering part of Italy and northern Europe (spde_data, the corresponding SPDE field); df_monthly_reanalysis is the data for the centroids of a model (more than 2000 points, spde_era5 the correspondin SPDE field).

Do you have any idea/suggestion, why the first model with alpha does not converge (unless an high prior precision) and the second one convergences without [issues?]

Thanks
Guido

[finnlindgren]

Without looking into the details, my current best guess is what I wrote in my previous reply; having alpha * in one or the other model is a fundamentally different model (one alpha is the reciprocal of the other, essentially), so it may be just "luck" if one works and not the other.

Note: recent versions of inlabru will auto-detect which components are needed by each predictor, so there's no longer a need to supply the include = c(...) argument.