% example 4 in Fall 1997 : Student Version 5 in Mac
% curve fitting by the least squares method
% read Example 8.2 in the textbook by Prof Nakamura
% ______________________________________________________________________
%
% data for curve fitting
x=[0.15,0.4,0.6,1.01,1.5,2.2,2.4,2.7,2.9,3.5,3.8,4.4,4.6,5.1,6.6,7.6]';
y=[4.4964,5.1284,5.6931,6.2884,7.0989,7.5507,7.5106,8.0756,7.8708,...
8.2403,8.5303,8.7394,8.9981,9.1450,9.5070,9.9115]';
%
% curve fitting by using the Log transformation
%
coef=polyfit(log(x),log(y),1)
[coef(1),exp(coef(2))]
g=exp(coef(2))*x.^coef(1);
plot(x,y,'+',x,g)
xlabel('x')
ylabel('y')
title('Curve Fitting by polyfit and Log transformation')
pause
%
% least squares & Newton method
%
% symbolics for the Newton method
imax=size(x)
f=0;
syms a b
for i=1:imax
	f=f+(y(i)-b*x(i)^a)^2;
end
f
f1=diff(f,a);
f2=diff(f,b);
j11=diff(f1,a);
j21=diff(f2,a);
j22=diff(f2,b);
% numerics for the least squares
na=0.2;
nb=5;
tol=1;
tolallow=10^(-5);
iteration=0;
maxiteration=5;
while tol>tolallow
	iteration=iteration+1;
	nf1=double(subs(subs(f1,a,na),b,nb));
	nf2=double(subs(subs(f2,a,na),b,nb));
	nf=[nf1;nf2];
	nj11=double(subs(subs(j11,a,na),b,nb));
	nj21=double(subs(subs(j21,a,na),b,nb));
	nj22=double(subs(j22,a,na));
	nj=[nj11,nj21;nj21,nj22];
	d=-nj\nf;
	tol=norm(d);
	na=na+d(1);
	nb=nb+d(2);
	[iteration,tol,na,nb]
end
iteration
[na,nb]
g=nb*x.^na;
plot(x,y,'*',x,g)
xlabel('x')
ylabel('y')
title('LeastSquares Fitting by the Newton Method')

