%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
%% FILE: uu.m
%%
%% Stefan Stuidun puoliso Urho Uraohjus miettii marraskuussa 2020 mita tehda 
%% syyskuu 2021-elokuu 2022 
%%
%% Han on paattanyt analysoida Stefan Stuidun vaihtoehdot
%%
%% a1 = Opiskelijavaihto Kiinaan
%% a2 = Keskittyminen opintoihin
%% a3 = Opintoja ja ruokalahetin toita iltaisin
%% a4 = Vahan opintoja ja paatyo isossa kv. yrityksessa
%%
%% Seuraavat skenaariot olivat Stefanin mielesta keskeisia
%%
%% s1 = Korona riehuu, yliopisto etaopetuksessa ja maailma lamassa
%% s2 = Korona riehuu, yliopisto etaopetuksessa ja maailman talous elpyy
%% s3 = Korona riehuu, yliopisto normiopetuksessa ja maailma lamassa
%% s4 = Korona riehuu, yliopisto normiopetuksessa ja maailman talous elpyy
%% s5 = Korona kukistettu ja maailma lamassa
%% s6 = Korona kukistettu ja maailman talous elpyy
%%
%% Skenaarioiden todennakoisyydet ovat tuntemattoman koronatodennakoisyyden p0
%% funktioina
%% 
%% p1 = p0/4
%% p2 = p0/4
%% p3 = p0/4
%% p4 = p0/4
%% p5 = (1-p0)/2
%% p6 = (1-p0)/2
%%
%% Koronatodennakoisyys p0 on herkkyysparametri.
%%
%% Urho perustaa paatoksensa odotusarvosaannolle ja katumuksen kaihtajan 
%% saannolle seka palkkiolle kuukausipalkka kolmen vuoden paasta (syksy 2023)
%%
%% Palkkiot oletetaan (varsin eparealistisesti) tunnetuiksi.
%%
%% Odotusarvosaannon ja katumuksen kaihtajan saannon yhdistamispaino w on
%% herkkyysparametri.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Paatosmatriisi 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

R = [
3000  3500  2800  3500  2800  3500;  %% Kiinaan
2500  2700  3300  3800  3500  4000;  %% opintoja
3000  3200  2000  2500  2000  2500;  %% ruokalahetti
3000  3200  2500  4000  3000  3500;  %% paatyot
];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Pienehko taulukko p0 herkkyyksille.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% p0 ja w otteistukset
p0arvot = 0:0.25:1;
warvot = 0:0.25:1;

%% Paatosongelmien vakiokentat
po1.pmat = R;
po1.psaanto = "oa";
po2.pmat = R;
po2.psaanto = "kk";

%% Skaalattu paatossaanto V2
pr2 = pmatriisi(po2);
V2 = -pr2.parvot;      %% Maksimointimuoto
V2 = (V2-min(V2))/(max(V2)-min(V2));

%% Skaalattu paatokset V1 
for i = 1:length(p0arvot)
    p = [
         p0arvot(i)/4;
         p0arvot(i)/4;
         p0arvot(i)/4; 
         p0arvot(i)/4;
         (1-p0arvot(i))/2;
         (1-p0arvot(i))/2
        ];
    po1.tn = p;
    pr1 = pmatriisi(po1);
    V1(:,i) = pr1.parvot;
    V1(:,i) = (V1(:,i) - min(V1(:,i)))/(max(V1(:,i))-min(V1(:,i)));
end
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Visualisointi
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Otteistustaajuus.  Mita isompi sen parempi, mutta laskeminen kestaa.
%% Jos N=500 laskeminen kestaa noin 20 sekuntia.
N=500;      

%% Laskut
p0arvot = linspace(0,1,N);
warvot = linspace(0,1,N);
z = zeros(N);

for i=1:length(p0arvot)
    for j=1:length(warvot)
        p0 = p0arvot(i);
        w = warvot(j);
        q0 = 1-p0;
        po1.tn = [p0/4 p0/4 p0/4 p0/4 q0/2 q0/2]';
        pr1 = pmatriisi(po1);
        v1 = pr1.parvot;
        v1 = (v1-min(v1))/(max(v1)-min(v1));
        v = w*v1 + (1-w)*V2;
        [roskaa, z(j,i)] = max(v);
    end    
end
    
%% Itse kuva    
contourf(p0arvot, warvot, z, 1:4)
set(gca, "fontsize", 24);
xlabel("p_0")
ylabel("w")
puna  = [1 0 0];
viher = [0 1 0];
colormap([puna; 0.6*puna; viher; 0.6*puna])  %% Likainen jippo!!!!
%% a1 = puna
%% a2 = tumma puna
%% a3 = viher
%% a4 = tumma viher
text(0.80,0.20,"a_1", "fontsize",24)
%text(0.20,0.10,"a_2", "fontsize",24)
%text(0.85,0.60,"a_3", "fontsize",24)
text(0.20,0.60,"a_2", "fontsize",24)
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Herkkyyskuvan piirto PDF-tiedostoon.  Poista kommentit, jos haluat PDF:n
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%orient("landscape");
%print("uu_herkkyyskuva.pdf");