fillPeaks has zeros ?

August 19, 2011, 10:55:59 AM

fillPeaks() will add intensities for peaks not observed in a certain sample, but in some/most of the other samples of a smple class.
Some people have observed, that even after fillPeaks(), some values are still "0" or even "NaN" (NotANumber).
fillPeaks() is trying to (blindly) extract/integrate intensities from those raw files, where the peaks are "missing" after the peak picking step (matchedFilter or centWave).
The information, where to look, is taken from the other samples where the peak was observed, i.e. the rtmin/rtmax and mzmin/mzmax after group()ing. Two cases where this will fail are most prominent:

If the peak is near the beginning (or end) of the LC run, this time range could be "before" or "beyond" the runtime of the sample in question, especially if retcor() was used to correct/modify any RT shifts. In this case, fillPeaks() will report "NaN". There is nothing xcms could sensibly do otherwise, and this non-number should flag a warning that there was something fishy. You might workaround by removing the first / last seconds (depending on the gradient they usually don't contain important/reliable peaks anyway) using the scanrange=c(min,max) parameter to xcmsSet() and findPeaks().
If your m/z calibration is very good for most runs, the group will have a very small mzin/mzmax window. For centroided data it is then quite easy to miss the actual raw points by just a tiny bit, as seen in the graph. In that case, the sum of intensities will be zero.

[attachment=1:2z4nlwmv]300px-ZeroPeaks-26.png[/attachment:2z4nlwmv] raw data just above raw data fillPeaks looking "beyond" raw data: The red EIC ends, because the peak is located at the end of the run, and the retention time has been corrected previously.

[attachment=1:2z4nlwmv]300px-ZeroPeaks-26.png[/attachment:2z4nlwmv] fillPeaks looking "beyond" raw data: The red EIC ends, because the peak is located at the end of the run, and the retention time has been corrected previously.

This is the code used:

Code: [Select]

library(xcms)
load("Anal10909_SP.Rdata")
xs <- fill.ret.Anal10909_SP

# Shortcut: aliases for groups
gr <- groups(xs)
gv <- groupval(xs, value="into")

# Shortcut: get offending group(s):
nanidx <- which(is.na(gv), arr.ind=TRUE)
zeridx <- which(gv==0, arr.ind=TRUE)

# Plot offending group(s):
e <- getEIC(xs, groupidx=nanidx[1], )
# Bad NaN guy in red:
col=rep("grey", length(sampnames(xs)))
col[nanidx[2]] <- "red"
plot(e, xs, col=col, sleep=1)

Code: [Select]

## PlotRaws
xr <- xcmsRaw(filepaths(xs)[nanidx[2]])
rtcor <- xs@rt$corrected
xr@scantime <- rtcor[[nanidx[2]]]

pdf(file="zeroPeaks.pdf")
par(mfrow=c(2,1))

for (i in 1:length(zeridx[,1])) {
    # EICs
    e <- getEIC(xs, groupidx=zeridx[i,1], )
    plot(e, xs, col=as.character(col[zeridx[i,1],]), sleep=1)

    ## PlotRaws
    xr <- xcmsRaw(filepaths(xs)[zeridx[i,2]])
    currentGroup <- gr[zeridx[i,1], ]

    ## Disable warnings for empty regions

    suppressWarnings(plotRaw(xr,
                             massrange=currentGroup[c("mzmin","mzmax")]+c(-0.1,+0.1),
                             timerange=currentGroup[c("rtmin","rtmax")]+c(-5,+5),
                             log=TRUE,
                             title=paste("Group/Samp/mz/RT: ", paste(zeridx[i,], collapse=" "), paste(round(currentGroup[c("mzmin","mzmax","rtmin","rtmax")], digits=2), collapse=" "))))

    rect(currentGroup["rtmin"],
         currentGroup["mzmin"],
         currentGroup["rtmax"],
         currentGroup["mzmax"], border="grey")
}
dev.off()

[attachment deleted by admin]