% Janne Koljonen
% Vaasa yliopisto
% 15.5.2008
%
% Soft Computing
%
% Harjoitus 9, tehtävä 2
% Ennustetaan Vaconin päätöskurssia 10 edellisen päivän päätöskurssilla.
% Yksikerroksinen Perceptron, opetus pseudoinversillä

% Päätöskursseja
% Opetusdata
prices_cal=[28.6;28.5;27.81;28.69;29.76;28.74;29.33;29.11;29.10;28.29;28.25;27.60;27.80;27.52;27.30;27.10;28.00;28.00;28.26;28.6;28.5;27.81;28.69;29.76;28.74;29.33;29.11;29.10;28.29;28.25;27.60;27.80;27.52;27.30;27.10;28.00;28.00;28.26];
% Testidata
prices_test=[28.29;28.25;27.60;27.80;27.52;27.30;27.10;28.00;28.00;28.26;28.90;29.26;29.50;29.99;29.35];

% Muodostetaan piirrevektoreista matriisit X=[x1'; x2'; ...], jossa vaakavektori xi' sisältää kymmenen 
% peräkkäisen päivän päätöskurssit alkaen indeksistä i (siis indeksit i:i+9). Piirrevektorit ovat siis riveillä.
% Lisäksi muodostetaan halutuista vasteista vektori D=[d1; d2; ...], jossa
% di on xi:n haluttu vaste eli päätöskurssi indeksistä i+10.
X_cal=[];
D_cal=[];
for i=1:size(prices_cal,1)-10
    X_cal(i,:)=prices_cal(i:i+9,1);
    D_cal(i,1)=prices_cal(i+10,1);
end
X_test=[];
D_test=[];
for i=1:size(prices_test,1)-10
    X_test(i,:)=prices_test(i:i+9,1);
    D_test(i,1)=prices_test(i+10,1);
end

% Opetetaan yksikerroksinen Perceptron, jolla lineaarinen neuroni, pseudoinverssimenetelmällä
W=inv(X_cal'*X_cal)*X_cal'*D_cal;

% Lasketaan hermoverkon ulostulot opetus- ja testijoukoille
Y_cal=X_cal*W;
Y_test=X_test*W;

% Lasketaan virheet kullekin näytteelle
E_cal=D_cal-Y_cal;
E_test=D_test-Y_test;

% Lasketaan neliövirheiden summat
CRESS=E_cal'*E_cal; % Calibration Residual Error Sum of Squares
PRESS=E_test'*E_test; % Prediction Residual Error Sum of Squares

% Lasketaan RMS-virhe näytettä kohden
RMSEC=sqrt(CRESS/size(E_cal,1)) % Root-Mean-Square Error of Calibration
RMSEP=sqrt(PRESS/size(E_test,1)) % Root-Mean-Square Error of Prediction


% Piirretään kurssikäyrät ja ennustetut kurssit
figure
plot([1:10],prices_cal(1:10,1),'b')
hold on
plot([10+1:10+size(D_cal,1)],D_cal,'k')
plot([10+1:10+size(D_cal,1)],Y_cal,'k--')
plot([10+size(D_cal,1)+1:10+size(D_cal,1)+size(D_test,1)],D_test,'r')
plot([10+size(D_cal,1)+1:10+size(D_cal,1)+size(D_test,1)],Y_test,'r--')
XLabel('Päivä')
YLabel('Päätöskurssi [€]')

% Piirretään virhekäyrät
figure
hold on
plot([10+1:10+size(D_cal,1)],E_cal,'k')
plot([10+size(D_cal,1)+1:10+size(D_cal,1)+size(D_test,1)],E_test,'r')
XLabel('Päivä')
YLabel('Virhe [€]')

% Piirretään %-virhekäyrät
figure
hold on
plot([10+1:10+size(D_cal,1)], 100*(E_cal./D_cal),'k')
plot([10+size(D_cal,1)+1:10+size(D_cal,1)+size(D_test,1)], 100*(E_test./D_test),'r')
XLabel('Päivä')
YLabel('Virhe [%]')

% Päivävolatiliteetti eli päivätuoton keskihajonta
% Päivätuotot
EPD=prices_cal(2:size(prices_cal,1))-prices_cal(1:size(prices_cal,1)-1);
Vola=std(EPD); % €

% Ennustusvirhe suhteessa volatiliteettiin
RMSEC/Vola
RMSEP/Vola

figure
hold on
plot([10+1:10+size(D_cal,1)], E_cal./Vola,'k')
plot([10+size(D_cal,1)+1:10+size(D_cal,1)+size(D_test,1)], E_test./Vola,'r')
XLabel('Päivä')
YLabel('Virhe [volatiliteettiyksikköä]')
