Problem 1

we want to find

\(V=\int_0^\infty \log(1+x^2)e^{-x}dx\)

use simple simulation and antithetic variables. Use simulation to find their respective variances.

\(V=E\{\log(1+X^2)\}\), where \(X\sim Exp(1)\), so

u <- runif(1e4)
x <- -log(u)
c(mean(log(1+x^2)), var(log(1+x^2)))
## [1] 0.6906371 0.5938963
y <- (-log(u)-log(1-u))/2
c(mean(log(1+y^2)), var(log(1+y^2)))
## [1] 0.6863131 0.1489408

Problem 2

Let \(X_1,..,X_{10}\sim Pois(1)\) and independent. Find \(P(\max\{X_i\}>10)\)

  1. analytically (you can use the R command ppois)

\[ \begin{aligned} &P(\max\{X_i\}>10) =1-P(\max\{X_i\}\le10) = \\ &1-P(X_1\le30,..,X_{10}\le10) \\ &1-\prod_{i=1}^{10} P(X_i\le 10) = \\ & 1-P(X_1\le 10)^{10} \\ \end{aligned} \]

1-ppois(10, 1)^10
## [1] 1.004777e-07
  1. via importance sampling

Let’s use \(Y\sim Pois(\lambda)\), and play around a bit with the \(\lambda\). We find

\[ \begin{aligned} &F_{\max\{Y_i\}}(k) =P(\max\{Y_i\}\le k) = P(Y_1\le k)^{10} = ppois(k, \lambda)^{10}\\ &f_{\max\{Y_i\}}(k) = P(Y_1\le k)^{10}-P(Y_1\le k-1)^{10} ; k>0 \end{aligned} \]

\[w =\frac{f_{\max\{X_i\}}(k)}{f_{\max\{Y_i\}}(k)}\]

hw12p2 <- function(lambda=1, B=1e5) {
     y=matrix(rpois(10*B, lambda), ncol=10)
     ymax = apply(y, 1, max)
     I=c(1:B)[ymax>10]
     k=ymax[I]
     w = (ppois(k, 1)^10-ppois(k-1, 1)^10)/
         (ppois(k, lambda)^10-ppois(k-1, lambda)^10)
     c(length(I)/B, sum(w)/B)
}
hw12p2(6)
## [1] 3.54780e-01 1.00513e-07
hw12p2(7)
## [1] 6.431000e-01 1.006184e-07
hw12p2(6.5)
## [1] 4.978900e-01 1.004634e-07

so here we have the event about 50% of the times, and so we use this result, which is very close to the true answer.