===== Cvičení 10 =====

[[https://cw.felk.cvut.cz/lib/exe/fetch.php/courses/a4b33opt/cviceni/markowitz.pdf?id=courses%3Aa4b33opt%3Acviceni%3Astart&cache=cache|Zadání (Markowitzův model)]]

  * bacha na to, zadání se změnilo, stáhněte si ze stránek předmětu [[http://cw.felk.cvut.cz/lib/exe/fetch.php/courses/a4b33opt/cviceni/markowitz_24nov09.pdf?id=courses%3Aa4b33opt%3Acviceni%3Astart&cache=cache|aktuální verzi]]

Ke spolecnemu vypracovani:
<code matlab>
clc;
clear;
close all;
%% Markowitzùv model
% 1 Nacteni dat
load akcie.mat;
% init
[T,I] = size(AkcieOdhad);
% vynosy akcii za cas
r = zeros(I,T);
vyberR = zeros(I,1);
p = 1/(T-1);
p2 = 1/(T-2);
for i=1:I
    for t=2:T
        r(i,t)= (AkcieOdhad(t,i) - AkcieOdhad(t-1,i))/(AkcieOdhad(t-1,i));
    end;
    vyberR(i) = p*sum(r(i,[2:T]));
end;

V = zeros (I);

for i=1:I
    for j=1:I
        V(i,j) = p2 * ( r(i,:) - vyberR(i) )*(( r(j,:) - vyberR(j) )');
    end;
end;

%% 3
% Interpretujte význam omezujících podmínek problému
%
% sum(X(i))=1   ... X(i) ~ vahy akcii i = jaky objem z celkove sumy investovat
%                                         do akcie i
%               ... mame predem danou sumu - tu pak delime podle x, takze
%               cela suma je 1*SUMA = SUMA, rozdeleni na 2 je pak X(1)*SUMA
%               + X(2)*SUMA = SUMA ... analogicky pro N
%
% X(i) >= 0     ... nelze nakoupit zaporny objem akcii...
%
% sum( X(i) * R(i) ) >= Rp ... Rp je pozadovany minimalni vynos, soucet
% vynosu z jednotlivych akcii musi dat minimalne Rp
%
%% 4
% Jak by se zmenila formulace problemu (1), pokud bychom chteli investovat do
% kazde akcie maximálne 30% z celkove hodnoty portfolia ?
%
% pribyla by jeste podminka X(i) <= 0.3
%
%% 5
MinRp = min(vyberR);
MaxRp = max(vyberR);
r_p = [MinRp : (MaxRp-MinRp)/50 : MaxRp ];
cnt = size(r_p);
rp = 0.02;

riziko = zeros(cnt);
vynos = zeros(cnt);

H = V;
f = zeros (I,1);
A = -vyberR';
Aeq = ones(1,I);
beq = 1;
lb = zeros(1,I);
ub = ones(1,I);

for i=1:cnt(2)
    rp = r_p(i);
    b = -rp;
    [x,fval] = quadprog(H,f,A,b,Aeq,beq,lb,ub);
    vynos(i) = x'*vyberR;
    riziko(i) = 2*fval;
end
plot(riziko, vynos);
</code>
Prosim editujte a komentujte...

Tak moje řešení včetně šestky. Ale ta se asi počítá poněkud špatně. Každopádně u nás na cviku to prošlo. Zřejmě je tam naopak to rovnoměrné rozložení a váhové rozložení, které by mělo být úplně nahoře (nejlepší výsledky). --- //[[zacham3@fel.cvut.cz|Martin Zachar]] 2009/11/26 18:38//
<code matlab>
clear vars;
close all;
clear all;
clc;
load akcie.mat;

%% 1 Vypocet parametru Markowitzova modelu

[aoH, aoW] = size(AkcieOdhad);

r = (AkcieOdhad((aoH - 1) : -1: 1, :) - AkcieOdhad(aoH : -1 : 2, :)) ./ AkcieOdhad(aoH : -1 : 2, :);
rprum = 1 / (aoH - 1) * sum (r);

r = r';
rprum = rprum';

V = zeros(aoW);

for i = 1 : aoW
    for j = 1 : aoW        
        V(i,j) = 1 / (aoH - 2) * sum ((r(i, :) - rprum(i)) * (r(j,:) - rprum(j))');        
    end;
end;

%% 5 Vypocet, kolik koupit akcii od ktere spolecnosti

r_p = min(rprum) : (max(rprum) - min(rprum)) / 50 : max(rprum);

Aeq = ones(1, aoW);
beq = 1;

lb = zeros(1, aoW);
A = -rprum';
b = - r_p;
f = zeros(1, aoW);

rx = zeros(size(r_p));
sigma = zeros(size(r_p));

for i = 1 : size(r_p,2)
    [x, fval] = quadprog(V, f, A, b(1,i), Aeq, beq, lb, []);
    rx(i) = x' * rprum;
    sigma(i) = 2*fval;
end;


figure;
hold on;
title('Zavislost r(x) na sigma');
xlabel('risk \sigma^2(x)');
ylabel('r(x)');
plot(sigma, rx, '-');
hold off;


%% 6 a) Od kazde spolecnosti nakoupim stejne akcii
clc; %close all;
x = ones(1,aoW) / aoW;
vklad = AkcieOdhad(1, :) .* x;

vynos = zeros(size(AkcieTest));

for i = size(AkcieTest,1) : -1 : 1
    vynos(i, :) = AkcieTest(i, :).*x - vklad;
end;

vynos = flipud(vynos);

figure;
hold on;
title('Celkovy vynos za jednotlive mesice');
xlabel('Mesice');
ylabel('Vynos');
plot(sum(vynos, 2), '-ro');

%% 6 quadprog - vypoctu vahy jako v uloze 5 a podle toho nakoupim akcie

r_p = 0.037;
b = - r_p;

x = quadprog(V, f, A, b, Aeq, beq, lb, []);
vklad = AkcieOdhad(1, :) .* x';
vynos = zeros(size(AkcieTest));

for i = size(AkcieTest,1) : -1 : 1
    vynos(i, :) = AkcieTest(i, :).*x' - vklad;
end;

vynos = flipud(vynos);

xlabel('Mesice');
ylabel('Vynos');
plot(sum(vynos, 2), '-go');


%% 6 b) Najdu nejsilnejsi spolecnot a do te dam vsechno

maximum = max(x);
for i = 1 : size(x, 1)
    if x(i) == maximum
        x(i) = 1;
    else
        x(i) = 0;
    end;
end;

vklad = AkcieOdhad(1, :) .* x';
vynos = zeros(size(AkcieTest));

for i = size(AkcieTest,1) : -1 : 1
    vynos(i, :) = AkcieTest(i, :).*x' - vklad;
end;

vynos = flipud(vynos);

plot(sum(vynos, 2), '-bo');

h = legend('rovnomerne rozlozeni', 'investice podle vah', 'investovani do nejocekavanejsiho', 'Location', 'Best');
set(h,'Interpreter', 'none');
hold off;
</code>
~~DISCUSSION~~