Inferencias de Wald: IC y CH

Sean \(X_1,\dots, X_n\) v.a.i.i.d. \(P_\theta, \theta\in \Theta\subset\mathbb{R}; f(x,\theta).\)

Como hemos visto, en situaciones regulares el EMV (\(\hat{\theta}_n\)) es CAN y AE para \(\theta_0\): \[\begin{equation} \sqrt{n}(\hat{\theta}-\theta)\xrightarrow{\mathcal{L}} N(0,\frac{1}{I_1(\theta)})\equiv N(0,V^2(\theta)), \quad (Eq. 1) \end{equation}\] donde \(V^2(\theta)=\frac{1}{I_1(\theta)}.\)

Las inferencias acerca de \(\theta\) basadas en esta distribución asintótica se conocen como inferencias de Wald. Por un lado, para obtener un IC para \(\theta,\) se necesita estimar \(I_1(\theta).\) Tal y como vimos la IF observada, que definimos como \(\widehat{I_1(\theta)}=-\frac{\partial^2}{\partial \theta^2}\log f(X,\theta)\vert_{\theta=\hat{\theta}}\), nos permite estimar esta cantidad, y concretamente: \[\widehat{V^2(\theta)}=\hat{V}^2(\theta)=\frac{1}{\widehat{I_1(\theta)}},\]

es decir, \(V^2(\theta)\) se estima como la inversa de la IF observada, y de donde también se sigue que:

\[\sqrt{n}(\theta-\hat{\theta})\xrightarrow{\mathcal{L}} N(0,\widehat{V^2(\theta)})\text{, o equivalentemente, } \hat\theta\cong N(\theta,\frac{\widehat{V^2(\theta)}}{n})\] NOTA: Cuidado con estas dos expresiones, soléis cometer fallos con la \(n.\)

Por tanto, \(Z_0=\frac{\sqrt n|(\hat{\theta}-\theta)|}{\sqrt{\hat{V}^2(\theta)}}\cong N(0,1).\) Así, fijado \(\alpha\), un IC un Wald \((1-\alpha)\%\) para \(\theta\) queda definido por: \[P_{\theta}( \frac{\sqrt n|(\hat{\theta}-\theta)|}{\sqrt{\hat{V}^2(\theta)}}\leq\lambda_{1-\alpha/2})\simeq 1-\alpha,\] donde \(\lambda_{1-\alpha/2}\) es el cuantil \(1-\alpha/2\) de una normal estandar (\(qnorm(1-\alpha/2)\) en {R}). O equivalentemente, despejando \(\theta\) en la expresión anterior se obtiene el IC de Wald pedido: \[P_{\theta}(\theta\in[\hat{\theta}\pm qnorm(1-\alpha/2)\frac{\sqrt{\hat{V}^2(\theta)}}{\sqrt{n}}])\simeq 1-\alpha\] Gráficamente:

El resultado anterior de la Eq. 1 también es clave para resolver CH para \(\theta\): \[\begin{equation} H_0:\theta=\theta_0 \quad \textit{vs}\quad H_1:\theta\neq \theta_0. \end{equation}\].

Tal y como vimos el estadístico de Wald: \[Q_W=\frac{n(\hat{\theta}-\theta_0)^2}{\widehat{V^2(\theta)}}= n(\hat{\theta}-\theta_0)^2 \widehat{I_1(\theta)}\overset{H_0}{\cong}\chi^2_{1}\] por tratarse del cuadrado de una distribución normal estándar. Esta distribución nos permitirá calcular el p-valor como \(P_{\theta_0}(Q_W>q_{W,obs}),\) donde \(q_{W,obs}\) es el valor observado para el estadístico en la muestra dada.

Inferencias basadas en el estadístico RV: IC y CH

Sean \(X_1,\dots, X_n\) v.a.i.i.d. \(P_\theta, \theta\in \Theta\subset\mathbb{R}; f(x;\theta).\) Dado un nivel de significación \(\alpha\) se quiere resolver el CH para \(\theta\): \[\begin{equation} H_0:\theta=\theta_0 \quad \textit{vs}\quad H_1:\theta\neq \theta_0. \end{equation}\] Según lo estudiado, el estadístico de la RV se define como: \[Q_L(\boldsymbol X)=2[\log L(\hat{\theta},\boldsymbol X)-\log L(\theta_0,\boldsymbol X)]\] y es asintóticamente equivalente a \(Q_w.\) Luego el test de RV de nivel \(\alpha\) para resolver CH para \(\theta\) tiene una región crítica de tamaño \(\alpha\) aproximadamente que viene dada por \((Q_L(X)>\chi^2_{1,1-\alpha})\simeq (Q_L(X)>qchisq(1-\alpha,1)),\) donde \(\chi^2_{1-\alpha,1}\) es el cuantil \(1-\alpha\) de una \(\chi^2_1.\) El p-valor del test vendrá dado por \(P_{\theta_0}(Q_L>q_{L,obs}),\) donde \(q_{L,obs}\) es el valor observado para el estadístico en la muestra dada.

