% ***** GRADIENT DESCENT QUADRATIC REGRESSION ***** %

% m:  Number of data samples
% n:  Number of factors or dependent variables
% p:  Index value of a single sample within the set of m samples
% q:  Index value of an individual feature within the 
%     ith sample's feature (or variable) vector
% x:  Feature vector with n+1 features, i.e. x=[x0 x1 x2 x3] & x0=1
% theta:  Parameter vector with at least n+1 parameters, 
%         i.e. theta=[theta0 theta1 theta2 theta3]
% h:  "Hypothesis" or prediction model, i.e. h=theta*transpose(x)
% J:  Cost function dependent on the parameter values in vector theta 

clear all

% Data Points
x=[0.1 1.2 3.0 4.6 4.9 5.7 6.2 7.0 8.4 10.1 10.3 10.9 12];
y=[1.3 4.9 3.5 3.9 5.6 4.3 6.7 4.4 7.4 7.8 9.2 11.4 10.4];

% Initialize parameters
theta1(1)=0; % This value changes each iteration
theta2(1)=0.02; % This value changes each iteration
theta3(1)=0.9; % This value changes each iteration

% Initialize "momentum", alpha
alpha=[0.01 0.001 0.0001];  % These values do not change during interation

m=length(y);  % Number of data points

h=theta1(1)+theta2(1).*x+theta3(1).*(x.^2);  % Hypothesized model

J(1)=(1/m)*sum((h-y).^2);  % Cost function

e=1;

count=2;  % count is the number of iterations that have occured plus one

while (e>0.001)

    theta1(count)=theta1(count-1)-alpha(1).*(2/m).*sum(h-y);
    theta2(count)=theta2(count-1)-alpha(2).*(2/m).*sum((h-y).*x);
    theta3(count)=theta3(count-1)-alpha(3).*(2/m).*sum((h-y).*(x.^2));
    h=theta1(count)+theta2(count).*x+theta3(count).*(x.^2);
    J(count)=(1/m)*sum((h-y).^2);
        
    e=abs(J(count)-J(count-1));
    
    count=count+1;
    
    plot(x,y,'bo')
    axis([0 12.5 -1 12.5])
    grid
    hold on
    plot(x,h,'r')
    pause(0.1)
    hold off
        
end

figure
loglog([1:count-1],J)
grid
hold off