clear; echo off;

% PARAMETERS:
% ==========

N = 100;                    % Number of data.
d = 2;                      % Input dimension.
%kernelFunction = 'expK';   % Kernel function.
%kernelFunction = 'sigmoidK';   % Kernel function.
kernelFunction = 'polyK';   % Kernel function.
bias = 0;                   
order = 2;
kernelPar =[order bias];    % Polynomial kernel parameters.



% GENERATE THE DATA:
% =================

theta = [2 3]';
x = rand(N,d);
y = (x.^2)*theta;  % quadratic regression.

% CREATE THE TRAIN AND TEST SETS:
% ==============================
trainSize = round(N/2);        % Number of training examples.
xTrain = x(1:trainSize,:);     % Training input data.
yTrain = y(1:trainSize,:);     % Training ouput data.
xTest = x(trainSize+1:N,:);    % Test input data.
yTest = y(trainSize+1:N,:);    % Test output data.


% TRAIN KERNEL MACHINE:
% =====================

KTrain = feval(kernelFunction,xTrain,xTrain,kernelPar);  % Kernel matrix.
deltasq = 0.1;  % Ridge parameter.
alpha = pinv(KTrain+deltasq*eye(trainSize,trainSize))*yTrain;  


% TEST THE LINEAR MODEL:
% =====================

KTest = feval(kernelFunction,xTest,xTrain,kernelPar);  % Kernel matrix.
yPredTrain = KTrain * alpha;  % Generate prediction.
yPredTest  = KTest * alpha;   % Generate prediction.
trainError = mean( (yTrain-yPredTrain).^2 ); % RMS error.
testError = mean( (yTest-yPredTest).^2 );    % RMS error.
disp(' ');
disp('Errors');
disp('------');
disp(' ');
disp(['Train error  = ' num2str(trainError)]);
disp(['Test error   = ' num2str(testError)]);
disp(' ');

% PLOT THE TRAIN ERROR AND THE TEST ERROR:
% =======================================
figure(2)
clf;
subplot(211)
plot(1:trainSize,yTrain,'ro',1:trainSize,yPredTrain,'b')
title('Training set');
zoom on;
subplot(212)
plot(trainSize+1:N,yTest,'ro',trainSize+1:N,yPredTest,'b')
title('Test set');
legend('True value','Prediction');
