Thanks for the proposed solution! It looks promising and I would try it!
If I do this with MCMC in PyMC, it seems that PyMC would automatically draw sample from the distribution to impute those missing values making it much easier to write the code.
The idea of taking expectations w.r.t W and using those as the missing values applies to VI, at least in theory (e.g., Blei, Kucukelbir, McAuliffe, 2017), when the variational distributions are all chosen as exponential family.