Take Home Key

Andy Long
Spring, 2020

Let’s start by creating the data:

An experiment on the relation between the velocity of an enzymatic reaction and the substrate concentration was conducted. The concentrations and observed reaction velocities for an experiment in which the enzyme was treated with puromycin are given by

concentration={.02,.02,.06,.06,.11,.11,.22,.22,.56,.56,1.1,1.1};
velocity={76,47,97,107,123,139,159,152,191,201,207,200};
data=Transpose[{concentration,velocity}];

a. Let’s build a model:

I’ll actually try the “linear model”: not a good idea if one has actually seen the data, but it does provide us with a benchmark (MSError: 954.71).

lm=LinearModelFit[data,x,x]
lm["ParameterTable"]
lm["AdjustedRSquared"]
lm["ANOVATable"]
lm["ParameterConfidenceIntervals"]
p1=ListPlot[data,PlotLabel->"velocity of an enzymatic reaction"];
p2=Plot[lm[x],{x,0,1.2}];
Show[p1,p2]

Fit a Michaelis-Menton model to this data,

Now let’s give the “Michaelis-Menton model” a try -- the one I suggested in the exam. We’ll use NonlinearModelFit to find it. I don’t need to give it any guidance, it turns out.

nlm=NonlinearModelFit[data,a c/(b+c),{a,b},c]
nlm["ParameterTable"]
nlm["AdjustedRSquared"]
nlm["ANOVATable"]
nlm["ParameterConfidenceIntervals"]
firstModel=Plot[nlm[x],{x,0.01,1.2},PlotRange->{30,210}];
Show[p1,firstModel,PlotRange->All]
firstRes=ListPlot[Transpose[{concentration,nlm["FitResiduals"]}]]

b. Evaluate the fit. Include a plot of the data with the fitted model, residuals, and
confidence intervals for the parameters.

MSError: 119.545

Both of a and b are significantly different from 0 -- no surprise there! 

The R^2 is essentially useless, as we’ve seen before with non-linear models.  

I’m not really happy with those residuals.  They’re fairly balanced near 0, but high at the far end. 
I’m wondering if an exponent of 1 is the right one?
It’s easy to try a different exponent -- just introduce a new parameter, the exponent on the c terms.

nlm=NonlinearModelFit[data,a c^ex/(b+c^ex),{{a,212},b,{ex,1}},c]
nlm["ParameterTable"]
nlm["AdjustedRSquared"]
nlm["ANOVATable"]
nlm["ParameterConfidenceIntervals"]
p2=Plot[nlm[x],{x,0.01,1.2},PlotRange->{30,210}];
Show[p1,p2,PlotRange->All]
ListPlot[Transpose[{concentration,nlm["FitResiduals"]}]]

MSError: 89.9846

I don’t like that the b term turns up non-significant: if I set b to zero, then I’m left with a constant!:)

Let’s take a slightly different approach: I’m also thinking “maybe this chemical process has a 3/4 power on it”. So I do a shift about .75, and try fitting b to the same power. 

Let’s see what happens....

nlm=NonlinearModelFit[data,a c^(.75-ex)/(b^.75+c^(.75-ex)),{{a,212},b,{ex,.01}},c]
nlm["ParameterTable"]
nlm["AdjustedRSquared"]
nlm["ANOVATable"]
nlm["ParameterConfidenceIntervals"]
p2=Plot[nlm[x],{x,0.01,1.2},PlotRange->{30,210}];
Show[p1,p2,PlotRange->All]
ListPlot[Transpose[{concentration,nlm["FitResiduals"]}]]

MSError: 89.9846

Okay, the p-value of 0.932969 says that ex may as well be set to zero. 

So I’m going to try a model where I assign the exponent to .75, and play with a and b.

nlm=NonlinearModelFit[data,a c^.75/(b^.75+c^.75),{{a,212},{b,.08}},c]
nlm["ParameterTable"]
nlm["AdjustedRSquared"]
nlm["ANOVATable"]
nlm["ParameterConfidenceIntervals"]
finalModel=Plot[nlm[x],{x,0.01,1.2},PlotRange->{30,210},PlotStyle->Orange];
Show[p1,finalModel,PlotRange->All]
finalRes=ListPlot[Transpose[{concentration,nlm["FitResiduals"]}],PlotStyle->Orange]

MSError: 81.0561

I’m happier now -- we eliminated that really ugly first residual, and I think that they look a little more balanced:

Show[firstRes,finalRes,PlotRange->All]
Show[p1,firstModel,finalModel,PlotRange->All]

c. Interpret α.

The parameter α (a, for me) is the asymptotic value of the velocity, as the concentration gets large (goes to infinity, but I don’t think that it can!).  It takes a value of 212 or so (95% Confidence interval {197.205,228.163}).

((212.7 c)/(0.064 +c))