Metabolomics Society Forum

## Load mzR library
library(mzR)

## Open the mzML file that you want to convert
dat<-openMSfile('myMSdata.mzML')

## Get spectra data
pks <- spectra(dat)
## Get file header
hdr <- header(dat)

## Remove all scans but ms2 scans
pks<-pks[hdr$msLevel==2]
hdr<-hdr[hdr$msLevel==2,]

## Take a quick look and make sure everything looks ok
head(hdr)

## Provide new scan/acquisition numbers
hdr$seqNum<-hdr$acquisitionNum<-seq(nrow(hdr))

## Clear all precursor charge/Intensity/MZ/ScanNum columns
## Not sure if this is required and how precursor ion scans will differ
hdr$precursorCharge<-hdr$precursorIntensity<-hdr$precursorMZ<-hdr$precursorScanNum<-0

## Overwrite msLevel to 'pretend' to be MS1 data
hdr$msLevel<-1

## Write out the new 'MS1' data
writeMSData(pks,'myAlteredMSdata.mzML',header=hdr)

mass <- seq(peaks[1, "mz"], peaks[nrow(peaks), "mz"] + binSize, by = binSize / 2)

mz<-c(153.03,153.06)
rtr_1<-c(155, 175)
binSize<-0.01

## Get peaks that you want to find in a feature
specificPeaks<-chromPeaks(xdata,rt=rtr_1,mz=mzr)

## Plot the peak locations
plot(specificPeaks[,"rt"],specificPeaks[,"mz"])

## Create approximate grouping windows (not the values that 'groupChromPeaks' will generate)
mass<-seq(min(specificPeaks[, "mz"]), max(specificPeaks[, "mz"]) + binSize, by = binSize / 2)
abline(h=mass,col=rep(c('red','green')))
## Note these are overlapping windows. So look at groupings between the coloured lines

## Save peak and feature definitions
Pks<-chromPeaks(xdata)
fDef<-featureDefinitions(xdata)

## Get peaks that you want to find in a feature
specificPeaks<-chromPeaks(xdata,rt=rtr_1,mz=mzr,type='all',msLevel=1L)

## Identify the row number of each peak from Pks
rows<-which(row.names(Pks)%in%row.names(specificPeaks))

## Search each feature definition peak list and identify the feature/s that contain every peak from above
fDef_row<-which(sapply(fDef$peakidx,function(x) all(rows%in%x)))

## Show the identified feature/s
fDef[fDef_row,]

# Read data file
raw_data<-readMSData(files = files[1], pdata = new("NAnnotatedDataFrame", pd), msLevel. = 1)
  
# Extract out spectra data
spec<-spectra(raw_data)

# Iterate through the spectra list and obtain relevant information
specl<-lapply(spec,function(x) c(x@msLevel,x@rt,x@scanIndex))

# Compact list to matrix
specl<-do.call(rbind,specl)

devtools::install_github("sneumann/xcms", ref = "RELEASE_3_8")

plotXIC_MSnExp_ex <- function(x, ...) {
  ## Restrict to MS level 1
  x <- filterMsLevel(x, 1L)
  if (!length(x))
    stop("No MS1 data available")
  fns <- basename(fileNames(x))
  if (isMSnbaseVerbose())
    message("Retrieving data ...", appendLF = FALSE)
  x <- as(x, "data.frame")
  x <- split(x, x$file)
  if (isMSnbaseVerbose())
    message("OK")
  ## Check if we are greedy and plot a too large area
  if (any(unlist(lapply(x, nrow)) > 20000))
    warning("The MS area to be plotted seems rather large. It is suggested",
            " to restrict the data first using 'filterRt' and 'filterMz'. ",
            "See also ?chromatogram and ?Chromatogram for more efficient ",
            "functions to plot a total ion chromatogram or base peak ",
            "chromatogram.",
            immediate = TRUE, call = FALSE)
  ## Define the layout.
  dots <- list(...)
  if (any(names(dots) == "layout")) {
    if (!is.null(dots$layout))
      layout(layout)
    dots$layout <- NULL
  } else
    layout(.vertical_sub_layout_ex(length(x)))
  tmp <- mapply(x, fns, FUN = function(z, main, ...) {
    .plotXIC_ex(x = z, main = main, layout = NULL, ...)
  }, MoreArgs = dots)
}


.plotXIC_ex <- function(x, main = "", col = "grey", colramp = topo.colors,
                     grid.color = "lightgrey", pch = 21,
                     layout = matrix(1:2, ncol = 1), ...) {
  # start edit
  plot_strip<-function(..., v, h) plot(...)
  dots <- list(...)
  # end edit
  
  print(list(...))
  if (is.matrix(layout))
    layout(layout)
  ## Chromatogram.
  bpi <- unlist(lapply(split(x$i, x$rt), max, na.rm = TRUE))
  brks <- lattice::do.breaks(range(x$i), nint = 256)
  par(mar = c(0, 4, 2, 1))
  
  # start edit
  plot_strip(as.numeric(names(bpi)), bpi, xaxt = "n", col = col, main = main, 
       bg = lattice::level.colors(bpi, at = brks, col.regions = colramp), xlab = "",
       pch = pch, ylab = "", las = 2, ...)
  # end edit
  
  if(!is.null(dots$v)) abline(v=dots$v,col='red',lty=2)
  
  mtext(side = 4, line = 0, "Intensity", cex = par("cex.lab"))
  grid(col = grid.color)
  par(mar = c(3.5, 4, 0, 1))
  
  # start edit
  plot_strip(x$rt, x$mz, main = "", pch = pch, col = col, xlab = "", ylab = "",
       yaxt = "n", bg = lattice::level.colors(x$i, at = brks, col.regions = colramp),
       ...)
  
  if(!is.null(dots$h)) abline(h=dots$h,col='red',lty=2)
  if(!is.null(dots$v)) abline(v=dots$v,col='red',lty=2)
  # end edit
  
  axis(side = 2, las = 2)
  grid(col = grid.color)
  mtext(side = 1, line = 2.5, "Retention time", cex = par("cex.lab"))
  mtext(side = 4, line = 0, "m/z", cex = par("cex.lab"))
}

.vertical_sub_layout_ex <- function(x, sub_plot = 2) {
  sqrt_x <- sqrt(x)
  ncol <- ceiling(sqrt_x)
  nrow <- round(sqrt_x)
  rws <- split(1:(ncol * nrow * sub_plot), f = rep(1:nrow,
                                                   each = sub_plot * ncol))
  do.call(rbind, lapply(rws, matrix, ncol = ncol))
}