clear all; close all;
clf;

% Part (a)

% Generate 200 samples (X_1 to X_200) of i.i.d. zero
% mean, unit variance Gaussian random variables.
% Hint: Use randn and cumsum.

n = 1:200;

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

X = randn(1,200);
S = cumsum(X,2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Do the divide by n part.
S = S./n;

subplot( 4, 1, 1 );
plot( n, S );
xlabel( 'n' );
ylabel( 'Sn' );
title( '7(a) Sample average sequence' );
grid on;

% Part (b)

% Now generate 5000 such sequences.
% Hint: use randn and cumsum (be careful!) again

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

X = randn(5000,200);
S = cumsum(X,2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

S = S./repmat( n, 5000, 1 );

% Part (c)  Strong Law of Large Numbers (this loop will run for a minute)

E_m = zeros(200,1);
for m = 0:199,

% N_m should be the number of rows in S that have an entry whose
% absolute value is > 0.1 in columns m+1 through 200.

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    Nm = 0; % Number of sequences with |S_(i,j)|>0.1 for some j > m
    for i = 1:5000 % Look at each sequence
        for j = (m+1):200 % We check every j > m
            if(abs(S(i,j))>0.1)
                Nm = Nm + 1;
                break;
            end
        end
    end
    E_m(m+1) = Nm/5000;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

end;
subplot( 4, 1, 2 );
plot( E_m );
xlabel( 'n' );
ylabel( 'E_m' );
title( '7(c) Strong Law of Large Numbers' );

%) Part (d) Convergence in Mean Square

% S_squared is the square of the S matrix.  EX = 0.
S_squared = S.^2;

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
M = 1/5000*sum(S_squared);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

subplot( 4, 1, 3 );
plot( n, M );
xlabel( 'n' );
ylabel( 'Mn' );
title( '7(d) Mean square convergence' );

% Part (e) Weak Law of Large Numbers

% Count the number of times |S_i,n| > 0.1 in each column.

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N = sum(abs(S)>0.1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Find E by dividing by 5000.
E = N/5000;

% Plot Pn and En vs. n.  Since EX = 0, Pn is the probability
% that |Sn| > 0.1.  Sn has zero mean and a standard deviation
% sigma_Sn.  Hint: sigma_Sn is a function of n.
% Therefore P( |Sn| > 0.1 ) = 2*Q( 0.1/sigma_Sn ).
% erfc( x ) = 2/sqrt( pi ) * integral from x to inf of exp( -t^2 ) dt
% So Q( x ) = 1/2*erfc( x/sqrt( 2 ) );

% Find sigma_Sn and Pn.  Hint: Pn and En look quite similar.

% WRITE MATLAB CODE HERE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
sigma_Sn = 1./sqrt(n);
P = erfc((0.1./sigma_Sn)/sqrt(2));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

subplot( 4, 1, 4 );
plot( n, E );
hold on;
plot( n, P, 'r--' );
axis( [ 0 200 0 1 ] );
xlabel( 'n' );
ylabel( 'En (solid), Pn (dashed)' );
title( '7(e) Convergence in probability' );

% Produce hardcopy
orient tall
print -dpdf converge_fig