\documentclass[12pt]{article}
\usepackage[spanish]{babel}
\newcommand{\soluciOn}{\sffamily}
\newcommand{\implica}{ $\quad\Rightarrow\quad$ }
\usepackage{Sweave}
\begin{document}
\title{Ejercicio resuelto}
\author{Inferencia}
\maketitle
\raggedright Sea $X$ la variable aleatoria que mide el tiempo, en
meses, entre interrupciones consecutivas de suministro ocurridas en
cierta instalación eléctrica.  Se supone que $X$ sigue una
distribución con densidad
\begin{Schunk}
\begin{Sinput}
> f = function (x, mu) mu^(-1) * exp(-x/mu) * (x>0)
\end{Sinput}
\end{Schunk}
y que las interrupciones son independientes entre sí.  Para estimar la
distribución de $X$ se tomó una muestra aleatoria simple cuyos
resultados fueron los siguientes:
\begin{Schunk}
\begin{Sinput}
> X = c(3.35, 1.30, 1.23, 2.13, 7.70, 11.61, 8.02, 1.44, 3.26, 
+  2.06, 6.06, 0.42, 5.58, 1.40, 6.26, 0.48, 1.46, 7.33, 7.55,
+  1.37, 5.41, 0.23, 3.96, 0.86, 0.85, 1.76, 1.3, 18.36, 2.76,
+  4.22)
\end{Sinput}
\end{Schunk}
\begin{enumerate}
\item Calcula el estimador máximo-verosímil de
  $\mu$ y el valor de la estimación para la muestra obtenida.

  {\soluciOn densidad $f=e^{-x/\mu}/\mu$ \implica verosimilitud $L=e^{-\sum x_i/\mu}/\mu^n$
    \implica $\ln L=-\sum x_i/\mu-n\ln\mu$ \implica $\partial\ln L/\partial\mu=\sum x_i/\mu^2-n/\mu$
    \implica $\bigl\{\quad0=\partial\ln L/\partial\mu$ \implica $\hat\mu_{\rm MV}=\sum x_i/n=\bar x\quad\bigr\}$

    o bien, sabiendo que $X$ es exponencial($1/\mu$) y que
    $\hat\lambda_{\rm MV}=1/\bar x$, por la propiedad de invariancia
    funcional del EMV, $\hat\mu_{\rm MV}=1/\hat\lambda_{\rm MV}=\bar
    x$}

\begin{Schunk}
\begin{Sinput}
> mean(X) # estimación 
\end{Sinput}
\begin{Soutput}
[1] 3.990667
\end{Soutput}
\end{Schunk}

\item Halla la esperanza y la varianza del estimador máximo-verosímil de 
  $\mu$.

  {\soluciOn
    $X_i$ exponencial($1/\mu$), independientes \implica
    $\sum X_i$ gamma($n,1/\mu$) \implica $\hat\mu_{\rm MV}$ gamma($n,n/\mu$) \implica
    esperanza $=$ $n/(n/\mu)$ $=$ $\mu$, \quad varianza $=$ $n/(n/\mu)^2$ $=$ $\mu^2/n$
  
  o bien, sin usar que $X$ es exponencial, 
  $$E(\hat\mu_{\rm MV}) = E\left(\frac{\sum X_i}{n}\right) = \frac{\sum E(X_i)}n
  = \frac{\sum E(X)}n = \frac{n\cdot E(X)}n = E(X)$$
  $$V(\hat\mu_{\rm MV}) = V\left(\frac{\sum X_i}{n^2}\right) = \frac{\sum V(X_i)}{n^2}
  = \frac{\sum V(X)}{n^2} = \frac{n\cdot V(X)}{n^2} = \frac{V(X)}n$$

  \begin{eqnarray*}
    E(X)&=&\int_0^{\infty}x\cdot\frac1\mu\cdot e^{-\frac x\mu}\cdot dx=
  \left<
  \matrix{u=x/\mu&dv=e^{-x/\mu}dx\cr du=dx/\mu&v=-\mu e^{-x/\mu}}
  \right>\\
    &=&\left[-x^2e^{-x/\mu}\right]_{x=0}^{\infty}
        -\int_0^{\infty}-e^{-x/\mu}dx=0+\mu=\mu\\
    &&\\
    E(X^2)&=&\int_0^{\infty}x^2\cdot\frac1\mu\cdot e^{-\frac x\mu}\cdot dx=
  \left<
  \matrix{u=x^2/\mu&dv=e^{-x/\mu}dx\cr du=2x\,dx/\mu&v=-\mu e^{-x/\mu}}
  \right>\\
    &=&\left[-x^2e^{-x/\mu}\right]_{x=0}^{\infty}
        -\int_0^{\infty}-2xe^{-x/\mu}dx=0+2\int_0^{\infty}xe^{-x/\mu}dx\\
    &=&2\mu\int_0^{\infty}x\cdot\frac1\mu\cdot e^{-x/\mu}dx=
        2\mu\cdot E(X) = 2\mu^2\\
    &&\\
    V(X)&=&E(X^2)-\bigl[E(X)\bigr]^2=2\mu^2-\mu^2=\mu^2
  \end{eqnarray*}
  }
\item Obtén un intervalo de confianza con coeficiente $1-\alpha=0.95$
  para la media y la mediana de $X$.

