Question on PMCAR with multiple likelihoods

[BasMatsuura]

Dear Professor Gómez-Rubio,

I am trying to implement your multivariate effects, particularly the PMCAR, for multiple likelihoods using INLA. Your book helped me figure out how to structure univariate data for multiple likelihoods, but I fail to understand how to do so for multivariate data. My aim is to then fit the PMCAR effect to it, using the copy() function to retain the hyperparameters across the likelihoods.

Have you any idea whether this would be possible? And if so, any suggestions on how to structure the data when passing it to INLA and the inla.PMCAR effect? Right now I am stuck performing trial and error to see if it runs with different configurations. Was hoping you had already given this a try!

Yours sincerely,

Bas

[becarioprecario]

Dear Bas,

I am trying to implement your multivariate effects, particularly the PMCAR, for multiple likelihoods using INLA. Your book helped me figure out how to structure univariate data for multiple likelihoods, but I fail to understand how to do so for multivariate data. My aim is to then fit the PMCAR effect to it, using the copy() function to retain the hyperparameters across the likelihoods.

I am sorry but I do not understand your question. The examples provided in the paper are for multivariate data. Depending on the structure of the data and the model, data can be modeled using a single likelihood or several likelihoods. If you could be more specific , then I think that I will be able to help you.

Best,

Virgilio

[BasMatsuura]

Dear Virgilio,

Thank you for getting back to me. I did not have time to reply yet. Just started a new job. Apologies.

The copy-function might be unnecessary. I played with the regular IAR spatial effect, and it appears that without copy() R-INLA fits the model I am trying to use. That should make things easier.

For context: I am using ideas from your paper on modal Gibbs with INLA to fit a mixture model. Earlier I tried using the INLA with AMIS approach you published with dr. Berild, but it appears that does not collect samples efficiently enough. So I abandoned that. Specifically, I am trying to model each multivariate response, y_i, as a sum of an intercept and an MCAR spatial effect, where I assume the MCAR is shared by the mixture components. Basically, the model is given by:

y_i | z_i = k \sim MVN (\mu_k + \Phi_i, \Sigma),

where k is the specific mixture component, \mu_k are component specific intercepts, and the \Phi term is a general mutlivariate spatial effect. (I am not quite sure how INLA estimates \Sigma if at all.) This is how I am currently trying to pass the data to get this done with INLA:

G = 3
n = 625
p = 3

alpha.min = 0
alpha.max = 1
W <- as(W, "Matrix")
mcar_data <- c(data$Y)
intercept_reg <- rep(c("y1", "y2", "y3"), each = n)
z_rep <- rep(data$z, G)
idx_reg <- 1:length(mcar_data)

model_mcar <- inla.MCAR.model(k = p, W = W, alpha.min = alpha.min, alpha.max = alpha.max)

%% Fit the model %%

Y_mix <- matrix(NA, ncol = G, nrow = n*p)
intercept_mix <- matrix(NA, ncol = G, nrow = n*p)

for(i in 1:G) {
    Y_mix[which(z_rep == i), i] <- mcar_data[which(z_rep == i)]
    intercept_mix[which(z_rep == i), i] <- intercept_reg[which(z_rep == i)]
  }

r <- inla(Y ~ -1 + Intercept + f(idx, model = model_mcar),
            data = list(Y = Y_mix, Intercept = intercept_mix, idx = idx_reg), family = rep("gaussian", G),
            control.compute = list(config = TRUE),
            control.predictor = list(compute = FALSE),
            verbose = FALSE)

Here data$z is just a vector of integers (e.g. c(1,2,2,1,3,1,3) etc.) that the simulation is based on, so I try to fit the model under the true underlying classes. I suspect I am doing something wrong in passing the data. The data has been simulated from the generating process, so I don't think that's causing issues. Using this code I get the error:

Error in inla.core.safe(formula = formula, family = family, contrasts = contrasts,  : 
  replacement has 5625 rows, data has 1875

Here 1875 is derived from n=625 and p=3. Thank you for looking into my issue!

Regards,

Bas

[BasMatsuura]

Making something reproducible would involve sending more information. Hopefully this is enough for you to get an idea of what I am running into. Please let me know if not.