% Skript zur Lösung von Aufgabe 3 - Stochastische Prozesse
clc, clear all

% Setze den 'random seed', um bei jedem Aufruf des Skripts diesselben Ergebnisse zu erhalten
rng(1);

%% ---------------------------------------------------------------
%% Teil a)
%% ---------------------------------------------------------------

% Anzahl Pfade
m = 5;

% Laufzeit des Prozesses
N = 100;

% Wertebereich für p
p_range = [0.1, 0.25, 0.5];
n_p = length(p_range);

figure(1), clf
for j = 1:n_p
	p = p_range(j);

	% Erzeuge N * m gleichverteilte Zufallsvariable auf [0,1]
	xi = rand(N,m);

	% Initialisiere Matrix für Zufallsvariable X
	X = ones(N,m);

	% Initialisiere Matrix für Pfade
	M = zeros(N+1,m);

	% Generiere aus xi die Zufallsvariable X mit Verteilung P(X_i = 1) = p, P(X_i = -1) = 1-p
	X(xi > p) = -1;

	M(2:end,:) = cumsum(X);
	subplot(1,n_p,j)
	plot(M)
	title(sprintf('Pfade von M (p = %3.3f)',p))
	xlabel('j')
	ylabel('M_j')

end

%% ---------------------------------------------------------------
%% Teil b) 
%% ---------------------------------------------------------------

% Anzahl Pfade
m = 10000;

% Laufzeit des Prozesses
N = 100;

figure(2), clf

for j = 1:n_p
	p = p_range(j);

	% Erzeuge N * m gleichverteilte Zufallsvariable auf [0,1]
	xi = rand(N,m);

	% Initialisiere Matrix für Zufallsvariable X
	X = ones(N,m);

	% Initialisiere Matrix für Pfade
	M = zeros(N+1,m);

	% Generiere aus xi die Zufallsvariable X mit Verteilung P(X_i = 1) = p, P(X_i = -1) = 1-p
	X(xi > p) = -1;

	M(2:end,:) = cumsum(X);
	subplot(1,n_p,j)
	M_100 = M(end,:);
	histogram(M_100, -N:2:N, 'Normalization', 'pdf')
	title(sprintf('Histogramm von M_{%i} (p = %3.3f)',N,p))

	% Gegen welche Verteilung konvergiert M_N?
	normpdf = @(x,mu,sigma_squared)  1./ sqrt(2*pi*sigma_squared) * exp(- ((x-mu).^2 ) /(2*sigma_squared) );
	mu_i = 2*p - 1;
	mu = N * mu_i;
	sigma_i_squared = 4 * p * (1-p);
	sigma_squared = N * sigma_i_squared;
	hold on
	x = linspace(-N,N,2*N);
	fx = normpdf(x,mu,sigma_squared);
	plot(x,fx,'r-')

end % for j = 1:n_p


%% ---------------------------------------------------------------
%% Teil c) 
%% ---------------------------------------------------------------

% Anzahl Pfade
m = 10000;

% Laufzeit des Prozesses
N = 100;

figure(3), clf

% Erzeuge N * m standard-normalverteilte Zufallsvariable
X = randn(N,m);

% Initialisiere Matrix für Pfade
B = zeros(N+1,m);

B(2:end,:) = cumsum(X);

B_100 = B(end,:);
histogram(B_100, -N:2:N, 'Normalization', 'pdf')
title(sprintf('Histogramm von B_{%i}', N))

% Gegen welche Verteilung konvergiert B_{100}?
normpdf = @(x,mu,sigma_squared)  1./ sqrt(2*pi*sigma_squared) * exp(- ((x-mu).^2 ) /(2*sigma_squared) );
mu_i = 0;
mu = N * mu_i;
sigma_i_squared = 1;
sigma_squared = N * sigma_i_squared;
hold on
x = linspace(-N,N,2*N);
fx = normpdf(x,mu,sigma_squared);
plot(x,fx,'r-')
