%simulace heteroscedasticita a jeji testovani

clear all
close all
clc
 
N = [50,100,500,1000]; k=3; %pocet pozorovani, pocet vysvetlujicich promennych

for i=1:length(N)
    rand('seed',10); %vyber generatoru
    
    x1 = randn(N(i),1)*5+6;   %1. vysvetlujici promenna ~N(6,100)
    x2 = randn(N(i),1)*10+3;    %2. vysvetlujici promenna ~N(3,16)
    
    %volba typu nahodne slozky
    e = randn(N(i),1).*x1;  % heteroscedastic error term (podle x1)
    
    b = ones(k,1);
    iota = ones(N(i),1);
    x = [iota x1 x2];

    % generate y-data
    y = x*b + e;

    % do ols regression
    reso = ols(y,x);
    figure(1)
    subplot(2,2,i);
    scatter(x1,reso.resid,'x');
    xlabel('x1');
    ylabel('e');
    title(num2str(N(i)));
end

for i=1:length(N)
    rand('seed',10); %vyber generatoru
    
    x1 = randn(N(i),1)*10+6;   %1. vysvetlujici promenna ~N(6,100)
    x2 = randn(N(i),1)*4+3;    %2. vysvetlujici promenna ~N(3,16)
    
    %volba typu nahodne slozky
    e = randn(N(i),1);        % homoscedastic error term
    
    b = ones(k,1);
    iota = ones(N(i),1);
    x = [iota x1 x2];

    % generate y-data
    y = x*b + e;

    % do ols regression
    figure(2)
    reso = ols(y,x);
    subplot(2,2,i);
    scatter(x1,reso.resid,'x');
    xlabel('x1');
    ylabel('e');
    title(num2str(N(i)));
end
