\documentstyle[12pt,fullpage,times]{article}

\begin{document}

\frenchspacing

\def\Var{\mbox{Var}}
\def\Cov{\mbox{Cov}}
\let\ul\underline

\begin{center}
\large\bf
Finger exercises for Week 13 --- To be discussed Monday Nov 25, 2002
\end{center}

\begin{enumerate}

% \item In the handout discussing elementary integration methods a
%    reasonably efficient ``automatic'' implementation of the trapezoid
%    rule in S was provided ``trap-int.s'').  What is the corresponding
%    implementation of Simpson's rule?  Give several examples showing
%    that it works as expected.

\item We continue with the data $Y$
\begin{verbatim}
1.2 2.4 1.3 1.3 0.0 1.0 1.8 0.8 4.6 1.4 ,
\end{verbatim}
modeled as iid observations $y$ from a Cauchy likelihood
\[
   f(y|\theta,s) 
   = \frac{1}{\pi s} \frac{1}{1+\left(\frac{y-\theta}{s}\right)^2}
\]
with $s$ fixed.  Assume now that $s=0.1$ and that $\theta$ has a prior
distribution, $\theta$ $\sim$ $N(1,1)$.

\begin{enumerate}

\item Graph the log-posterior density $\log p(\theta|Y)$ (or something
  log-proportional to it) for the Cauchy data above.

\item Compute an approximate 95\% credible interval for $\theta$,
using the posterior mode $\tilde\theta$ and
$-\frac{d^2}{d\theta^2}\log p(\tilde\theta|Y)$, using whatever tools
are convenient (optimization functions in R or Splus, numerical or
symbolic differentiation, etc.)
\label{approx}

\item Use the R/Splus function {\tt integrate} to calculate the coverage
  probability for the interval calculated in part \ref{approx}.  In
  other words, calculate $P(\theta\in[a,b]|Y)$, where $[a,b]$ is the
  interval you calculated in part \ref{approx}.

\item Try calculating $E[\theta|Y]$ and $\sqrt{\Var(\theta|Y)}$ using
  SPLUS's {\tt integrate} function.  Is the approximate 95\% credible
  interval based on $E[\theta|Y]$ and $\sqrt{\Var(\theta|Y)}$ much
  different from the interval in part \ref{approx}?  Is the coverage
  probability much different?

\item To explore what happens to numerical integrals when the
likelihood becomes more and more peaked, recalculate $E[\theta|Y]$ and 
$\sqrt{\Var(\theta|Y)}$, when the dataset $Y$ is replaced with 
\begin{itemize}
\item The same dataset, with each observation repeated 2 times.
\item The same dataset, with each observation repeated 3 times.
\item The same dataset, with each observation repeated 5 times.
\item The same dataset, with each observation repeated 10 times.
\item The same dataset, with each observation repeated 20 times.
\item The same dataset, with each observation repeated 30 times.
\end{itemize}
and describe/explain your findings.

\end{enumerate}

% \item  Verify the following methods for generating common random
%    variables:
% \begin{itemize}
% \item $U \sim Unif[0,1]$ imples $X = \log(U)
% \sim Exp(1)$.
% \item $U \sim Unif[0,1]$ imples  $X=\tan[\pi (U-1/2)]$ $\sim$ Cauchy.
% \end{itemize}

\item The inverse-CDF method says that if $U \sim Unif[0,1]$ then
  $X=F^{-1}(U)$ will have CDF $P[X\leq x]=F(x)$.  If $X$ is a discrete
  random variable taking on the values $0,1,2,\ldots$ then $F(x)$ is a
  step function: $F(x)\equiv F(k)$ for all $x\in[k,k+1)$.  This
  suggests a ``search'' method of generating $X$: First generate
  $U\sim Unif[0,1]$; then find the $k$ such that $F(k)\leq U< F(k+1)$,
  and return $X=k$.
\begin{enumerate}
\item Show that this ``search'' method produces a random draw from
  $X\sim F(x)$.
\item Implement the method in an Splus or R function to generate
  random draws from a Poisson($\lambda$) random variable.  Show that
  it works.
\item Discuss at least two drawbacks of the method.
\end{enumerate}

\item In the rejection methods, on average, how many pairs $(Y,U)$
  have to be drawn, before $Y$ is accepted as an observation from the
  ``target'' density $f(x)$?  What does this say about what the
  relationship between the ``target'' density $f(x)$ and the
  ``convenient'' density $g(x)$ should be?

\item Show that you can generate $X \sim N(0,1)$ by applying the
  rejection method to Cauchy observations $Y=\tan(\pi(U-1/2))$.  What
  is (the \ul{\em best}) $a$?  Show (analytically) that it does
  \ul{\em not} work to generate $X \sim Cauchy$ by applying the
  rejection method to $N(0,1)$ observations $Y$.  {\em What goes
    wrong???} 

\item  Can you develop a rejection method for discrete random variables
  (the method for continuous rv's was discussed in class).  Does the
  ``convenient'' density $g(x)$ have to be for a discrete random
  variable?

\item In the notes and in class we discussed two versions of
  importance sampling, one involving
\[
        \tilde I_n = \frac1n\sum_1^n  h(Y_i)f(Y_i)/g(Y_i) \rightarrow
        E_f[h(X)] \mbox{ ,~~~ as $n\rightarrow\infty$}
\]
and the other involving $\hat I_n = \tilde I_n/\tilde O_n$ where
\[
        \tilde O_n = \frac1n\sum_1^n 1\cdot f(Y_i)/g(Y_i) \rightarrow
        1  \mbox{ ,~~~ as $n\rightarrow\infty$}
\]
Apply the delta method to the function $g(x,y) = x/y$, where
$x=\frac1n\tilde I_n$ and $y=\frac1n\tilde O_n$, to show that 
\[
           \Var(\hat I_n) \approx \Var(\tilde I_n)
                          - 2E_f[h(X)] \Cov(\tilde I_n,\tilde O_n)
                          + (E_f[h(X)])^2 \Var(\tilde O_n)
\]
Under what conditions could we expect $\Var(\hat I_n) \leq \Var(\tilde
I_n)$?



\end{enumerate}

\end{document}
