Si faltan los datos de proceso es conocido y es MNAR, es posible obtener una estimación insesgada de parámetro?

Question

Estaba leyendo sobre los datos que faltan el manejo de las técnicas de "Análisis Estadístico con Los Datos que faltan de la Segunda Edición de" Little y Rubin. En el ejemplo 6.22 consideran que $y_1, y_2, ... y_n$ son yo.yo.d. exponencial de las variables aleatorias. es decir,$f(y_i) = \theta^{-1}e^{-\dfrac y \theta}$. Más $y_i$ sólo se registran si es menos de un corte de $c$. Supongamos $r$ $n$ observaciones se registran y el resto que falta. Correspondiente missingness indicador de $M_i$ es 0 si no se encuentra y 1 si se observa. Ellos muestran que la estimación basada en la articulación de la probabilidad de falta de proceso de datos y los datos observados da

$\hat{\theta} = \dfrac {\Sigma_{i=1}^{r} y_i + (n-r)c} r $

Mi pregunta es si esta estimación es imparcial?

Aquí es lo que yo he probado hasta ahora:

Desde la observó $y_i$ pertenecen a un tronco de una distribución exponencial con el apoyo $[0, c)$, $E(y_i) = \dfrac \theta {1-e^{-\dfrac c \theta}} (1-e^{-\dfrac c \theta}(1+\dfrac c \theta))$.

$\Rightarrow E(\hat{\theta}) = \dfrac \theta {1-e^{-\dfrac c \theta}} (1-e^{-\dfrac c \theta}(1+\dfrac c \theta)) + nc * E(\dfrac {1} {r}) - c$

Desde $r$, el número de datos observados que son menos de $c$ sigue una distribución binomial con probabilidad de éxito $p(y < c)$, $\dfrac {1} {r}$ es la inversa de la distribución binomial. Traté de calcular este valor esperado pero terminé concluyendo que es infinito. Así que me parece que el sesgo también infinito?

Enlace: Inversa Binomio valor Esperado

Las simulaciones que he realizado para comprobar unbiasedness del estimador de

También he realizado algunos simualations para tener una idea aproximada acerca de cómo el estimador es sesgado, y a partir de las simulaciones a mí me parece que puede ser imparcial y su varianza en inversamente proporcional a la frecuencia de corte $c$. Sin embargo yo no puedo concluir a partir de simulaciones. Sin embargo, aquí está mi código R:

library(ggplot2)

numObs = 100
len = 100
cutoff = 30
origEstimator = vector("numeric", len)
newEstimator = vector("numeric", len)
for(i in 1:len){
  sample = rexp(n = numObs, rate = 0.05)
  origEstimator[i] = mean(sample)

  truncSample = sample[sample<cutoff]
  r = length(truncSample)
  newEstimator[i] = (sum(truncSample) + (numObs-r)*cutoff)/r
}

ggplot(data = data.frame(mean=c(origEstimator, newEstimator), type=c(rep("Complete data", len), rep("Adjusted", len)))) + 
  geom_density(aes(x=mean, color=type))

Captura de pantalla del libro:

Answer 1

Esta es una pregunta interesante.

Primero voy a mostrar que $r>0$. Si $r=0$, entonces no hay datos observados, y la probabilidad de una función no es cóncava, por lo que este problema estadístico no está bien definida. Dado $r>0$, $E(1/r)$ debe ser finito.

Deje $p=\mathrm{Pr}(y_i<c) = 1-\mathrm{exp}(-c/\theta)$. Tenemos $$E(1/r) = \sum_{r=1}^n \binom{n}{r} p^r (1-p)^{n-r}/r = n p (1-p)^{n-1} F(1, 1, 1-n;2, 2;p/(p-1)),$$ donde $F(\cdot)$ es la función hipergeométrica generalizada. $E(\hat{\theta})$ puede ser calculado por el taponamiento $E(1/r)$ en la ecuación en la pregunta.

El siguiente computación numérica muestra que $E(\hat{\theta}) = \theta$, es decir, $\hat\theta$ que considera que el MNAR (falta no al azar) mecanismo es imparcial. Tenga en cuenta que genhypergeo() se utiliza para calcular la función hipergeométrica generalizada numérico de error, pero por encima de la suma puede ser calculado exactamente.

library(hypergeo)

theta_hat = function(theta, n = 100) {
  c = .5*theta # arbitrary c can be used
  p = 1 - exp(-c/theta)
  theta*(1-exp(-c/theta)*(c/theta+1))/(1-exp(-c/theta)) + (n *c) *
    (genhypergeo(U=c(1,1,1-n), L=c(2,2), z=p/(p-1)) * n * p * (1-p)^(n-1)) - c
}

theta = 1:10
plot(theta, sapply(theta, theta_hat), ylab='theta_hat')
abline(a=0, b=1)