% Log-optimal portfolio
% jp vert, June 15, 2006
% Originally taken from Steve Boyd's course:
% http://www.stanford.edu/class/ee364/matlab/log_opt_invest.m


clear all

P = [3.5000    1.1100    1.1100    1.0400    1.0100;
     0.5000    0.9700    0.9800    1.0500    1.0100;
     0.5000    0.9900    0.9900    0.9900    1.0100;
     0.5000    1.0500    1.0600    0.9900    1.0100;
     0.5000    1.1600    0.9900    1.0700    1.0100;
     0.5000    0.9900    0.9900    1.0600    1.0100;
     0.5000    0.9200    1.0800    0.9900    1.0100;
     0.5000    1.1300    1.1000    0.9900    1.0100;
     0.5000    0.9300    0.9500    1.0400    1.0100;
     3.5000    0.9900    0.9700    0.9800    1.0100];

[m,n] = size(P);


% Solve the log-optimal portfolio problem
cvx_begin
    variable x_opt(n)
    maximize(geomean(P*x_opt))
    sum(x_opt) == 1
    x_opt >= 0
cvx_end

% Display some results:
disp('The log-optimal portforlio is:');
disp(x_opt);

x_unif = ones(n,1)/n;
R_opt = sum(log(P*x_opt))/m;
R_unif = sum(log(P*x_unif))/m;
disp(sprintf('The long-term growth rate with the log-optimal portfolio is %.3g percent',R_opt*100))
disp(sprintf('The long-term growth rate with the uniform portfolio is %.3g percent',R_unif*100))

% Generate random event sequences
rand('state',10);
N = 10;  % number of random trajectories 
T = 200; % time horizon
w_opt = []; w_unif = [];
for i = 1:N
    events = ceil(rand(1,T)*m);
    P_event = P(events,:);
    w_opt = [w_opt [1; cumprod(P_event*x_opt)]];
    w_unif = [w_unif [1; cumprod(P_event*x_unif)]];
end

% Plot wealth versus time
figure
semilogy(w_opt,'g')
hold on
semilogy(w_unif,'r--')
grid
axis tight
xlabel('time')
ylabel('wealth')
