Hi,
I am trying to use piecewise linear regression to approximate a
nonlinear function. Actually, I don't know how many linear functions I
need, therefore, I want build an array of regression models to automate
the approximation job. Could you please give me any clue?
Attached is ongoing code:
rawData = scan("c:/zyang/mass/data/A01/1.PRN",
what=list(numeric(),numeric()));
len = length(rawData[[1]]);
cuts = len*c(0.01, 0.03, 0.08, 0.18, 0.38, 0.69, 1);
cuts = as.integer(cuts);
mod1 = lm(rawData[[2]][1:cuts[1]]~rawData[[1]][1:cuts[1]]);
mod2 lm(rawData[[2]][cuts[1]:cuts[2]]~rawData[[1]][cuts[1]:cuts[2]]);
mod3 lm(rawData[[2]][cuts[2]:cuts[3]]~rawData[[1]][cuts[2]:cuts[3]]);
mod4 lm(rawData[[2]][cuts[3]:cuts[4]]~rawData[[1]][cuts[3]:cuts[4]]);
mod5 lm(rawData[[2]][cuts[4]:cuts[5]]~rawData[[1]][cuts[4]:cuts[5]]);
mod6 lm(rawData[[2]][cuts[5]:cuts[6]]~rawData[[1]][cuts[5]:cuts[6]]);
mod7 lm(rawData[[2]][cuts[6]:cuts[7]]~rawData[[1]][cuts[6]:cuts[7]]);
plot(rawData[[1]],rawData[[2]],type='l', col="green",
xlab="Da/z",
ylab="m/z");
abline(mod1, lty="1", col="red");
abline(mod2, lty="1", col="red");
abline(mod3, lty="1", col="red");
abline(mod4, lty="1", col="red");
abline(mod5, lty="1", col="red");
abline(mod6, lty="1", col="red");
abline(mod7, lty="1", col="red");
Thanks in advance,
On Wed, Dec 18, 2002 at 03:51:47PM -0500, Zhongming Yang wrote:> I am trying to use piecewise linear regression to approximate a > nonlinear function.Why not smooth regression, or non-linear regression?> Actually, I don't know how many linear functions I > need, therefore, I want build an array of regression models to automate > the approximation job. Could you please give me any clue?Clue 1) See above. You might be using the wrong tool. A smooth regression might be better here. help(loess), library(gss), and library(sm) are your friends. Clue 2) If you really want piecewise linear, a list makes more sense than a vector. R does handle vectors quite nicely, but I find its real strength is the way it handles complex lists with ease.> Attached is ongoing code: > > rawData = scan("c:/zyang/mass/data/A01/1.PRN", > what=list(numeric(),numeric())); > len = length(rawData[[1]]); > cuts = len*c(0.01, 0.03, 0.08, 0.18, 0.38, 0.69, 1); > cuts = as.integer(cuts);#change cuts to a matrix of values, col 1 is the lower #bound, col 2 is the upper bound for your segments. cuts <- cbind(c(1,cuts[1:(length(cuts)-1)], cuts) #make an empty list mod.list <- list() #fill that list with models for(ii in 1:dim(cuts)[1]) { start <- cuts[ii,1] end <- cuts[ii,2] mod.list[[ii]] <- lm(rawData[[2]][start,end] ~ rawData[[1]][start,end]) } #to extract coefficients lapply(mod.list,coef) #to extract coefficients, and confidence intervals lapply(mod.list,function(x,...){ coef(summary(x))} ) #to reproduce your ablines lapply(mod.list,abline,col="red") etc Cheers Jason -- Indigo Industrial Controls Ltd. 64-21-343-545 jasont at indigoindustrial.co.nz
And if you *really* want piecewise linear function (and most likely you want the pieces to be continuous, no?), there are better ways than yours. For "manual" fitting, use something like: library(splines) lm(y ~ bs(x, knots=..., deg=1)) For more automatic fitting, I believe bruto() or even mars() in the package `mda' will do. Andy> -----Original Message----- > From: Jason Turner [mailto:jasont at indigoindustrial.co.nz] > Sent: Wednesday, December 18, 2002 5:08 PM > To: Zhongming Yang > Cc: r-help at stat.math.ethz.ch > Subject: Re: [R] Can I build an array of regrssion model? > > > On Wed, Dec 18, 2002 at 03:51:47PM -0500, Zhongming Yang wrote: > > I am trying to use piecewise linear regression to approximate a > > nonlinear function. > > Why not smooth regression, or non-linear regression? > > > Actually, I don't know how many linear functions I > > need, therefore, I want build an array of regression models > to automate > > the approximation job. Could you please give me any clue? > > Clue 1) See above. You might be using the wrong tool. A smooth > regression might be better here. help(loess), library(gss), and > library(sm) are your friends. > > Clue 2) If you really want piecewise linear, a list makes more > sense than a vector. R does handle vectors quite nicely, but I > find its real strength is the way it handles complex lists with > ease. > > > Attached is ongoing code: > > > > rawData = scan("c:/zyang/mass/data/A01/1.PRN", > > what=list(numeric(),numeric())); > > len = length(rawData[[1]]); > > cuts = len*c(0.01, 0.03, 0.08, 0.18, 0.38, 0.69, 1); > > cuts = as.integer(cuts); > > #change cuts to a matrix of values, col 1 is the lower > #bound, col 2 is the upper bound for your segments. > > cuts <- cbind(c(1,cuts[1:(length(cuts)-1)], cuts) > > #make an empty list > mod.list <- list() > #fill that list with models > for(ii in 1:dim(cuts)[1]) { > start <- cuts[ii,1] > end <- cuts[ii,2] > mod.list[[ii]] <- lm(rawData[[2]][start,end] ~ > rawData[[1]][start,end]) > } > > #to extract coefficients > lapply(mod.list,coef) > > #to extract coefficients, and confidence intervals > lapply(mod.list,function(x,...){ coef(summary(x))} ) > > #to reproduce your ablines > lapply(mod.list,abline,col="red") > > etc > > Cheers > > Jason > -- > Indigo Industrial Controls Ltd. > 64-21-343-545 > jasont at indigoindustrial.co.nz > > ______________________________________________ > R-help at stat.math.ethz.ch mailing list > http://www.stat.math.ethz.ch/mailman/listinfo/r-help >------------------------------------------------------------------------------
Thank you for your time.
Actually, the purpose of this program is not smoothing but rather to
pick some peaks from the noisy exponential decay data. Piecewise linear
regression is ued to find the baseline. I don't know whether there is
other technique to solve this problem.
Could you please check the following code? The concern is how to draw
the raw data and piecewise linear functions on one diagram.
rawData = scan("c:/zyang/mass/data/A01/1.PRN",
what=list(numeric(),numeric()));
#write.table(rawData, 'rawdata.txt', quote=FALSE, row.names=FALSE,
col.names=FALSE);
len = length(rawData[[1]]);
cuts = len*c(0.01, 0.03, 0.08, 0.18, 0.38, 0.69, 1);
cuts = as.integer(cuts);
#make an empty list
mod.list <- list();
for(i in 1:length(cuts))
{
if (i==1)
{
start = 1;
}
else
{
start = cuts[i-1];
}
end = cuts[i];
mod.list[[i]] = lm(rawData[[2]][start:end] ~
rawData[[1]][start:end])
}
plot(rawData[[1]],rawData[[2]],type='l', col="green",
xlab="Da/z",
ylab="m/z");
for(i in 1:length(cuts))
{
if (i==1)
{
start = 1;
}
else
{
start = cuts[i-1];
}
end = cuts[i];
fitted.value[start:end] = predict(mod.list[[i]],
newdata=rawData[[1]][start:end], type='response' );
}
abline(rawData[[1]],fitted.value, col="red");
Thanks,
>>> Jason Turner <jasont at indigoindustrial.co.nz> 12/18/02
05:08PM >>>
On Wed, Dec 18, 2002 at 03:51:47PM -0500, Zhongming Yang
wrote:> I am trying to use piecewise linear regression to approximate a
> nonlinear function.
Why not smooth regression, or non-linear regression?
> Actually, I don't know how many linear functions I
> need, therefore, I want build an array of regression models to
automate> the approximation job. Could you please give me any clue?
Clue 1) See above. You might be using the wrong tool. A smooth
regression might be better here. help(loess), library(gss), and
library(sm) are your friends.
Clue 2) If you really want piecewise linear, a list makes more
sense than a vector. R does handle vectors quite nicely, but I
find its real strength is the way it handles complex lists with
ease.
> Attached is ongoing code:
>
> rawData = scan("c:/zyang/mass/data/A01/1.PRN",
> what=list(numeric(),numeric()));
> len = length(rawData[[1]]);
> cuts = len*c(0.01, 0.03, 0.08, 0.18, 0.38, 0.69, 1);
> cuts = as.integer(cuts);
#change cuts to a matrix of values, col 1 is the lower
#bound, col 2 is the upper bound for your segments.
cuts <- cbind(c(1,cuts[1:(length(cuts)-1)], cuts)
#make an empty list
mod.list <- list()
#fill that list with models
for(ii in 1:dim(cuts)[1]) {
start <- cuts[ii,1]
end <- cuts[ii,2]
mod.list[[ii]] <- lm(rawData[[2]][start,end] ~
rawData[[1]][start,end])
}
#to extract coefficients
lapply(mod.list,coef)
#to extract coefficients, and confidence intervals
lapply(mod.list,function(x,...){ coef(summary(x))} )
#to reproduce your ablines
lapply(mod.list,abline,col="red")
etc
Cheers
Jason
--
Indigo Industrial Controls Ltd.
64-21-343-545
jasont at indigoindustrial.co.nz
Thanks,
But why I can't draw regression line with the following code:
rawData = scan("c:/zyang/mass/data/A01/1.PRN",
what=list(numeric(),numeric()));
len = length(rawData[[1]]);
mod = lm(rawData[[2]]~bs(rawData[[1]],10,degree=1));
plot(rawData[[1]],rawData[[2]],type='l', col="green",
xlab="Da/z",
ylab="m/z");
abline(mod,col="red");
>>> Thomas Lumley <tlumley at u.washington.edu> 12/19/02 11:04AM
>>>
On Thu, 19 Dec 2002, Jason Turner wrote:
> On Wed, Dec 18, 2002 at 03:51:47PM -0500, Zhongming Yang wrote:
> > I am trying to use piecewise linear regression to approximate a
> > nonlinear function.
>
> Why not smooth regression, or non-linear regression?
>
> > Actually, I don't know how many linear functions I
> > need, therefore, I want build an array of regression models to
automate> > the approximation job. Could you please give me any clue?
>
> Clue 1) See above. You might be using the wrong tool. A smooth
> regression might be better here. help(loess), library(gss), and
> library(sm) are your friends.
>
> Clue 2) If you really want piecewise linear, a list makes more
> sense than a vector. R does handle vectors quite nicely, but I
> find its real strength is the way it handles complex lists with
> ease.
I don't see any problem with wanting to fit linear splines. It's
quite
easy, as well
eg
models <- lapply( 1:8, function(n) lm(y~bs(x, n, degree=1)))
fits piecewise linear functions with 1 to 8 pieces.
-thomas
______________________________________________
R-help at stat.math.ethz.ch mailing list
http://www.stat.math.ethz.ch/mailman/listinfo/r-help
On Thu, 19 Dec 2002, Zhongming Yang wrote:> Thanks, > > But why I can't draw regression line with the following code: > > rawData = scan("c:/zyang/mass/data/A01/1.PRN", > what=list(numeric(),numeric())); > len = length(rawData[[1]]); > mod = lm(rawData[[2]]~bs(rawData[[1]],10,degree=1)); > plot(rawData[[1]],rawData[[2]],type='l', col="green", xlab="Da/z", > ylab="m/z"); > abline(mod,col="red"); >Because that's not what abline does. It draws a straight line, not a piecewise linear regression curve. You might find termplot() useful. -thomas
Try lines(predict(mod), col="red"). (This is assuming your original x
data
are sorted.)
abline(mod) draws a single straight line, using intercept and slope from a
simple linear regression model ('mod'). If mod has more than one slope
(as
is the case with bs()), it won't make sense.
Andy
> -----Original Message-----
> From: Zhongming Yang [mailto:Zhongming.Yang at cchmc.org]
> Sent: Thursday, December 19, 2002 11:21 AM
> To: jasont at indigoindustrial.co.nz; tlumley at u.washington.edu
> Cc: r-help at stat.math.ethz.ch
> Subject: Re: [R] Can I build an array of regrssion model?
>
>
> Thanks,
>
> But why I can't draw regression line with the following code:
>
> rawData = scan("c:/zyang/mass/data/A01/1.PRN",
> what=list(numeric(),numeric()));
> len = length(rawData[[1]]);
> mod = lm(rawData[[2]]~bs(rawData[[1]],10,degree=1));
> plot(rawData[[1]],rawData[[2]],type='l', col="green",
xlab="Da/z",
> ylab="m/z");
> abline(mod,col="red");
>
>
>
>
>
>
>
> >>> Thomas Lumley <tlumley at u.washington.edu> 12/19/02
11:04AM >>>
> On Thu, 19 Dec 2002, Jason Turner wrote:
>
> > On Wed, Dec 18, 2002 at 03:51:47PM -0500, Zhongming Yang wrote:
> > > I am trying to use piecewise linear regression to approximate a
> > > nonlinear function.
> >
> > Why not smooth regression, or non-linear regression?
> >
> > > Actually, I don't know how many linear functions I
> > > need, therefore, I want build an array of regression models to
> automate
> > > the approximation job. Could you please give me any clue?
> >
> > Clue 1) See above. You might be using the wrong tool. A smooth
> > regression might be better here. help(loess), library(gss), and
> > library(sm) are your friends.
> >
> > Clue 2) If you really want piecewise linear, a list makes more
> > sense than a vector. R does handle vectors quite nicely, but I
> > find its real strength is the way it handles complex lists with
> > ease.
>
>
> I don't see any problem with wanting to fit linear splines. It's
> quite
> easy, as well
> eg
> models <- lapply( 1:8, function(n) lm(y~bs(x, n, degree=1)))
>
> fits piecewise linear functions with 1 to 8 pieces.
>
> -thomas
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> http://www.stat.math.ethz.ch/mailman/listinfo/r-help
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> http://www.stat.math.ethz.ch/mailman/listinfo/r-help
>
------------------------------------------------------------------------------