Missing Values • INLAvaan

INLAvaan handles missing data in one of two ways: listwise deletion (default, i.e. uses all complete cases) or Full Information Maximum Likelihood (FIML; missing = "ML").

Simulate Data

library(INLAvaan)
mod <- "
    # Latent variable definitions
    ind60 =~ x1 + x2 + x3
    dem60 =~ y1 + y2 + y3 + y4
    dem65 =~ y5 + y6 + y7 + y8

    # Latent regressions
    dem60 ~ ind60
    dem65 ~ ind60 + dem60

    # Residual correlations
    y1 ~~ y5
    y2 ~~ y4 + y6
    y3 ~~ y7
    y4 ~~ y8
    y6 ~~ y8
  "
dat <- lavaan::PoliticalDemocracy

# Simulate missingness (MCAR)
set.seed(221)
mis <- matrix(rbinom(prod(dim(dat)), 1, 0.99), nrow(dat), ncol(dat))
datmiss <- dat * mis
datmiss[datmiss == 0] <- NA

Listwise Deletion

fit1 <- asem(mod, datmiss, meanstructure = TRUE)
#> ℹ Mode finding and Hessian computation.
#> ✔ Posterior mode and Hessian. [295ms]
#> 
#> ℹ Performing VB correction.
#> ✔ VB correction; mean |δ| = 0.190σ. [937ms]
#> 
#> ⠙ Fitting 0/42 skew-normal marginals.
#> ⠹ Fitting 14/42 skew-normal marginals.
#> ✔ Fit 42/42 skew-normal marginals. [2.4s]
#> 
#> ⠙ Posterior sampling and summarising.
#> ✔ Summarise 1000 posterior draws. [1s]
#> 
#> ℹ Fit measures: PPP, DIC.
fit1@Data@nobs[[1]] == nrow(datmiss[complete.cases(datmiss), ])
#> [1] TRUE
print(fit1)
#> INLAvaan 0.3.2.9003 ended normally after 71 iterations
#> 
#>   Estimator                                      BAYES
#>   Optimization method                           NLMINB
#>   Number of model parameters                        42
#> 
#>                                                   Used       Total
#>   Number of observations                            35          75
#> 
#> Model Test (User Model):
#> 
#>    Marginal log-likelihood                    -795.666 
#>    PPP (Chi-square)                              0.509
coef(fit1)
#>    ind60=~x2    ind60=~x3    dem60=~y2    dem60=~y3    dem60=~y4    dem65=~y6 
#>        1.808        1.750        0.943        0.805        1.472        1.061 
#>    dem65=~y7    dem65=~y8  dem60~ind60  dem65~ind60  dem65~dem60       y1~~y5 
#>        0.721        1.354        0.913        0.567        1.080        0.470 
#>       y2~~y4       y2~~y6       y3~~y7       y4~~y8       y6~~y8       x1~~x1 
#>        1.844        3.596        0.607       -0.668        1.283        0.071 
#>       x2~~x2       x3~~x3       y1~~y1       y2~~y2       y3~~y3       y4~~y4 
#>        0.144        0.424        1.698        7.790        4.086        2.786 
#>       y5~~y5       y6~~y6       y7~~y7       y8~~y8 ind60~~ind60 dem60~~dem60 
#>        1.537        6.752        2.038        3.983        0.506        1.439 
#> dem65~~dem65         x1~1         x2~1         x3~1         y1~1         y2~1 
#>        0.139        5.424        5.534        4.107        7.250        6.635 
#>         y3~1         y4~1         y5~1         y6~1         y7~1         y8~1 
#>        8.307        6.834        6.639        5.224        8.267        6.157

Full Information Maximum Likelihood (FIML)

