R: Fit chromatographic peaks with a gaussian profile
fitpeaks
R Documentation
Fit chromatographic peaks with a gaussian profile
Description
Find chromatographic peaks, and fit peak parameters using a gaussian
profile. The algorithm is extremely simple and could be replaced by a
more sophisticated algorithm. In particular one can expect bad fits if
peaks are overlapping significantly.
Usage
findpeaks(y, span = NULL)
fitpeaks(y, pos)
Arguments
y
response (numerical vector)
span
number of points used in the definition of what
constitutes a "local" maximum. If not given, a default value of 20
percent of the number of time points is used.
pos
locations of local maxima in vector y
Details
Finding peaks with function findpeaks is based on the position
of local maxima within a window of width span.
Peak parameters are calculated using fitpeaks, assuming a
normal distribution. Peak width is given as a standard deviation,
calculated from the full width at half maximum (FWHM); the peak area
is given by the ratio of the peak height and the density.
Value
Function findpeaks simply returns the locations of the local
maxima, expressed as indices.
Function fitpeaks returns a matrix, whose columns contain the
following information:
rt
location of the maximum of the peak (x)
sd
width of the peak (x)
FWHM
full width at half maximum (x)
height
height of the peak (y)
area
peak area
Again, the first three elements (rt, sd and FWHM) are expressed as
indices, so not in terms of the real retention times. The
transformation to "real" time is done in function getAllPeaks.
Note
Function findpeaks was modelled after code suggested by
Brian Ripley on the R help list.
Author(s)
Ron Wehrens
See Also
getAllPeaks
Examples
data(tea)
new.lambdas <- seq(260, 500, by = 2)
tea <- lapply(tea.raw, preprocess, dim2 = new.lambdas)
tea.split <- splitTimeWindow(tea, c(12, 14), overlap = 10)
Xl <- tea.split[[2]]
Xl.opa <- opa(Xl, 4)
Xl.als <- doALS(Xl, Xl.opa)
tpoints <- getTime(Xl.als)
plot(tpoints, Xl.als$CList[[2]][,2], type = "l", col = "gray")
pk.pos <- findpeaks(Xl.als$CList[[2]][,2], span = 11)
abline(v = tpoints[pk.pos], col = 4)
pks <- fitpeaks(Xl.als$CList[[2]][,2], pk.pos)
apply(pks, 1,
function(pkmodel) {
lines(tpoints,
dnorm(1:length(tpoints), pkmodel["rt"], pkmodel["sd"]) *
pkmodel["area"],
col = 2)
invisible()
})
## reasonably close fit, apart from the small peak in the middle...
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(alsace)
Loading required package: ALS
Loading required package: nnls
Loading required package: Iso
Iso 0.0-17
Loading required package: ptw
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/alsace/fitpeaks.Rd_%03d_medium.png", width=480, height=480)
> ### Name: fitpeaks
> ### Title: Fit chromatographic peaks with a gaussian profile
> ### Aliases: findpeaks fitpeaks
> ### Keywords: manip
>
> ### ** Examples
>
> data(tea)
> new.lambdas <- seq(260, 500, by = 2)
> tea <- lapply(tea.raw, preprocess, dim2 = new.lambdas)
> tea.split <- splitTimeWindow(tea, c(12, 14), overlap = 10)
>
> Xl <- tea.split[[2]]
> Xl.opa <- opa(Xl, 4)
>
> Xl.als <- doALS(Xl, Xl.opa)
>
> tpoints <- getTime(Xl.als)
> plot(tpoints, Xl.als$CList[[2]][,2], type = "l", col = "gray")
> pk.pos <- findpeaks(Xl.als$CList[[2]][,2], span = 11)
> abline(v = tpoints[pk.pos], col = 4)
>
> pks <- fitpeaks(Xl.als$CList[[2]][,2], pk.pos)
> apply(pks, 1,
+ function(pkmodel) {
+ lines(tpoints,
+ dnorm(1:length(tpoints), pkmodel["rt"], pkmodel["sd"]) *
+ pkmodel["area"],
+ col = 2)
+ invisible()
+ })
NULL
> ## reasonably close fit, apart from the small peak in the middle...
>
>
>
>
>
> dev.off()
null device
1
>