R help - Mar 2013 - multiple peak fit

Hi

I went through some extensive search to find suitable method (package, function)
to fit multiple peaks. The best I found is ALS package but it requires rather
complicated input structure probably resulting from GC-MS experimental data and
seems to be an overkill to my problem.

I have basically simple two column data frame with columns time and sig 

dput(temp)
structure(list(time = c(33, 34, 37.33333, 44.66667, 56, 70, 85, 
99.66667, 114, 127.33333, 140, 151.66667, 163, 174, 185, 196, 
206.66667, 217.33333, 228, 238.66667, 249.33333, 260, 271, 282, 
292.66667, 303.66667, 314.33333, 325.33333, 336, 346.66667, 357.33333, 
368, 379, 389.66667, 400.33333, 411, 421.66667, 432.33333, 443, 
454, 465, 476, 487, 497.66667, 508.33333, 519, 530, 541, 551.66667, 
562.33333, 573, 584, 595, 606, 616.66667, 627.33333, 638, 649, 
659.66667, 670.33333, 681, 692, 703, 713.66667, 724.33333, 735, 
746, 757, 768, 779, 789), sig = c(1.558, 1.5549, 1.5619, 1.5614, 
1.5618, 1.6044, 1.6161, 1.6287, 1.6432, 1.6925, 1.7273, 1.6932, 
1.669, 1.6863, 1.6962, 1.7186, 1.7513, 1.8325, 1.9181, 2.01, 
2.1276, 2.2821, 2.5596, 2.9844, 4.1272, 13.421, 15.422, 14.119, 
11.491, 8.8799, 6.7774, 5.6223, 4.8775, 4.3628, 4.0517, 3.9146, 
3.8704, 3.9162, 4.0372, 4.0948, 4.2054, 4.1221, 3.9145, 3.5724, 
3.2108, 2.8311, 2.4605, 2.1985, 1.9685, 1.8158, 1.7487, 1.692, 
1.6565, 1.6374, 1.609, 1.5927, 1.5401, 1.5366, 1.5614, 1.5314, 
1.4989, 1.5053, 1.4953, 1.4824, 1.4743, 1.4468, 1.4698, 1.4671, 
1.4675, 1.4704, 1.4966)), .Names = c("time", "sig"),
row.names = c(NA,
-71L), class = "data.frame")

When it is plotted there are clearly 2 peaks visible.

plot(temp, type="l")

Does anybody know how to decompose those data to two (or more) gaussian (or
other) peaks?

I can only think about nlme but before I start from scratch I try to ask others.

Best regards
Petr

Petr,

You could use nonlinear regression to fit two Guassian peaks.  For example,

fit <- nls(sig ~ d + a1*exp(-0.5*((time-c1)/b1)^2) +
a2*exp(-0.5*((time-c2)/b2)^2),
 start=list(a1=15, b1=20, c1=315, a2=4, b2=50, c2=465, d=1.5), data=temp)
plot(temp, type="l")
lines(temp$time, predict(fit), lty=2, col="red")

Jean



On Thu, Mar 21, 2013 at 8:17 AM, PIKAL Petr <petr.pikal@precheza.cz>
wrote:
> Hi
>
> I went through some extensive search to find suitable method (package,
> function) to fit multiple peaks. The best I found is ALS package but it
> requires rather complicated input structure probably resulting from GC-MS
> experimental data and seems to be an overkill to my problem.
>
> I have basically simple two column data frame with columns time and sig
>
> dput(temp)
> structure(list(time = c(33, 34, 37.33333, 44.66667, 56, 70, 85,
> 99.66667, 114, 127.33333, 140, 151.66667, 163, 174, 185, 196,
> 206.66667, 217.33333, 228, 238.66667, 249.33333, 260, 271, 282,
> 292.66667, 303.66667, 314.33333, 325.33333, 336, 346.66667, 357.33333,
> 368, 379, 389.66667, 400.33333, 411, 421.66667, 432.33333, 443,
> 454, 465, 476, 487, 497.66667, 508.33333, 519, 530, 541, 551.66667,
> 562.33333, 573, 584, 595, 606, 616.66667, 627.33333, 638, 649,
> 659.66667, 670.33333, 681, 692, 703, 713.66667, 724.33333, 735,
> 746, 757, 768, 779, 789), sig = c(1.558, 1.5549, 1.5619, 1.5614,
> 1.5618, 1.6044, 1.6161, 1.6287, 1.6432, 1.6925, 1.7273, 1.6932,
> 1.669, 1.6863, 1.6962, 1.7186, 1.7513, 1.8325, 1.9181, 2.01,
> 2.1276, 2.2821, 2.5596, 2.9844, 4.1272, 13.421, 15.422, 14.119,
> 11.491, 8.8799, 6.7774, 5.6223, 4.8775, 4.3628, 4.0517, 3.9146,
> 3.8704, 3.9162, 4.0372, 4.0948, 4.2054, 4.1221, 3.9145, 3.5724,
> 3.2108, 2.8311, 2.4605, 2.1985, 1.9685, 1.8158, 1.7487, 1.692,
> 1.6565, 1.6374, 1.609, 1.5927, 1.5401, 1.5366, 1.5614, 1.5314,
> 1.4989, 1.5053, 1.4953, 1.4824, 1.4743, 1.4468, 1.4698, 1.4671,
> 1.4675, 1.4704, 1.4966)), .Names = c("time", "sig"),
row.names = c(NA,
> -71L), class = "data.frame")
>
> When it is plotted there are clearly 2 peaks visible.
>
> plot(temp, type="l")
>
> Does anybody know how to decompose those data to two (or more) gaussian
> (or other) peaks?
>
> I can only think about nlme but before I start from scratch I try to ask
> others.
>
> Best regards
> Petr
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
	[[alternative HTML version deleted]]