De la relación entre test de nivel \(\alpha\) e IC \((1-\alpha)\), se tiene que el IC basado en el estadístico de la RV (ICRV\((1-\alpha)\)) se corresponde con el conjunto de valores de \(\theta\) para los que el test de la RV de nivel \(\alpha\) no rechaza la \(H_0\) es decir: \[\begin{align} ICRV(1-\alpha)&=\lbrace \theta_0: Q_L(X) \leq qchisq(1-\alpha,1) \rbrace\\ &= \lbrace \theta_0: 2[\log L(\hat{\theta},X)-\log L(\theta_0,X)] \leq qchisq(1-\alpha,1) \rbrace\\ &= \lbrace \theta_0: \log L(\theta_0,X) \geq \log L(\hat{\theta},X)-\frac{qchisq(1-\alpha,1)}{2} \rbrace\\ &= \lbrace \theta_0: \log L(\theta_0,X) \geq d_1 \rbrace \end{align}\]

Ejercicio 1: Modelo de Poisson.

Considere \(X_1,\dots,X_{50}\) m.a.s. de \(\mathcal{P}(\lambda),\) \(f(x,\lambda)=e^{-\lambda}\frac{\lambda^x}{x!},x=0,1,2,\dots.\) Para los datos observados se tiene \(S=\sum_{i=1}^{50} P_{\lambda}(X_i=x_i)=75.\)

  1. Obtener y representar la función de verosimilitud y log verosimilitud.

Verosimilitud: \(L(\lambda,X_1,\dots,X_n)\propto e^{-n\lambda}\frac{\lambda^{\sum_{i=1}^n X_i}}{\prod_{i=1}^n X_i!}\)

Log verosimilitud: \(\log L(\lambda,X_1,\dots,X_n) = -n\lambda +\sum_{i=1}^n X_i \log \lambda\)

  1. Calcular el EMV de forma analítica. Incluye una representación en las gráficas para la verosimilitud y logverosimilitud. Ecuación de verosimilitud: \(\frac{\partial}{\partial \lambda}\log L(\lambda,X_1,\dots,X_n)=-n+\frac{\sum_{i=1}^n X_i}{\lambda}=0,\) de donde se sigue \(\hat{\lambda}=\overline{X}=\frac{\sum_{i=1}^n X_i}{n}=\frac{75}{50}=1.5\)

  2. Obetener el EMV de forma numérica (aproximada) mediante la función {optim}. No siempre es posible calcular EMV de forma analítica y es común recurrir a métodos numéricos (por ejemplo Newton Raphson). En {R} podemos usar la función {optim} que minimiza funciones lineales y no lineales usando Nelder–Mead, quasi-Newton o de descenso de gradiente.

#Leer https://cran.r-project.org/doc/contrib/Optimizacion_Matematica_con_R_Volumen_I.pdf

#optim(par, fn, gr = NULL, ...,
#      method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN",
#                 "Brent"),
#      lower = -Inf, upper = Inf,
#      control = list(), hessian = FALSE)
#By default optim performs minimization, but it will maximize negative functions.


#par: Initial values for the parameters to be optimized over.

#fn: A function to be minimized, with first argument the vector of parameters over which minimization is to take place. It should return a scalar result.

#method: The default method is an implementation of that of Nelder and Mead. Method "Brent" is for one-dimensional problems only

#lower, upper: Bounds bounds in which to search for method "Brent".

#hessian: Logical. Should a numerically differentiated Hessian matrix be returned? Valor de la segunda derivada de fn
  
#Caso uniparamétrico basta con:

NOTA: Estamos trabajando con {mlv}, de forma que la inversa de {minimo$hessian} es un estimador de la varianza asintótica.

  1. Calcular IC de Wald para \(\lambda\) con confianza aproximada del \(95\%.\)

Podemos llegar de dos maneras diferentes.

Distribución asintótica del EMV (TCL): \(\sqrt{n}(\overline{X}-\lambda)\xrightarrow{\mathcal{L}}N(0,\lambda)\), puesto que \(E_\lambda(X_i)=\lambda\) y \(Var_\lambda(X_i)=\lambda \forall i.\)

