load('../Data/StudData2012.txt')

kjonn=StudData2012(:,1);
hogde=StudData2012(:,4);
mor=StudData2012(:,5);
far=StudData2012(:,6);

f=find(1-isnan(mor) & (1-isnan(far)));

%Hogde barn, jenter justert til gutar
hogdeBarn=hogde(f).*(1+0.08*(kjonn(f)-1));
%Gj.snitt hogda til foreldra. Mor*1.08
hogdeForeldre=0.5*far(f)+0.5*1.08*mor(f);

plot(hogdeForeldre, hogdeBarn,'.')
xlabel('Høgde foreldre','FontSize',16)
ylabel('Høgde barn','FontSize',16)
%estimerer regresjonsparametra
[b,bint,r,rint,stats]=regress(hogdeBarn,[ones(length(f),1),hogdeForeldre]);
hold
plot([160,200],[b(1)+b(2)*160, b(1)+b(2)*200],'r','LineWidth',2)

%plotter hogdebarn=hogdeforeldre linje
%plot([160,200],[160,200],'k','LineWidth',2)

a=b(1)
bb=b(2)
sigmae=stats(4)
sqrt(sigmae)



%Predikerer voksenhogde for eit barn med gj.snitt foreldrehogde 182.5
x0=182.5;
ypred=a+bb*x0


%Finn Pred.varians
Sxx=sum((hogdeForeldre-hogdeBarn).^2)
predvar=sigmae*(1+1/length(hogdeForeldre)+(x0-mean(hogdeBarn)).^2/Sxx)
sqrt(predvar)

%Finn 95% pred.intervall
ypred-1.96*sqrt(predvar)
ypred+1.96*sqrt(predvar)


%Finn KI for barn med gj.snitt foreldrehogde 182.5
muvar=sigmae*(1/length(hogdeForeldre)+(x0-mean(hogdeBarn)).^2/Sxx)
sqrt(muvar)
%Finn 95% KI
ypred-1.96*sqrt(muvar)
ypred+1.96*sqrt(muvar)


%Finn KI for mangen x0
x0=160:1:200;
ypred=a+bb*x0;
muvar=sigmae*(1/length(hogdeForeldre)+(x0-mean(hogdeBarn)).^2/Sxx);
KInedre=ypred-1.96*sqrt(muvar);
KIovre=ypred+1.96*sqrt(muvar);
plot(x0,KIovre,'r--', 'LineWidth',2)

plot(x0,KInedre,'r--','LineWidth',2)


%Finn PI for mangen x0
x0=160:1:200;
ypred=a+bb*x0;
predvar=sigmae*(1+1/length(hogdeForeldre)+(x0-mean(hogdeBarn)).^2/Sxx);
PInedre=ypred-1.96*sqrt(predvar);
PIovre=ypred+1.96*sqrt(predvar);
plot(x0,PIovre,':g','LineWidth',2)
plot(x0,PInedre,':g','LineWidth',2)


%Hogde far og hogde mor

mean(mor(f))
mean(far(f))
std(mor(f))
std(far(f))

figure
plot(far(f),mor(f),'.')

xlabel('Høgde far','FontSize',16)
ylabel('Høgde mor','FontSize',16)

[b,bint,r,rint,stats]=regress(mor(f),[ones(length(f),1),far(f)]);
hold
plot([150,200],[b(1)+b(2)*150, b(1)+b(2)*200],'r','LineWidth',2)

[b,bint,r,rint,stats]=regress(far(f),[ones(length(f),1),mor(f)]);