Problem 7

In this problem we solve a combined linear and
non-linear least squares regression problem.
We are trying to look at the growth of bacterial
cells over time.  There are two species present
in our nutrient broth, each with different growth
rates and initial concentrations, however we can't
easily distinguish between them.  Instead we measure
the total concentration of cells as a function of
time, and try to infer the growth rate of each of
the two species.  Thus:

   Total = A * exp (ka*t) + B * exp (kb*t)

We have the observed data on total cell concentration
given below.  We wish to use least-squares regression
to calculate the two growth rate constants.  Note that
you cannot linearize the above equation!  Solve this
problem by formulating it in the least squares sense
(e.g., the sum of the square of the deviation between
the model prediction for A, B, ka, and kb and the
observed data points).  Since the problem is linear in
the parameters A and B, write a function which returns
the optimum values of these parameters and the sum of
squares for fixed ka and kb.  This will be a linear least
squares problem.  Then write an outer program which
does the non-linear optimization for ka and kb.  This
is much faster than solving the problem as a four parameter
non-linear optimization problem.  You may use the matlab
routine fminsearch for the two-parameter non-linear optimization
problem.  Type help fminsearch to learn more about its use.

In addition to getting the best fit values of the parameters,
determine the matrix of covariance.  This is most easily
done by 1) using the deviation between the data points and
the fitted model to calculate the variance in the data
measurements themselves, 2) determining the change in the
best fit parameter values with changes in each of the data
points (e.g., take the derivative grad f where f is the
array of fitted parameter values and the gradient is with
respect to each of the n data points.  Grad f is thus a
4 x n matrix, which is computed numerically), and 3) use
the dependence matrix and the data variance to calculate the
matrix of covariance of the parameters.  This will only work
if the variance and / or the dependence is small (the penalty
of ignoring higher order terms when calculating the variance).

Using all this, I want you to do the following:

a.  Determine the best fit model parameters for the data: the
two growth rates and the two initial concentrations.

b.  Determine the matrix of covariance of the fitting parameters
from the scatter in the data.

c.  Calculate the 95% confidence interval for the two growth rates.

d.  Plot up the data, the model, and the 95% confidence interval of
the model prediction.  This last bit is easy to do once you have the matrix
of covariance of the fitting parameters - since the cell concentrations
predicted by the model are just a function of the fitting parameters, 
the standard deviation of the model prediction at each time can be 
calculated via the standard error propagation formula.  Although
the model prediction is non-linear in the modelling parameters, it
is simple to take the necessary gradient analytically.  The gradient will
thus be an nx4 array, where each row corresponds to evaluating the
gradient at a different modelling time.


A couple of useful tips:

1. Start close to the correct value for the growth rates. One of the 
growth constants will be negative (some of the cells are dying off).

2. When computing the dependence of the fitting parameters
on the data, don't vary the data by too small amount, since the
routine fminsearch has a preset target precision.  If you make too small
a change, then you can't compute the derivative you need!  I found
that a change in data scaled with the standard deviation of the
data was effective in getting grad f.



The data is:

       t         c
         0    4.1077
    0.0500    3.8447
    0.1000    3.5597
    0.1500    3.2212
    0.2000    3.0120
    0.2500    2.8910
    0.3000    2.8429
    0.3500    2.3134
    0.4000    2.3515
    0.4500    2.5100
    0.5000    2.4260
    0.5500    2.3457
    0.6000    2.4561
    0.6500    2.3981
    0.7000    2.5175
    0.7500    2.8218
    0.8000    3.1638
    0.8500    3.3671
    0.9000    3.5864
    0.9500    4.1444
    1.0000    4.5336