2
$\begingroup$

Use MC integration to estimate the probability that X * exp(X) < 2.5, assuming that X ~ Gamma(1.2,3.7)

data mcprobdata;
/* generate the sample */
call streaminit(23891);
count = 0;
do i=1 to 2000;
    p = rand('uniform');
    x = quantile('Gamma', p, 1.2, 3.7);
    count = count + (x*exp(x) < 2.5);
end;

/* calculate the estimate */
t = count/2000;
se = sqrt(t*(1-t)/2000);

run;

Is my t variable correctly calculating this interval probability?

$\endgroup$
  • 2
    $\begingroup$ This probability can be found using other techniques of integration: why don't you compare those results with yours? It might help to know that $x\exp(x)\lt 2.5$ is equivalent to $x\lt 0.958586$ (approximately). $\endgroup$ – whuber Apr 18 '14 at 14:13
1
$\begingroup$

Solution in R for reference.

N <- 10^6
x <- rgamma(N, shape =  1.2, scale = 3.7)
sum(x * exp(x) < 2.5) / N
[1] 0.155834

Check the SAS documentation to see how they parameterize the gamma distribution (shape/scale or shape/rate?).

$\endgroup$
0
$\begingroup$

Yes, it is, as long as x is successfully being drawn as intended.

2000 is quite a low number of MC draws to use.

You can test it by doing a numerical integration in Excel or checking the Gamma cdf by some other means.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.