%% Problem Sheet 2
% Methods of Monte Carlo Simulation

clear all; close all; clc;

%% Exercise 2 (inverse transform)

N = 10^3;
X = zeros(N,1);

for i = 1:N
    U = rand();
    if U <= 0.25
        X(i) = 0;
    else
        X(i) = -(1/3)*log((4/3)*(1-U));
    end
end

figure(1);
cdfplot(X);


%% Exercise 3 (MC integration)

clear all; close all;

N = 10^3;
lambda = [0.5, 1, 2];

% function to integrate and density of exponential distribution
f = @(x) 2 * x .* exp(-2*x);
g = @(x, lambda) lambda * exp(-lambda*x);

for l = lambda
    % generate exponentials using inverse transform
    U = rand(N,1);
    Y = -(1/l)*log(U);
    % plug exponentials into the function to integrate
    Z = f(Y) ./ g(Y, l);
    fprintf('Estimate using lambda = %.1f: %f\n', l, mean(Z));
    fprintf('The standard deviation was: %f\n\n', std(Z)/sqrt(N));
end