[Case Study]Reject Sampling using Matlab

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

In mathematics, rejection sampling is a basic technique used to generate observations from a distribution. It is also commonly called the acceptance-rejection method or "accept-reject algorithm" and is a type of Monte Carlo method. Rejection sampling is based on the observation that to sample a random variable one can sample uniformly from the region under the graph of its density function.


Description
To visualize the motivation behind rejection sampling, imagine graphing the density function of a random variable onto a large rectangular board and throwing darts at it. Assume that the darts are uniformly distributed around the board. Now take off (reject) all of the darts that are outside the area under the curve. The remaining darts will be distributed uniformly within the area under the curve, and the x-positions of these darts will be distributed according to the random variable's density. This is because there is the most room for the darts to land where the curve is highest and thus the probability density is greatest.

The visualization as just described is equivalent to a particular form of rejection sampling where the proposal distribution is uniform (hence its graph is a rectangle). The general form of rejection sampling assumes that the board is not necessarily rectangular but is shaped according to some distribution that we know how to sample from (e.g. using inversion sampling), and which is at least as high at every point as the distribution we want to sample from, so that the former completely encloses the latter. (Otherwise, there will be parts of the curved area we want to sample from that can never be reached.) Rejection sampling works as follows:

Sample a point (an x-position) from the proposal distribution.

• Draw a vertical line at this x-position, up to the curve of the proposal distribution.
• Sample uniformly along this line (i.e. uniformly from 0 to the value of the proposal distribution (maximum of the probability density function)). If the sampled value is greater than the value of the desired distribution at this vertical line, return to step 1.
• Note also that this algorithm can be used to sample from the area under any curve, regardless of whether the function integrates to 1. In fact, scaling a function by a constant has no effect on the sampled x-positions. This means that the algorithm can be used to sample from a distribution whose probability density function is only known up to a constant (i.e. whose normalizing constant is unknown), which is common in computational statistics.
  

Now we perform the rejection sampling algorithm as follows:
1. Sample x from g(x) and u from U(0,1).
2. Check whether u<f(x)/Mg(x) or not.
3. If this holds, accept x as a realization of f(x).
4. If not, reject the value of x and repeat the sampling step.
Let us choose the uniform distribution as g(x). We see that we can choose M=8/5 since f (x) < 8/5 g(x) . The following code implements the rejection sampling:

1. M=8/5;N=10000;
   
2. f=@(x) 4/5*((x>0 & x<1)+(x>3/8&x<5/8));
   
3. g=@(x) (x<1 & x>0);
   
4. xsamples=[];
   
5. for i=1:N
   
6. u=rand(1);
   
7. x=rand(1);
   
8. if u<f(x)/(M*g(x))
   
9. xsamples=[xsamples x];
   
10. end
    
11. end
    
12. hist(xsamples);
    
13. We can also implement the same algorithm using vectorization. We can sample N number of us and xs and accept
    
14. or reject them in one shot.
    
15. u=rand(N,1);x=rand(N,1);
    
16. accept=u<f(x)./(M*g(x));
    
17. xsamples=x(accept);
    
18. hist(xsamples)

复制代码

Reference
Numerical Methods using Matlab Abhishek Gupta Apress

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：Case study Sampling MATLAB Reject matla motivation generate called method