%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% MATLAB script to predict population in the US in 2017, based on census
% data from 1900 to 2010.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Load census data
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

load census.dat;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Create new dates, up to 2017
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

dates = [census(:,1); 2017];
n2017 = max(size(dates));

predict = zeros(5,2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Start with polynomial of degree 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p2 = polyfit(census(:,1),census(:,2),2);
y2 = polyval(p2,dates);

predict(1,1) = 2;
predict(1,2) = y2(n2017);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Repeat with polynomial of degree 4
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p4 = polyfit(census(:,1),census(:,2),4);
y4 = polyval(p4,dates);

predict(2,1) = 4;
predict(2,2) = y4(n2017);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Repeat with polynomial of degree 6
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p6 = polyfit(census(:,1),census(:,2),6);
y6 = polyval(p6,dates);

predict(3,1) = 6;
predict(3,2) = y6(n2017);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Repeat with polynomial of degree 8
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p8 = polyfit(census(:,1),census(:,2),8);
y8 = polyval(p8,dates);

predict(4,1) = 8;
predict(4,2) = y8(n2017);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Repeat with polynomial of degree 10
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p10 = polyfit(census(:,1),census(:,2),10);
y10 = polyval(p10,dates);

predict(5,1) = 10;
predict(5,2) = y10(n2017);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Plot results
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

figure
plot(census(:,1),census(:,2),'ok');
hold on
plot(dates,y2,'-r');
plot(dates,y4,'-b');
plot(dates,y6,'-g');
plot(dates,y8,'-c');
plot(dates,y10,'-m');
h2=plot(dates(n2017),y2(n2017),'o');
set(h2,'MarkerEdgeColor','none','MarkerFaceColor','r')
h4=plot(dates(n2017),y4(n2017),'o');
set(h4,'MarkerEdgeColor','none','MarkerFaceColor','b')
h6=plot(dates(n2017),y6(n2017),'o');
set(h6,'MarkerEdgeColor','none','MarkerFaceColor','g')
h8=plot(dates(n2017),y8(n2017),'o');
set(h8,'MarkerEdgeColor','none','MarkerFaceColor','c')
h10=plot(dates(n2017),y10(n2017),'o');
set(h10,'MarkerEdgeColor','none','MarkerFaceColor','m')
legend('census','d=2','d=4','d=6','d=8','d=10','Location','NorthWest');
