Explore BrainMass
Share

Matlab normal equations without polyfit

This content was STOLEN from BrainMass.com - View the original, and get the already-completed solution here!

See attached file for full problem description.

Does anyone know how to do this without using Matlab functions such as polyfit and etc? The regression and error analysis needs to be performed by solving the normal equations. I also attached an example of how it is supposed to be modeled.

---
%Tossing a lead weight: example of linear regression error from BME221 lecture

%matrix of modeling functions: first column of ones, second column = t,
%third column = -0.5*t^2 (data is modeled by h = h0 +V0*t +0.5g*t^2).
%We shall denote our model by the vector equation A*x=b
A=[ones(1,5);0:4;-0.5*(0:4).^2]'
%height measurements in meters
h=[0.3 14.8 20.7 15.9 2]'
K=inv(A'*A)*A'
% our modeling parameters, (h0,V0,g)
x=K*h
%variance in data
Shsq=sum((h-A*x).^2)/(5-3)
%variance in first modeling parameter, for demonstration only.
%other variances and covariances can be extracted from the symmetric SIGMAxsq.
Sx1sq=K(1,:)*K(1,:)'*Shsq
%SIGMAbsq: matrix of covariance of the observations b.
%If the observations are independent then SIGMAbsq will be diagonal only.
%We further assume that all observations have the same variance.
%(a common assumption)
SIGMAbsq=eye(5)*Shsq
%matrix of covariance of the regression parameters
SIGMAxsq=K*SIGMAbsq*K'
%Suppose we wish to determine the model value at 100 points ranging from 0 to 4.
%Thus, we now have the vector equation: bm = Am*x
Am=[ones(1,100);linspace(0,4,100);-0.5*linspace(0,4,100).^2]';
%Finally, we have the matrix of covariance of bm.
%In general, only the diagonal elements of SIGMAbmsq are of interest.
sigmabmsq=diag(Am*SIGMAxsq*Am');
plot(A(:,2),h,'o',Am(:,2),Am*x)
hold on
plot(Am(:,2),Am*x+sqrt(sigmabmsq),'r--',Am(:,2),Am*x-sqrt(sigmabmsq),'r--')
hold off
legend('data','best-fit model','69% confidence (1sigma) interval')
xlabel('t (sec)')
ylabel('h (m)')
---

https://brainmass.com/computer-science/matlab/matlab-normal-equations-without-polyfit-51678

Solution Preview

The example you gave did indeed help.

The attached file follows the example you gave.

function Michaelis3()

% Data for Rxn rate with progesterone - case 1: with progesterone
dog1 = [ 5 10 20 25 30 40 50 60 ] ; % dog is C in the formula
rxn1 = [ 3.25 2.4 5.1 6 7.2 7.7 9.1 9.05 ] ; % Rxn is R in the formula

% The variables in which we do linear regression:
% x = 1/C and y = 1/R
x1 = 1./dog1' ;
y1 = 1./rxn1' ;

% model: y1 = (Km/Rmax)*x1 + (1/Rmax)
% matrix of modeling functions:
A1 = [ ones(1,length(x1)) ; x1' ]' ;
K1 = inv(A1'*A1)*A1'
% Modeling parameters: ModPars = [ (1/Rmax), (Km/Rmax) ]
ModPars1 = K1 * y1

% Now we know what Km and Rmax are in the fit:
Rmax1 = 1/ModPars1(1)
Km1 = ModPars1(2)/ModPars1(1)
% This answers the 1st question
% but we also have to look for confidence intervals

% Data Variance: subtract 2 as here we have 2 modeling parameters
% and so used 2 degrees of freedom in the ...

Solution Summary

The expert uses Matlab normal equations without using the polyfit function.

\$2.19

Matlab model using normal equations model to analyze the body's response to a bolus injection of antibiotic

(See attached file for full problem description)

---
This model is used to analyze the body's response to a bolus injection of this antibiotic:
B(t) = B0e(-k*t)
Using this model, and data provided for the average plasma concentration data in ug/ml with time in hours :
a) perform linear regression on this penicillin clearance using the normal equations
b) plot the concentration data and re-plot the log of the concentration data to show that it appears linear.
c) Determine initial concentration, B0, and elimination rate constant, k

Concentration (ug/ml) Time (min)
89 22
60 44
30 88
15 132
7 176

a) for linear regression
ln( B ) = ln( B0 )- Kt

I'm not sure how to use the normal equations in MATLAB to plot this? We're not supposed to use polyfit function... I think you have to make a 2 column matrix with a column of 1's and a column of the time?
---

(See attached file for full problem description)

View Full Posting Details