Study with the several resources on Docsity
Earn points by helping other students or get them with a premium plan
Prepare for your exams
Study with the several resources on Docsity
Earn points to download
Earn points by helping other students or get them with a premium plan
Lecture notes on monte carlo simulation, direct sampling techniques, and exponential distribution sampling from the ne-628 course by dr. Nasir m mirza. The notes cover specific algorithms for direct sampling, exponential sampling methods, and matlab programs for monte carlo simulations.
Typology: Lecture notes
1 / 11
This page cannot be seen from the preview
Don't miss anything!
Lecture 5: Nasir M Mirza
6/23/2012 1
Dr. Nasir M Mirza
Lecture 8:
Direct Sampling Techniques
Course: NE-
Lecture 5: Nasir M Mirza
6/23/2012 2
Specific Algorithms: For univariant distributions, there is a general inversion technique that we use now.
1 , 1
, 0 1
0 , 0
( )
u
u u
u
F u
Let y(x) be an increasing function of x. Suppose that u is uniform variable and its cumulative distribution is
There for on (0, 1) range the cumulative distribution function for y is determined by solving the following equation:
F ( y ) , where , isuniform
And y is determined from following:
y F , where , is uniform
Lecture 5: Nasir M Mirza
6/23/2012 4
% MC Sampling for f(y) = lambdaexp(-lambday) % Exponential Distribution sampling program in MATLAB N = 10000 ; % No of histories to be generated max_bins = 500; % Max. no of bins xlambda = 2.0; rand('state', 0) % initialize to zero** for j=1:max_bins % initialize ibin(j) = 0; xmid(j) = 0; ntheory(j)=0; end for i= 1:N % Start Monte Carlo loop yy = rand; % choose random F ff(i) = -log(yy)/xlambda; % calculate function end max_ff = 5; % find max value of func estimated size_bin = max_ff/max_bins; % find size of bins for kk=1: N % put all in bins l_limit=0; u_limit=0; for ii=1: max_bins l_limit = (ii-1)size_bin; u_limit = iisize_bin; if((ff(kk)>l_limit)&(ff(kk)<=u_limit)) jj=ii; end end ibin(jj) = ibin(jj) + 1; % one more for the bin end for k = 1:max_bins xmid(k) = (k - 0.5)size_bin;* % compute mid-point % compute the expected pdf fx = xlambdaexp(-xlambdaxmid(k)); ntheory(k) = Nfxsize_bin;** % no by theory end plot(xmid, ibin,'r+',xmid, ntheory,'b.')
Lecture 5: Nasir M Mirza
6/23/2012 5
0 1 2 3 4 5
0
20
40
60
80
100
120
140
N = 10000 No. of Bins = 500 lambda = 1. f(y) = lambdaexp(-lambday)
f(y)
y
Monte Carlo Method Theory
0 1 2 3 4 5
0
50
100
150
200
250
No. of Bins = 500, Lambda = 2.
Lecture 5: Nasir M Mirza
6/23/2012 7
N = 90000 ; max_bins = 200; a = 10.0; b = 2.0; nm = 2; ss = a./b; rand('state', 0) % initialize the generator to zero for j=1:max_bins % initialize ibin(j) = 0; xmid(j) = 0; ntheory(j)=0; end for i= 1:N % Start Monte Carlo loop yy = rand; % choose random F power = -1./(nm - 1.0); ff(i) = ss.(yy.^power - 1.); % calculate function end max_ff = 5; % find max value of func estimated size_bin = max_ff/max_bins; % find size of bins for kk=1: N l_limit=0; u_limit=0; for ii=1: max_bins l_limit = (ii-1)size_bin; u_limit = iisize_bin; if((ff(kk)>l_limit)&(ff(kk)<=u_limit)) jj=ii; end end ibin(jj) = ibin(jj) + 1; % one more for the bin end for k = 1:max_bins % compute mid-point xmid(k) = (k - 0.5)size_bin; % compute the expected pdf ddd =(a + bxmid(k)).^nm ; fx = (nm - 1.).b.(a.^(nm -1.))./ddd ntheory(k) = 2Nfxsize_bin; % no of point by theory end plot(xmid, ibin,'r+',xmid, ntheory,'b.')
The program in MATLAB
Lecture 5: Nasir M Mirza
6/23/2012 8
0 1 2 3 4 5
0
100
200
300
400
500
600
700
800
900
1000
f(x)
x
Parameters Used: a =10, b = 1, n = 2 N = 10000, No. of bins = 200 Range for x == (0, 5)
Direct Sampling
; 1 , 0.
( 1 ) ( )
n x a bx
n ba f x
Lecture 5: Nasir M Mirza
6/23/2012 10
Transforming back,
Sampling a Gaussian Distribution
2 2 ; ,
2 e rdr e du where u r
r u
x ln 1 cos( 2 2 )
with r u ; 2
and 0 1 1 ; 0 2 1 ;
x r cos ; x r sin
The sampling functions are then given as
This method is called Box-Muller method and gives two independent Gaussian variables. It is slow method.
Lecture 5: Nasir M Mirza
6/23/2012 11
0 1 2 3 4 5
0
100
200
300
400
Results od Sampling
Using a MATLAB program, we have following:
Parameters Used: N = 10000, No. of bins = 200 Range for x == (0, 5)