\begin{Schunk}
\begin{Sinput}
> (n = length(X))
\end{Sinput}
\begin{Soutput}
[1] 30
\end{Soutput}
\begin{Sinput}
> alfa = 1 - 0.95
> (a = qgamma(c(alfa/2, 1-alfa/2), n, n))
\end{Sinput}
\begin{Soutput}
[1] 0.6746958 1.3882946
\end{Soutput}
\begin{Sinput}
> (muMV = mean(X))
\end{Sinput}
\begin{Soutput}
[1] 3.990667
\end{Soutput}
\end{Schunk}

  {\soluciOn media de $X$ = $\mu$ ; $\hat\mu_{\rm MV}$ gamma($n,n/\mu$) \implica  $\hat\mu_{\rm MV}/\mu$ gamma($n,n$)
    \implica $0.95=\Pr(0.674695800714031 < \hat\mu_{\rm MV}/\mu < 1.38829458128622)$
    \implica intervalo de confianza para la media $=$ [$\hat\mu_{\rm MV}/1.38829458128622$; $\hat\mu_{\rm MV}/0.674695800714031$]
    = [$2.87451000706883$; $5.91476434630739$]

    mediana de $X$ = $m$ \implica $1/2=F(m)=1-e^{-m/\mu}$ \implica $\ln(1/2)=-\ln2=-m/\mu$ \implica $\hat m_{\rm MV}=\ln2\cdot\hat\mu_{\rm MV}$
    por equivariancia \implica intervalo de confianza para $m$ = $\ln2$ $\cdot$ intervalo de confianza para $\mu$ =
    [$1.99245850689111$; $4.09980223031946$]

    intervalo para $\mu$ mediante \emph{bootstrap} :

\begin{Schunk}
\begin{Sinput}
>     B = 10000
>     ## bootstrap percentil
>     distri = replicate (B,
+     {
+         Xb = sample(X,,TRUE)
+         mean(Xb)
+     })
>     quantile (distri, c(alfa/2, 1-alfa/2))
\end{Sinput}
\begin{Soutput}
    2.5%    97.5% 
2.705333 5.491675 
\end{Soutput}
\begin{Sinput}
>     ## bootstrap básico
>     distri = replicate (B,
+     {
+         Xb = sample(X,,TRUE)
+         mean(Xb) - muMV
+     })
>     muMV - quantile (distri, c(1-alfa/2, alfa/2))
\end{Sinput}
\begin{Soutput}
  97.5%    2.5% 
2.49625 5.27935 
\end{Soutput}
\end{Schunk}

intervalo exacto de longitud óptima:
\begin{Schunk}
\begin{Sinput}
> ## alfa == alfa1 + alfa2
> intervalo = function (alfa1) 
+     muMV / qgamma(c(1-alfa+alfa1,alfa1), n, n)
> longitud = function (alfa1) diff (intervalo (alfa1))
> (o = optimize (longitud, c(0,alfa)))
\end{Sinput}
\begin{Soutput}
$minimum
[1] 0.03802793

$objective
[1] 2.953237
\end{Soutput}
\begin{Sinput}
> intervalo (o$minimum) # intervalo óptimo
\end{Sinput}
\begin{Soutput}
[1] 2.739056 5.692293
\end{Soutput}
\end{Schunk}
  }

% \item Contrasta al nivel de significación $\alpha=0.01$, las hipótesis
% $$H_0: X \equiv \; Exponencial \qquad H_1:X \ne \; Exponencial$$
% \item Contrasta al nivel de significación $\alpha=0.01$, las hipótesis
% $$H0: \mu \le 3 \qquad H_1:\mu > 3$$
\item Calcula la estimación máximo-verosímil de $\mu$, suponiendo que
  los datos de la muestra estuvieran agrupados en las siguientes
  clases:
%$$[0,\, 1),\, [1,\, 1.5),\, [1.5, \,3),\, [3, 4.5),\, [4.5,\,  7.5), (\ge 7.5) $$

 \begin{center}
 \begin{tabular}{c|c|c|c|c|c|c}
 Intervalo & $[0;1)$ & $[1;1.5)$ & $[1.5;3)$ & $[3;4.5)$ & $[4.5;7.5)$ & $\ge 7.5$\\
 \hline Frecuencia & $5$ & $7$ & $4$ & $4$ & $5$ & $5$
 \end{tabular}
 \end{center}
\begin{Schunk}
\begin{Sinput}
> corte = c(0, 1, 1.5, 3, 4.5, 7.5, Inf)
> frecu = c(5, 7, 4, 4, 5, 5)
> L = function (mu) prod (diff (pexp (corte, 1/mu)) ^ frecu)
> lnL = function (mu)
+     sum (log (diff (pexp(corte, 1/mu))) * frecu)
> (muMVa <- optimize (L, c(0, 100), maximum=TRUE) $ maximum)
\end{Sinput}
\begin{Soutput}
[1] 4.053007
\end{Soutput}
\begin{Sinput}
> optimize (lnL, c(0, 100), maximum=TRUE) $ maximum
\end{Sinput}
\begin{Soutput}
[1] 4.053015
\end{Soutput}
\end{Schunk}
 {\soluciOn
salvo por el error de la optimización numérica, el máximo se alcanza
en el mismo lugar, en torno a $\mu=4.05301$
}
\end{enumerate}
\end{document}