fit2 <- asem(mod, datmiss, missing = "ML", meanstructure = TRUE)
#> ℹ Mode finding and Hessian computation.
#> ✔ Posterior mode and Hessian. [593ms]
#> 
#> ℹ Performing VB correction.
#> ✔ VB correction; mean |δ| = 0.164σ. [1.4s]
#> 
#> ⠙ Fitting 0/42 skew-normal marginals.
#> ⠹ Fitting 9/42 skew-normal marginals.
#> ⠸ Fitting 28/42 skew-normal marginals.
#> ✔ Fit 42/42 skew-normal marginals. [6.7s]
#> 
#> ⠙ Posterior sampling and summarising.
#> ⠹ Computing fit indices (PPP/DIC).
#> ✔ Summarise 1000 posterior draws. [1.8s]
#> 
#> ℹ Fit measures: PPP, DIC.
print(fit2)
#> INLAvaan 0.3.2.9003 ended normally after 91 iterations
#> 
#>   Estimator                                      BAYES
#>   Optimization method                           NLMINB
#>   Number of model parameters                        42
#> 
#>   Number of observations                            75
#>   Number of missing patterns                        19
#> 
#> Model Test (User Model):
#> 
#>    Marginal log-likelihood                   -1393.186 
#>    PPP (Chi-square)                              0.047
coef(fit2)
#>    ind60=~x2    ind60=~x3    dem60=~y2    dem60=~y3    dem60=~y4    dem65=~y6 
#>        2.214        1.835        0.650        0.791        0.954        1.026 
#>    dem65=~y7    dem65=~y8  dem60~ind60  dem65~ind60  dem65~dem60       y1~~y5 
#>        1.050        1.287        1.299        0.524        0.756        0.482 
#>       y2~~y4       y2~~y6       y3~~y7       y4~~y8       y6~~y8       x1~~x1 
#>        0.898        3.316        0.301        0.365        1.294        0.084 
#>       x2~~x2       x3~~x3       y1~~y1       y2~~y2       y3~~y3       y4~~y4 
#>        0.136        0.511        1.699        7.448        3.408        2.945 
#>       y5~~y5       y6~~y6       y7~~y7       y8~~y8 ind60~~ind60 dem60~~dem60 
#>        1.794        5.910        2.073        3.682        0.462        4.558 
#> dem65~~dem65         x1~1         x2~1         x3~1         y1~1         y2~1 
#>        0.206        5.060        4.791        3.557        5.462        5.787 
#>         y3~1         y4~1         y5~1         y6~1         y7~1         y8~1 
#>        7.158        5.253        5.356        4.092        6.855        4.422

plot(
  coef(fit1),
  coef(fit2),
  xlab = "Listwise Deletion Estimates",
  ylab = "FIML Estimates"
)
abline(0, 1)

Model criteria under FIML

loo() and waic() work directly on a FIML fit. Each unit is scored on the entries it actually has – the observed-data predictive $\log p(y_{i,\text{obs}} \mid D_{-i})$ , with the full row deleted from the conditioning set – so a case with more missing entries contributes a smaller log-likelihood term and a smaller score, self-weighting in the expected log predictive density. The missing-at-random assumption that justifies FIML estimation also justifies this predictive score.

loo(fit2)
#> ── Leave-one-subject-out ────────────────────────── 75 subjects, second-order ──
#> 
#>          Estimate   SE
#> elpd_loo  -1287.1 36.6
#> p_loo        38.6  3.2
#> looic      2574.2 73.2
#> 
#> ── Curvature check ─────────────────────────────────────────────────────────────
#> 
#>   first-to-second-order gap        23.4
#>   pD/2 (trace)                     19.2
#>   excess over pD/2 (trace)       +22.1%
#> 
#> ℹ The gap approaches pD/2 (trace) from above. A large excess says the
#>   second-order expansion has not settled over the sample.

Comparing two FIML fits with compare(..., loo = TRUE) is valid only when they share the same observed entries – the same data and the same missingness pattern – since each unit is scored on the entries it has. See the cross-validation article for the Taylor case-deletion method itself.

Two-level FIML fits are also supported: they are scored per cluster (leave-one-cluster-out), each cluster contributing its observed-data marginal likelihood. The per-row deletion diagnostic (type = "loso") is not available under missing data.