2

Given $I=\int_{3}^{+\infty}\dfrac{1}{\sqrt{2\pi}}e^{-x^2/2} \, dx$, I want to use a Control variate function to reduce the variance of a simple Monte Carlo simulation to compute this integral. In order to change this integral domain to a real one, I've made a variable substitution of the type $y=\dfrac{1}{x}$, which results in $I=\int_0^{1/3}\dfrac{1}{\sqrt{2\pi}}y^{-2}e^{-1/(2y^2)}$ if I didn't do any error.
For the simple Monte Carlo this is what i did in R:

hx<-function(y){1/sqrt(2*pi)*y^(-2)*exp(-1/(2*y^2))} # Variable substitution
u<-runif(m,0,1/3)
I_MC<-(1/3-0)*mean(hx(u)); I_MC
# Variance
V_MC<-(1/3-0)^2*var(hx(u))/m;V_MC

And I've used as control variate the function $f(x)=x$:

fx<-function(x){x}
m=10000
u<-runif(m,0,1/3)
c_star<- -cov(hx(u),fx(u))/var(fx(u));c_star
Efx<-mean(fx(u));Efx
I_MC_CV<-(1/3)*mean(hx(u)+c_star*(fx(u)-Efx));I_MC_CV
# Variance
V_MC_CV<-(1/3)^2*var(hx(u)-c_star*(fx(u)-Efx))/m;V_MC_CV

I can't understand what I'm doing wrong here, as the result for both integrals is the correct one, but the variance reduction I'm getting is: -132%, which means that it increased variance instead.
Can someone help me out finding the error? Thank you!

tomas-silveira
  • 141
  • 1
    See https://stats.stackexchange.com/a/146983/919 for an account of the conditions you need to check. In particular, consider the correlation between $f$ and the integrand you have computed. – whuber Nov 06 '20 at 18:25
  • Thank you, I think i check all the conditions, the control variate is correlated with my integrand (f), and the expressions seem to be alright... Can you help me out finding what I'm missing? – tomas-silveira Nov 06 '20 at 19:02

0 Answers0