%% Rejection Sampling Matlab Implementation:By Shoumik Roychoudhury
%%Sampling a Gaussian G(0,1) distribution using Uniform distribution
%%U(-5,5)
clear all
close all 
clc
%% 
N =1000; % No. of samples drawn
M = 4;
accept = 1;
reject = 1;
a = -5;
b = 5;
%sample from uniform distribution U(-5,5)(Proposal Distribution)
x = a + (b-a)*rand(N,1);
x = x';
u = rand(N,1);% Uniform Distribution u(0,1)
u = u';
x1 = -5:0.1:5;
p = @(x) 1/(sqrt(2*pi)) *exp(-x.^2/2);%Gaussian distribution p(x) = G(0,1) as target distribution
q = @(x) 1/(b-a);%Unifrom distribution q(x) = U(-5,5)

hold on
%rejection sampling
for i = 1:N
   if(u(i)<(p(x(i))/(M*(q(x(i))))))
        samples(accept) = x(i);
        accept = accept + 1;
   else 
        reject = reject + 1;
   end
end
%normalization of histogram
[hs, hx] = hist(samples,x1); 
hs=(hs/length(samples))*(length(x1)/10); 
figure(1)
bar(x1,hs); %plots normalized histogram
hold on
plot(x1,p(x1),'r','linewidth',4)
title('Rejection Sampling');
% rejection percentage
rejection = ((reject-1)/N)*100