# Multivariate meta-analysis cross-validation with missing data imputed during model fitting

**URL:** <https://discourse.mc-stan.org/t/multivariate-meta-analysis-cross-validation-with-missing-data-imputed-during-model-fitting/26596>\
**Category:** brms\
**Tags:** loo, meta-analysis, brms\
**Created:** [February 28, 2022, 4:25pm UTC](https://discourse.mc-stan.org/t/multivariate-meta-analysis-cross-validation-with-missing-data-imputed-during-model-fitting/26596 "2022-02-28T16:25:47Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![fborszcz](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/fborszcz/32/6946_2.png) [@fborszcz](https://discourse.mc-stan.org/u/fborszcz)\
**Post date:** [February 28, 2022, 4:25pm UTC](https://discourse.mc-stan.org/t/multivariate-meta-analysis-cross-validation-with-missing-data-imputed-during-model-fitting/26596/1 "2022-02-28T16:25:47Z")

</div>

Hello all,

I’m performing a meta-analysis with correlated data (various effect sizes from the same subjects), predictors with missing data and with measurement error (i.e., sd), and a combination of both. As shown in the example below:

```nohighlight
library(brms)

#covariance matrix of outcomes
cor_matrix <- function(x, r, v = rep(1, length(x)), na.rm = FALSE) {
  mat <- diag(v)
  se <- sqrt(v)
  se[is.na(se)] <- 0
  if (length(x) > 1L) {
    for (i in 2:nrow(mat)) {
      for (j in 1:(i-1)) {
        if (x[i] == x[j]) {
          mat[i, j] <- mat[j, i] <- r * se[i] * se[j] 
        }
      }
    }
  }
  dimnames(mat) <- list(1:nrow(mat), 1:ncol(mat))
  if (na.rm) {
    keep <- !is.na(diag(mat))
    mat <- mat[keep, keep]
  }
  mat
}

#data
dat <- data.frame(
  yi = rnorm(10),
  x_1 = c(1, 1, NA, NA, NA, 1, NA, 3, 5, 1),
  x_2 = rnorm(10),
  x_2_sd = runif(n = 10, min = 0.5, max = 1.5),
  x_3 = rnorm(10),
  effect_ID = factor(1:10),
  study_ID = factor(c(1:4,5,5,6,6,7:8)),
  vi = runif(n = 10, min = .25, max = .85)
)

#matrix
V <- cor_matrix(dat$study_ID, r = 0.9, v = dat$vi)

#model
m1_f <- 
  bf(yi|mi() ~ 1 + mi(x_1) + mi(x_2) + x_3 +
       (1|effect_ID) +
       fcor(V), 
     sigma = 1) +
  bf(x_1|mi() ~ 1) +
  bf(x_2|mi(x_2_sd) ~ 1) +
  set_rescor(FALSE)

m1 <- 
  brm(m1_f,
      family = gaussian(),
      data2 = list(V = V),
      data = dat,
      chains = 2,
      save_pars = save_pars())
        

```

Since I have a combination of other predictors, I want to check the best-fitting model. So, to cross-validate, I used only the complete data (excluded the missing data):

```nohighlight
#loo CV
m1_loo <- add_criterion(m1, newdata = na.omit(m1$data), criterion = "loo",
                        resp = "yi")

```

But as the variance-covariance matrix was inserted into the model when adding the “loo” criterion, the following error message appears:

```nohighlight
Error: Dimensions of 'M' for FCOR terms must be equal to the number of observations.

```

My question is if is there a way to cross-validate these models, given the circumstances described above?

- Operating System: Windows 10
- brms Version: 2.16.8
