Page 20
ht=[75 92 108 121 130 142 155];
age=[1 3 5 7 9 11 13];
sumy=sum(ht); sumx=sum(age);
age2=age.*age;
sumx2=sum(age2);
ageht=age.*ht;
sumxy=sum(ageht);
m=(7*sumxy-sumx*sumy)/... %3 dots mean continue to next line
(7*sumx2-sumx^2)
b=(sumx2*sumy-sumx*sumxy)/...
(7*sumx2-sumx^2)
Page 21
% Now the matrix solution
% matrix of independent variable
% experimental values as columns
MT=[1 3 5 7 9 11 13; 1 1 1 1 1 1 1]; % two rows
M=MT'; % transpose to columns
% M = transpose of MT
% dependent variable data next, as col. vec.
Y=[75; 92; 108; 121; 130; 142; 155];
M\Y % Matlab syntax for leastsquare
% plot data and fit for comparison
plot(age,ht,'o') % point plot ht vs age with circles
hold on
fit=m*age+b;
plot(age,fit); xlabel('age'); ylabel('Height (cm)');
Page 22/23
% year by year cases; note ellipses to continue the line
AIDS=[97 206 406 700 1289 1654 2576 3392 4922 ...
6343 8359 9968 12990 14397 16604 17124 19585 ...
19707 21392 20846 23690 24610 26228 22768];
CAC=cumsum(AIDS)/1000; %cumumlative (scaled down 1000)
% housekeeping to get the sequence 0,0.5,1,1.5, ...
s=size(AIDS); % number of half years
count=[0:s(2)-1];
time =1981+count/2;
% shifted and scaled time
scaledTime=(time-1980)/10
% log the data to do a log-log plot
lnCAC=log(CAC)
lnTime=log(scaledTime)
% form the coefficient matrix
% for lnCAC = k*lnTime + b fit
MT=[lnTime; ones(1,24)] %second row is ones
M=MT';
params=M\(lnCAC') % do the leastsquares
k=params(1)
A=exp(params(2))
% now compare the fit to the data
% in log-log space
plot(lnTime,lnCAC,'o')
lnFit= params(1).*lnTime+params(2)
plot(lnTime,lnFit)
Page 24
% and compare in regular space
hold off
plot(time,CAC)
CACFit=exp(params(2)).*scaledTime.^params(1)
plot(time,CACFit)
Exercises/Experiments
1.
ht=[62,63,64,65,66,67,68,69,70,71,72,73,74];
wt=[128,131,135,139,142,146,150,154,158,162,167,172,177];
MT=[ht.^3; ht.^2; ht; ones(1,13)];
params=MT'\wt';
plot(ht,wt,'x') hold on
fit=params(1)*ht.^3+params(2)*ht.^2+params(3)*ht+params(4);
plot(ht,fit)
errorLinear=sum((4.04*ht-124.14-wt).^2)
errorcubic=sum((fit-wt).^2)
2.
age60=[0,10,20,30,40,50,60,80,100];
percent60=[100,98.5,98,96.5,95,92.5,79,34,4];
MT=[age60.^4; age60.^3; age60.^2; age60; ones(size(age60))];
parms=MT'\percent60'
fit=parms(1)*age60.^4+parms(2)*age60.^3+...
parms(3)*age60.^2+parms(4)*age60+parms(5);
plot(age60,percent60,age60,fit)
3.
AIDS=[97, 206, 406, 700, 1289, 1654,2576, 3392, 4922, 6343,...
8359, 9968,12990,14397, 16604, 17124, 19585, 19707, 21392,...
20846, 23690, 24610, 26228, 22768];
CAC=cumsum(AIDS)/1000;
s=size(AIDS); % number of half years
count=[0:s(2)-1];
Time =1981+count/2;
pts=[Time' CAC'];
plot(pts(:,1),pts(:,2)); hold on
Times=(Time-1980)/10; LnCAC=log(CAC);
MT=[Times; ones(1,s(2))];
params=MT'\LnCAC'
k=params(1); A=params(2);
y=exp(A)*exp(k.*Times);
plot(10*Times+1980,y)
4.
CW=[24.2,22.9,27.,21.5,23.5,22.4,...
23.8,25.5,24.5,25.5,22.,24.5];
GSW=[38.5,26.,34.,25.5,37.,30.,...
34.,43.5,30.5,36.,29.,32];
MT=[CW; ones(size(CW))];
parmsW=MT'\GSW';
plot(CW,GSW,'x'); hold on
x=21:28; plot(x,parmsW(1)*x+parmsW(2))
%%%
CM=[28.5,24.5,26.5,28.25,28.2,...
29.5,24.5,26.9,28.2,25.6,28.1,...
27.8,29.5,29.5,29];
GSM=[45.8,47.5,50.8,51.5,55.0,...
51.,47.5,45.,56.0,49.5,57.5,...
51.,59.5,58.,68.25];
MT=[CM; ones(size(CM))];
parmsM=MT'\GSM'
plot(CM,GSM,'o')
x=24:30;
plot(x,parmsM(1)*x+parmsM(2))}