To Start R: Click on the R icon
If no R icon, go to Start and then All Programs to find R
R function: nls
> help(nls)
What is the default algorithm?
Here is an example:
# fitting Michaelis and Menten's original data
conc <- c(0.3330, 0.1670, 0.0833, 0.0416, 0.0208, 0.0104, 0.0052)
vel <- c(3.636, 3.636, 3.236, 2.666, 2.114, 1.466, 0.866)
Micmen <- data.frame(conc=conc, vel=vel)
Let's plot the data.
> plot(conc,vel, pch="+")
You must give starting values for the optimization.
fit <- nls(vel~Vm*conc/(K+conc),Micmen,start=list(Vm=3.7,K=0.02))
What is Vm? What is K?
What's in an object of class "nls" ?
There is a summary and a plot function! What other diagnostic plots would be useful?
There are also functions resid, fitted, predict, and AIC function.
Data from a biochemical experiment where the intial velocity of a reaction measured for different concentations of the substrate are given in the data frame Puromycin . For more information use
> help(Puromycin)
Plot velocity against concentration separately for the two levels of treatment.
> plot(Puromycin$conc, Puromycin$rate, type="n")
> text(Puromycin$conc, Puromycin$rate, ifelse(Puromycin$state== "treated", "T", "U"))
Fit nonlinear model to treated data.
> Treated <- Puromycin[Puromycin$state == "treated",]
> Purfit1 <- nls(rate ~ Vm*conc/(K + conc), Treated, list(Vm = 200, K= 0.1), trace=T)
What does the trace=T option do?
Use the summary command and make some diagonistic plots for the model fit.
What is nls.control? To see, use the help function.
Weight versus age of chicks on different diets For more information use
> help(ChickWeight)
What does the coplot function do? Copy and paste from the (ChickWeight) help file to see.
Plot chick weight over time by chicks.
> plot(ChickWeight$Time, ChickWeight$weight, type="n")
> text(ChickWeight$Time, ChickWeight$weight, ChickWeight$Chick)
Let's pick out information for Chick #1, a representative Chick
> Chick.1 <- ChickWeight[ChickWeight$Chick == 1, ]
> lines(Chick.1$Time, Chick.1$weight)
Let's fit a logistic growth model:
> fit1 <- nls(weight~Asym/(1+exp((xmid-Time)/scal)),Chick.1,start=list(Asym=368,xmid=14,scal=6),trace=T)
What happens if you change scal=6 to scal=0.1?
What happens if you change xmid=14 to xmid=0.1?
How did I come up with starting values? There is a getInitial and SSlogis function in R.
> getInitial(weight ~ SSlogis(Time, Asym, xmid, scal), data = Chick.1)
Note: Initial values are in fact the converged values.
Here is a place to start: myplotnls1.r
Don't look at this one, yet! Try to make your own, first.
Here is mine: myplotnls.r (Feel free to use this one!)