Distribución asintótica del EMV (Eq. 1): \(\sqrt{n}(\overline{X}-\lambda)\xrightarrow{\mathcal{L}}N(0,\frac{1}{I_1(\lambda)}),\) donde \(I_1(\lambda)=-E_\lambda(\frac{\partial^2}{\partial \lambda^2}\log f(X,\lambda))=-E_\lambda(\frac{-X}{\lambda^2 })=\frac{1}{\lambda}.\)

Para poder hacer inferencias, necesitamos estimar \(I_n(\lambda)\) mediante la IF observada: Calculamos \(\frac{\partial^2}{\partial\lambda^2}\log L(\lambda,X_1,\dots,X_n)=-\frac{\sum_{i=1}^n X_i}{\lambda^2}\), \(\widehat{I_n(\lambda)}=-\frac{\partial^2}{\partial \lambda}\log f(\boldsymbol{X},\lambda)\vert_{\theta=\hat{\lambda}}=\frac{\sum_{i=1}^n X_i}{\lambda^2}\vert_{\lambda=\overline{X}}=\frac{n^2}{\sum_{i=1}^n X_i}=\frac{2500}{75}.\) Nótese que \(\widehat{I_1(\lambda)}=\frac{\widehat{I_n(\lambda)}}{n}=\frac{50}{75}\), luego \(\widehat{V^2(\lambda)}=\frac{1}{\widehat{I_1(\lambda)}}=\frac{75}{50}\).

De donde se sigue: \((\overline{X}-\lambda)\xrightarrow{\mathcal{L}}N(0,\frac{75}{2500}),\) o \(\overline{X}\xrightarrow{\mathcal{L}}N(\lambda,\frac{75}{2500}).\) Nótese que \(se(\hat{\lambda})=\sqrt{\frac{75}{2500}}=\sqrt{\frac{\widehat{V^2(\lambda)}}{50}}\)

Por tanto, el IC \((1-\alpha)\%\) de Wald para \(\lambda\) es: \[(\overline{X}\pm qnorm(0.975)\sqrt{\frac{75}{2500}})\simeq(1.5\pm qnorm(0.975)\sqrt{0.03})\simeq(1.16,1.84)\]

  1. Calcular el estadístico de Wald para contrastar \(H_0:\lambda=1 \quad \textit{vs}\quad H_1:\lambda\neq 1\). ¿Cuál es el p-valor del test?

\(Q_W=\frac{n(\hat{\lambda}-\lambda_0)^2}{\widehat{V^2(\lambda)}}= n(\hat{\lambda}-\lambda_0)^2 \widehat{I_1(\lambda)}=\frac{n(\overline{X}-1)^2}{\overline{X}}\overset{H_0}{\cong}\chi^2_{1}\)

Región crítica: \((Q_W>C_{0.05}),\) donde \(C_{0.05}=qchisq(0.95,1)=3.841,\) pues \(Q_W\sim_{H_0} \chi^2_1.\)

P-valor: \(P_{\lambda=1}(Q_W>q_{W,obs})\cong 1-P_{\lambda=1}(\chi^2_1\leq q_{W,obs})=1-pchisq(q_{W,obs},1),\) donde \(q_{W,obs}=\frac{n(\overline{X}-1)^2}{\overline{X}}=\frac{50(1.5-1)}{1.5^2}\)

  1. Calcular el estadístico test de la RV para contrastar \(H_0:\lambda=1 \quad \textit{vs}\quad H_1:\lambda\neq 1\). ¿Cuál es el p-valor del test? Estadístico RV: \[\begin{align} Q_L&=2[\log L(\hat\lambda,X_1,\dots,X_n)-\log L (\lambda_0,X_1,\dots,X_n)]\\ &=2[\log L(\overline{X},X_1,\dots,X_n)-\log L(1,X_1,\dots,X_n)]\\ &=2[-n\overline{X}+n\overline{X}\log\overline{X}+n]\overset{H_0}{\cong}\chi^2_{1} \end{align}\]

