rprot_ex.R

library(magrittr)
library(ggplot2)
library(rpx)

rpx:::apply_fix_issue_5(FALSE)

## https://www.ebi.ac.uk/pride/archive/projects/PXD022816
## RawBeans: A Simple, Vendor-Independent, Raw-Data Quality-Control
## Tool (10.1021/acs.jproteome.0c00956)

PXD022816 <- PXDataset("PXD022816")
PXD022816
pxtax(PXD022816)
pxfiles(PXD022816)

pxurl(PXD022816)
(mzids <- pxget(PXD022816, grep("mzID", pxfiles(PXD022816))[1:3]))
(mzmls <- pxget(PXD022816, grep("mzML", pxfiles(PXD022816))[1:3]))

library(PSM)
library(Spectra)

sp <- Spectra(mzmls) %>%
    filterMsLevel(2)
sp

table(basename(sp$dataOrigin))


id <- readPSMs(mzids)
table(id$idFile)

tidyr::as_tibble(id) %>%
    ggplot(aes(x = MetaMorpheus.score,
               colour = isDecoy)) +
    geom_density() +
    facet_wrap(~ spectrumFile)

id_filtered <- filterPSMs(id)
summary(id_filtered$PSM.level.q.value)



sp$by <-
    paste(sub("207e04bccd3dd_", "", basename(dataOrigin(sp))),
          sub("^.+=1 ", "", sp$spectrumId),
          sep = "_")
id$by <-
    paste(sub("^.+\\QEP", "QEP", id$spectrumFile),
          sub("^.+=1 ", "", id$spectrumID),
          sep = "_")

sp <- joinSpectraData(sp, id, by.x = "by")

table(is.na(sp$sequence))

sp <- sp[!is.na(sp$sequence)]

k <- sp$sequence %in%
    c("DGSDEPGTAACPNGSFHCTNTGYK",
      "DGQVINETSQHHDDLE",
      "DSYVGDEAQSK")

sp_k <- sp[k]
sp_k <- setBackend(sp_k, MsBackendDataFrame())
## normalise intensities
sp_k$intensity <- sp_k$intensity/max(sp_k$intensity)

cmat <- compareSpectra(sp_k)

rownames(cmat) <-
    colnames(cmat) <- strtrim(sp_k$sequence, 3)

pheatmap::pheatmap(cmat)

plotSpectra(sp_k[1:4], labels = addFragments,
            labelPos = 3, labelCol = "steelblue",
            main = sp_k$sequence[1:4])

plotSpectraMirror(sp_k[1], sp_k[2])
plotSpectraMirror(sp_k[1], sp_k[3])
plotSpectraMirror(sp_k[1], sp_k[4])

plotSpectraMirror(sp_k[5], sp_k[6])
plotSpectraMirror(sp_k[5], sp_k[10])

plotSpectraMirror(sp_k[7], sp_k[8])
plotSpectraMirror(sp_k[7], sp_k[9])