Re: [R] Fitting data from a spectrophotometer.

About this list Date view Thread view Subject view Author view Attachment view

From: Sundar Dorai-Raj (sundar.dorai-raj@pdf.com)
Date: Wed 12 May 2004 - 06:43:23 EST


Message-id: <40A13AEB.90407@pdf.com>


TAPO (Thomas Agersten Poulsen) wrote:
> Dear R-list,
>
> It is not uncommon for laboratory equipment (e.g.
spectrophotometers) to have a linear response in a certain interval and
then go into saturation. I wonder if there is an R-function that models
this; for instance by estimating the breakpoint and fitting a line below
the breakpoint and a constant above.
>
> Best regards
> Thomas Poulsen
>
> --
> Thomas Poulsen Research Scientist. PhD, MSc
> Novozymes A/S Protein Design / Bioinformatics
> Brudelysvej 26, 1US.24 Phone: +45 44 42 27 23
> DK-2880 Bagsværd.
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

I don't know if there is something this specific out there already, but
you could use ?optim to do this. Here's a quick example:

set.seed(1)
x <- rnorm(100)
y <- ifelse(x < 0, 1 + 4 * x, 1) + rnorm(100, sd = 0.2)

brk <- function(par, x, y) {
   a <- par[1] # intercept
   b <- par[2] # slope
   p <- par[3] # break point
   h <- par[4] # height of saturation
   yhat <- ifelse(x < p, a + b * x, h)
   sum((y - yhat)^2)
}

st <- coef(lm(y ~ x)) # starting values
fit <- optim(c(st, 0, 0), brk, x = x, y = y)
a <- fit$par[1]
b <- fit$par[2]
p <- fit$par[3]
h <- fit$par[4]

yhat <- ifelse(x < p, a + b * x, h)
plot(y ~ x)
lines(x[order(x)], yhat[order(x)])

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


About this list Date view Thread view Subject view Author view Attachment view

This archive was generated by hypermail 2.1.3 : Mon 31 May 2004 - 23:05:09 EST