filterMz,Spectrum not returning empty spectrum if mzlim outside of range

[jorainer]

If no values are found by the filterMz method for a given mzlim range the method warns about not finding data points for that range and returning an empty spectrum, it returns however a Spectrum with the mz corresponding to the mzlim and intensities of 0.

Wouldn't it be better to return a Spectrum with @mz and @intensity being numeric()?

Running example:

library(MSnbase)
mzf <- system.file("microtofq/MM14.mzML", package = "msdata")
od <- readMSData2(files = mzf, msLevel. = 1, centroided. = TRUE)
sp <- od[[1]]

empty_sp <- filterMz(sp, mz = c(1005, 1006))
Warning message:
In trimMz_Spectrum(object, mzlim = mz, msLevel. = msLevel., ...) :
  No data points between 1005 and 1006 for spectrum with acquisition
number 1. Returning empty spectrum.

## Spectrum is not so empty...
mz(empty_sp)
[1] 1005 1006

intensity(empty_sp)
[1] 0 0

[lgatto]

I don't remember if there was a good reason to implement it like this - may be @sgibb remembers.
One reason I could think of is that, after running filterMz, there is a guarantee that there is some data (possibly just two 0s; and these spectra's indices are reported in the warning message anyway), rather that having one or a few spectra with nothing. This could just facilitate systematic processing.

May be this call for a filterEmptySpectra function. And not that the isEmpty method is not in line with the current filterMz issue you describe, as it looks for an empty mz slot.

> getMethod("isEmpty", "Spectrum")
Method Definition:

function (x) 
length(x@mz) == 0
<environment: namespace:MSnbase>

Signatures:
        x         
target  "Spectrum"
defined "Spectrum"

Let's see what Seb has to say, otherwise, I would be in favour as long as we add an filterEmptySpectra,MSnExp method.

[jorainer]

Understand. It's no problem, I was just puzzled that I got some values although it was saying it returns an empty spectrum.

Currently I'm playing around in xcms with a method to extract a MS data subset defined by an rt and mz range - this should be as fast as possible, thus I'm not really keen in calling another method filterEmptySpectra after filterRt and filterMz. My workaround for now is to check if the peaksCount of the spectrum is 0. That seems to work. Would be nice if isEmpty could get fixed, so I could use that instead.

[lgatto]

Let's wait for Seb. Having thought a bit about it, I would prefer to implement your suggestion.

Would calling filterEmptySpectra be really slow? You wouldn't even have too, unless having empty spectra causes problems later. Could you clarify.

[jorainer]

No, filterEmptySpectra would not be slow - I just prefer to not call too many methods in this function I'm implementing to keep it fast. Also, I'd like to know if there was an empty spectrum but I don't want to remove them. Best is I share the code:

The function is thought to replace the getEIC method in xcms in the long run which extracts an ion chromatogram for a specified rt and mz range. It has some drawbacks (works on binned data only etc), but is quite fast. Now, my function is a little more generic - I'm extracting all of the MS data for a mz/rt range grouped per file so I can generate the chromatogram or whatever later.

extractMsData <- function(x, rtrange, mzrange) {
    ## Subset the OnDiskMSnExp
    fns <- fileNames(x)
    tmp <- filterMz(filterRt(x, rt = rtrange), mz = mzrange)
    fromF <- match(fileNames(tmp), fns)
    suppressWarnings(
        dfs <- spectrapply(tmp, function(z) {
            if (peaksCount(z))
                return(data.frame(rt = rep_len(rtime(z), length(z@mz)),
                                  as.data.frame(z)))
            else
                return(data.frame(rt = numeric(), mz = numeric(),
                                  i = integer()))
        })
    )
    res <- vector(mode = "list", length = length(fns))
    res[fromF] <- split(dfs, f = fromFile(tmp))
    return(lapply(res, do.call, what = rbind))
}

[sgibb]

I would prefer to return an empty Spectrum with length(mz) == 0 && length(intensity) == 0 && peaksCount(x) == 0). I can't remember why filterMz is implemented in the current way but I am not aware that I ever touched this function at all.
By looking at this function i was wondering why we don't have an [ operator for the Spectrum class and its subclasses (at least internally it would be useful, we don't have to export it).

[jorainer]

Previously the filterMz was called trimMz, so eventually there was some reason for this. In any rate, I'm fine with any solution as long as there is a way to determine that the Spectrum is empty.

[lgatto]

I would suggest to fix filterMz to return an empty spectrum with no mz, intensity and add a filterEmptySpectra method. I am happy to do this once I get a chance.

[lgatto]

Added a filterEmptySpectra method in commit f1ee470