P-valor: \(P_{\lambda=1}(Q_L>q_{L,obs})\cong 1-P_{\lambda=1}(\chi^2_1\leq q_{L,obs})=1-pchisq(q_{L,obs},1),\) donde \(q_{L,obs}=2[-75+75\log(1.5)+50]\)

  1. Calcular IC basado en la RV para \(\lambda\) con confianza aproximada del \(95\%.\) \[\begin{align} ICRV(0.95)&=\lbrace \lambda_0: 2[\log L (\overline{X},\boldsymbol{X})-\log L (\lambda_0,\boldsymbol{X})]\leq qchisq(0.95,1) \rbrace\\ &=\lbrace \lambda_0: \log L(\lambda_0,\boldsymbol{X})\geq\log L(\overline{X},\boldsymbol{X})-\frac{qchisq(0.95,1)}{2}\rbrace\\ &=\lbrace \lambda_0: \log L(\lambda_0,\boldsymbol{X})\geq d_1\rbrace \end{align}\]

  2. Calcular IC de Wald para \(g(\lambda)=P_\lambda(X=0)=e^{-\lambda}\) con confianza aproximada del \(95\%.\) Por un lado, el EMV es invariante para cualquier función \(g\) (con inversa), se tiene que \(\widehat{g(\lambda)}=g(\hat{\lambda})=g(\overline{X})=e^{-\overline{X}}.\) Por otro lado, el Delta Método permite a partir de la distribución asintótica del EMV calcular la distribución de cualquier función \(g\) (derivable) del EMV como: \[\sqrt{n}(g(\hat\lambda)-g(\lambda))\xrightarrow{\mathcal{L}}N(0,\widehat{V^2(\lambda)}(g^{'}(\hat\lambda))^2,\] donde \(g(\hat\lambda)=e^{-\overline{X}},\) \(g(\lambda)=e^{-\lambda},\) \(\widehat{V^2(\lambda)}=\overline{X}\) y \((g^{'}(\hat\lambda))^2=(-e^{-\hat\lambda})^2=e^{-2\overline{X}}.\) De donde se sigue: \[\frac{\sqrt{n}(e^{-\overline{X}}-e^{-\lambda})}{\sqrt{\overline{X}e^{-2\overline{X}}}}\cong N(0,1).\] El IC \((1-\alpha)\%\) de Wald viene dado por: \(P_{\lambda}(\frac{\sqrt{n}\vert e^{-\overline{X}}-e^{-\lambda} \vert}{\sqrt{\overline{X}e^{-2\overline{X}}}}\leq qnorm(1-\alpha/2))\simeq 1-\alpha.\) De forma que despejando \(\lambda\) de la expresión anterior se llega al resultado:\[(e^{-\overline{X}}\pm qnorm(1-\alpha/2)\sqrt{\frac{\overline{X}e^{-2\overline{X}}}{n}})\]

  3. Calcular el IC basado en la RV para \(g(\lambda)=P_\lambda(X=0)\) con confianza aproximada del \(95\%.\)

Los IC basados en el estadístico de la RV son invariantes frente a transformaciones. Para calcular un IC para una función del parámetro bastará con aplicar dicha función a los extremos del intervalo RV obtenidos para el parámetro. En este caso el IC pedido es: \[(e^{-1.19},e^{-1.86})\simeq (0.16,0.30)\]

Ejercicio 2: Modelo de Exponencial.

Considere \(X_1,\dots,X_{30}\) m.a.s. de \(\mathcal{Exp}(\lambda),\) \(f(x,\lambda)=\lambda e^{-\lambda x},x>0\), \(\lambda>0.\) Para los datos observados se tiene \(\overline{X}=0.9\).

  1. Obtener y representar la función de verosimilitud y log verosimilitud.
  2. Calcular el EMV de forma analítica.
  3. Obetener el EMV de forma numérica (aproximada) mediante la función {optim}.
  4. Calcular IC de Wald para \(\lambda\) con confianza aproximada del \(95\%.\)
  5. Calcular el estadístico test de Wald para contrastar \(H_0:\lambda=1 \quad \textit{vs}\quad H_1:\lambda\neq 1\). ¿Cuál es el p-valor del test?
  6. Calcular el estadístico test de la RV para contrastar \(H_0:\lambda=1 \quad \textit{vs}\quad H_1:\lambda\neq 1\). ¿Cuál es el p-valor del test?
  7. Calcular IC basado en la RV para \(\lambda\) con confianza aproximada del \(95\%.\)

Ejercicio 3:

Se seleccionan aleatoriamente 1000 muestras de suelo que se combinan en lotes de tamaño 10 muestras. Los lotes son examinados para evaluar la presencia de una toxina, es decir, las muestras individuales de suelo no se examinan. Concretamente, en 100 lotes, de 10 muestras de suelo cada uno, se detecto la toxina en 12 lotes.

  1. Calcular el EMV para \(p=\text{probabilidad de que la toxina esté presente en una muestra de suelo individual}.\)
  2. Detaerminar la distribución asintótica del EMV de \(p\) y calcular IC aproximados de Wald para \(p\) con confianza 0.95.
  3. Calcular ICRV con confianza aproximada de 0.95 para \(p\).