Page 166

```Matlab: Coin Toss Random Walk

nSteps = 30;
nTrials = 200000;
H=zeros(1,nTrials);
for j=1:nTrials
% do all nSteps coin tosses
T = 2*(rand(1,nSteps) < 0.5 ) -1;
H(j) = sum(T);
end
hist(H)
```

Page 172

```Matlab: 2-D Brownian Motion

sig=2; % step size std dev
nSteps = 2000; % # steps in each walk
nSteps = 2*fix(nSteps/2); % assure even
% box-muller for N(0,1) random variables
u1=rand(nSteps/2,1); %U[0,1) vector nSteps/2 long
u2=rand(nSteps/2,1);
x1 = cos(2*pi*u1).*sqrt(-2*log(u2)); %N(0,1) vector
x2 = sin(2*pi*u1).*sqrt(-2*log(u2));
% concatenate x1 and x2 into h
h(1:nSteps/2) = x1;
h(1+nSteps/2:nSteps) = x2;
x(1)=0; y(1)=0; %start at the origin
for i=2:nSteps
theta = 2*pi*rand; % random direction
x(i)=x(i-1)+sig*h(i)*cos(theta);
y(i)=y(i-1)+sig*h(i)*sin(theta);
end
%typical walk
figure(1)
subplot(2,2,1);
plot(x,y);
title('Typical walk for nSteps steps');
xlabel('X position');
ylabel('Y position');
% do the walk from subStart to subEnd
sub1Start = 500;
sub1End = 1500;
k = 1;
for i=sub1Start:sub1End
xsub1(k) = x(i);
ysub1(k) = y(i);
k = k+1;
end
figure(1)
subplot(2,2,2);
plot(xsub1,ysub1);
sub2Start = 1000;
sub2End = 1300;
k = 1;
for i=sub2Start:sub2End
xsub2(k) = x(i);
ysub2(k) = y(i);
k = k+1;
end
figure(1)
subplot(2,2,3);
plot(xsub2,ysub2);
nWalks = 400; % # walks for density
for j=1:nWalks
u1=rand(nSteps/2,1);
u2=rand(nSteps/2,1);
x1 = cos(2*pi*u1).*sqrt(-2*log(u2));
x2 = sin(2*pi*u1).*sqrt(-2*log(u2));
h(1:nSteps/2) = x1;
h(1+nSteps/2:nSteps) = x2;
x(1)=0; y(1)=0;
for i=2:nSteps
theta = 2*pi*rand;
x(i)=x(i-1)+sig*h(i)*cos(theta);
y(i)=y(i-1)+sig*h(i)*sin(theta);
end
Endx(j) = x(nSteps);
Endy(j) = y(nSteps);
end
%density plot
figure(1);
subplot(2,2,4);
scatter(Endx,Endy,3);
title('Density Plot for nSteps steps');
xlabel('X end');
ylabel('Y end');
```

Page 178

```Matlab: Option Price via Simulation

nDays = 90;
dt = 1/365.0; % use time in years
T = nDays*dt; % expiry in years
S0=20; % starting stock price
X=25; % strike price
r = .031; % risk-free rate
sigma = .6; % volatility
expTerm = r*dt;
stddev = sigma*sqrt(dt);
nTrials = 10000;
value=0;
for j=1:nTrials
n = randn(1,nDays); % vector of N(0,1) random numbers
S = S0;
for i=1:nDays
dS = S*(expTerm + stddev*n(i));
S = S + dS;
end
S90(j) = S; % to histogram
value = value + max(S-X,0);
end
value = value/nTrials;  % expected value
Price = exp(-r*T)*value % discount back
hist(S90,0:.5:65)
```

Page 188

```Matlab: Kelly Problem Evaluation via Simulation

function fitnesses = kellyEval(popSize,chromPop);
popSize = 1; % remove for use in GA
chromPop=0.05; % remove for use in GA
F=linspace(1,1,popSize); % initial fortune
N=1000; % iterations of the game
h=linspace(1,1,N); % histogram, remove for GA
for i=2:N
w = F.*chromPop; %term-by-term wager
% play 60/40 game
r = (rand(1,popSize) < 0.6); % 1 if <0.6 else 0
r = 2*r-1; % +1 if <0.6 else -1
w = r.*w; % +bet if <0.6, else -bet
h(i)=F; % remove for use in GA
F = F+w; % new fortune
end
fitnesses = exp((1/N)*log(F));
```

Page 195

```Matlab: Kelly Problem Solution via Genetic Algorithms

% a "chromosome" here is a fraction between 0 and 1
% a population is a vector of chromosomes
popSize=32;
chromPop= rand(1,popSize); %initalize randomly
% space for their fitnesses = exp(growth rate)
fitnesses=linspace(0,0,popSize);
newPop=linspace(0,0,popSize);
upper = linspace(1,1,popSize); % keep fractions < 1
lower = linspace(0,0,popSize); % keep fractions > 0
bp = linspace(0,0,popSize); % roulette wheel breakpts
pm = 0.1; % mutation probability
orelax = 1.4; % over relaxation factor for mating
epsilon=0.07 % mutation range
bestFitness=0;
%chromPop % for debugging
for g=1:100  % 100 generations
% mutate each w/prob pm
m = (rand(1,popSize)<pm); % mutation indices
% mutate each chromoPop(i) within epsilon
a = max(lower,chromPop-epsilon);
b = min(upper,chromPop+epsilon);
r = m.*(a+(b-a).*rand(1,popSize)); % new value
chromPop = (1-m).*chromPop; %null out old value
chromPop = chromPop+r; %replace by new value
%chromPop % for debugging
% mate each chromosome via a permutation of chromPop
permPop = chromPop;
for i=1:popSize-1
k = i+floor((popSize+1-i)*rand); %random in remaining
swap = permPop(i);
permPop(i) = permPop(k);
permPop(k)=swap;
end
%permPop % for debugging
newPop = chromPop +...
orelax*(permPop-chromPop).*rand(1,popSize);
% stay within 0 and 1
newPop = max(lower,newPop);
newPop = min(upper,newPop);
%newPop % for debugging
% evaluate newPop
fitnesses = kellyEval(popSize,newPop);
oldBest = bestFitness;
bestFitness = max(oldBest,max(fitnesses));
if bestFitness > oldBest
m = fitnesses>=bestFitness; % where is it
solution = m.*newPop;
end
% roulette wheel select, make up the wheel
sum = 0;
for i=1:popSize
sum = sum + fitnesses(i);
bp(i) = sum;
end
%bp % for debugging
% now pick popSize times
for i=1:popSize
U=rand; k=1;
while( U>= bp(k))
k=k+1;
end; %k is selected
chromPop(i) = newPop(k);
end
%chromPop % for debugging
end
bestFitness
solution
```