function EnzymeKinetEst
% Enzyme Kinetics Parameter Estimation
%
%   Author: HUANG Huajiang
%   Copyright 2003 UNILAB Research Center, 
%   East China University of Science and Technology, Shanghai, PRC
%   $Revision: 1.0 $  $Date: 2003/05/16 $
%
%  [Ref]: Englezos & Kalogerakis. Applied Parameter Estimation for 
%   Chemical Engineers. Marcel Dekker, 2000（p.323）

clear all
clc

S = [20   15   10    5    2.5]';
rate = [330  300  260   220  110]';

% 線性回歸
y = 1./rate;
X = [ones(size(y))  1./S];
[b,bint,r,rint,stats] = regress(y,X);
rmax_LinReg = 1/b(1);
Km_LinReg = rmax_LinReg * b(2);
R2_LinReg = stats(1);           % R2_LinReg：相關(guān)系數(shù)的平方值

% 非線性回歸
x = S;
y = rate;
b0 = [rmax_LinReg  Km_LinReg];  % 以線性回歸的結(jié)果作為非線性回歸的初值
[beta,resid,J] = nlinfit(x,y,@myfun, b0);
ci = nlparci(beta,resid,J);
ratec = myfun(beta, S);         % 速率擬合值(計算值)

% 非線性擬合效果圖(實驗與擬合的比較)
S4plot = linspace(S(1),S(end),300);
rate4plot = myfun(beta, S4plot);
plot(x,y,'o',S4plot,rate4plot,'-')
legend('Exp','Model')
xlabel('S（μM）')
ylabel('r（μM/min）')

% 殘差關(guān)于擬合值的殘差圖
figure
plot(ratec,resid,'*')
xlabel('擬合（μM/min）')
ylabel('殘差R (μM/min)')
refline(0,0)

% 參數(shù)辨識結(jié)果
fprintf('Estimated Parameters by Linear Regression:\n')
fprintf('\trmax = %.1f\n',rmax_LinReg)
fprintf('\tKm = %.3f\n',Km_LinReg)
fprintf('\tR2 = %.3f\n',R2_LinReg)

fprintf('Estimated Parameters by Nonlinear Regression:\n')
fprintf('\trmax = %.1f\n',beta(1))
fprintf('\tKm = %.3f\n',beta(2))

% ------------------------------------------------------------------
function yhat = myfun(beta, x)
yhat = beta(1)*x./(beta(2)+x);      % beta(1) = rmax, beta(2) = Km