clear all
close all
clc

%% Part I 

p=0.25; % take 0 < p < 1/2
N=1000;
M=1000;
X=2*(rand(N,M)<p)-1;
S=cumsum(X);

lambda=1.1;
Y=lambda.^S;
figure,plot(Y);

% we observe that when lambda<1 (and 0<p<1/2), the process Y diverges
% for lambda>1 (and 0<p<1/2), the process Y converges to 0 at infinity,
% but may go for large excursions before that, especially for large
% values of lambda

%% Part II

lambda=(1-p)/p;  % value of lambda for which Y is a martingale
Y=lambda.^S;
a=2;             % threshold for S

T=ones(1,M); % stopping time T
A=ones(1,M);
for n=1:N-1
    A=A.*(S(n,:)<a).*(S(n,:)>-a); % keeps track when S remains between -a and a
    T=T+A; % increments the time T as long as A=1
end
ST=zeros(M,1);
for m=1:M
    ST(m)=S(T(m),m); % value of S at time T
end
P=sum(ST==a)/M  % probability that S(T)=a
                % theoretical value: P=1/(1+lambda^a) with lambda=(1-p)/p

Tp=ones(1,M); % same thing as above, but this time for T'
Ap=ones(1,M);
for n=1:N-1
    Ap=Ap.*(S(n,:)<a);
    Tp=Tp+Ap;
end
STp=zeros(M,1);
for m=1:M
    STp(m)=S(Tp(m),m);
end
Pp=sum(STp==a)/M  % probability that S(T')=a
                  % theoretical value: P=1/lambda^a with lambda=(1-p)/p