Ignorability in Rubin's theory of missing data mechanisms I am trying to understand Rubin's theory of bayesian inference with
missing data, specifically how the missing data mechanism affects the
inference on a superpopulation parameter. The theory is exposed for
example in the chapter 7 of [1].
The complete $N$-dimensional data vector $y$ is splitted into an
observed and unobserved components: $y = (y_{\mathrm{obs}},
y_{\mathrm{nob}})$ and $I$ is the $N$-vector of indicators whose
components $I_{i}$ equal 0 when $y_{i}$ is unobserved and 1 when it is
observed. The sample space is thus the product of the $y$ space and
the $I$ space. For simplicity, I consider a missing data mechanism
$p(I ~|~ y_{\mathrm{obs}})$ which does not depend on a fully observed
covariate $x$ or an unknown parameter $\phi$.
The key result in Rubin is that when $p(I ~|~ y) = p(I ~|~
y_{\mathrm{obs}})$, the missing data mechanism is ignorable
and the analysis can proceed as usual. This can be seen in the
following development:
$$
\begin{align*}
  p(\theta ~|~ y_{\mathrm{obs}}, I) &~=~ \frac{p(\theta) \int p(y~|~\theta)\, p(I ~|~
    y)~d y_{\mathrm{nob}}~}{\int\int p(y~|~\theta)\,p(\theta) ~d\theta p(I ~|~ y)~d y_{\mathrm{nob}}~} \\
\end{align*}
$$
If $p(I ~|~ y) = p(I ~|~ y_{\mathrm{obs}})$, it can be taken out of
the integrals and cancelled out such that the posterior distribution
$p(\theta ~|~ y_{\mathrm{obs}}, I) = p(\theta ~|~ y_{\mathrm{obs}})$.
So far so good. But what I don't understand is the following
derivation that contradicts the preceding result:
$$
\begin{align*}
  p(\theta ~|~ y_{\mathrm{obs}}, I)  &~=~ \frac{p(y_{\mathrm{obs}}, I, \theta)}{p(y_{\mathrm{obs}}, I)} \\
  &~=~ \frac{p(y_{\mathrm{obs}} ~|~ \theta) p(\theta)
    p(I ~|~ y_{\mathrm{obs}})}{\int p(y_{\mathrm{obs}} ~|~\theta)p(\theta) ~d\theta~ p(I ~|~ y_{\mathrm{obs}})} \\
  &~=~ \frac{p(y_{\mathrm{obs}} ~|~ \theta)p(\theta)}{\int p(y_{\mathrm{obs}} ~|~\theta)p(\theta)
    ~d\theta} \\
  &~=~ p(\theta ~|~ y_{\mathrm{obs}}).
\end{align*}
$$
Here the sampling mechanism cancels out no matter what. It does not
matter if $p(I ~|~ y) = p(I ~|~ y_{\mathrm{obs}})$. I must be doing
something wrong, but I can't see what.
[1] Gelman, A., Carlin, J. B., Stern, H. S., & Rubin,
D. B. (2004). Bayesian Data Analysis (2nd ed.). Boca Raton: Chapman &
Hall/CRC.
 A: In the second row of your derivation, you seem to assume 
\begin{equation}
p(y_{obs},I,\theta) = p(y_{obs}|\theta)p(\theta)p(I|y_{obs}),
\end{equation}
but this is not in general true - the correct factorization would be 
\begin{equation}
p(y_{obs},I,\theta) = p(y_{obs}|\theta)p(\theta)p(I|y_{obs},\theta).
\end{equation}
Without the assumption $p(I|y) = p(I|y_{obs})$, $\theta$  'impacts' $I$ via $y_{nob}$ even when conditioning on $y_{obs}$, and thus $p(I|y_{obs},\theta)\neq p(I|y_{obs})$.
A: Juho Kokkala provided the right answer and I will just add an example
here, based on p. 121 of Little and Rubin, 2002, 2nd ed.
The correct factorization of the missingness mechanism is indeed $p(I
~|~ y_{\mathrm{obs}}, \theta)$. More specifically, it can be further
factored into the components depending on $y_{\mathrm{obs}}$ and those
which depend on $y_{\mathrm{nob}}$. Since the latter are unobserved
data, we rely on their predictive distribution for a given $\theta$:
$$
\begin{align*}
  p(I ~|~ y_{\mathrm{obs}}, \theta) &~=~ p(I_{\mathrm{obs}} ~|~
  y_{\mathrm{obs}}) p(I_{\mathrm{nob}} ~|~ \theta) \\
  &~=~\prod_{i=1}^{n} p(I_{i} ~|~ y_{i}) \prod_{i=n+1}^{N} p(I_{i} ~|~ \theta)
\end{align*}
$$
Note the $n$ first components of $y$ constitute $y_{\mathrm{obs}}$,
and the complete data has $N$ observations.
Now the example. Imagine that $y \sim \mathrm{Exp}(\theta)$ and our missing
data is due to the censoring mechanism
$$
\begin{align}
  p(I_{i} ~|~ y_{i}) ~=~
  \begin{cases}
    1, & \mathrm{if}~I_{i}= 1 \mathrm{~and~} y_{i} \ge c,
    \mathrm{~or~} I_{i} = 0 \mathrm{~and~} y_{i} < c \\
  0, & \mathrm{otherwise.}
  \end{cases}
\end{align}
$$
Then we have the likelihood
$$
\begin{align*}
  p(y_{\mathrm{obs}}, I ~|~ \theta) &~=~ p(y_{\mathrm{obs}} ~|~ \theta)
  p(I_{\mathrm{obs}} ~|~ y_{\mathrm{obs}}) p(I_{\mathrm{nob}} ~|~ \theta) \\
  % &~=~ \prod_{i=1}^{n} p(y_{i}, I_{i} ~|~ \theta) \prod_{i=n+1}^{N}
  %   p(I_{i} ~|~ \theta) \\
  &~=~ \prod_{i=1}^{n} p(y_{i} ~|~ \theta) \mathrm{Pr}(y_{i} < c ~|~
  y_{i}) \prod_{i=n+1}^{N} \mathrm{Pr}(y_{i} \ge c ~|~ \theta) \\
  &~=~ \theta^{-n} \exp\left (-\sum_{i=1}^{n} \frac{y_{i}}{\theta} \right)
  \exp\left (-\frac{(N-n)c}{\theta}\right).
\end{align*}
$$
Using the hypothesized missing data mechanism and the hypothesized
distribution of $y_{\mathrm{nob}}$, we compute the conditional
density of the unobserved components of $I$.
