## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
    collapse = TRUE,
    comment = "#>",
    fig.width = 10
)

## ----include = FALSE, eval=FALSE----------------------------------------------
# library(devtools)
# build_vignettes()
# rmarkdown::render("Score-function-test-of-PSM.Rmd", output_dir = "../doc/")

## ----setup--------------------------------------------------------------------
library(Aerith)
library(dplyr)
library(stringr)
library(ggplot2)
library(tidyr)

## -----------------------------------------------------------------------------
demo_file <- system.file("extdata", "107728.FT2", package = "Aerith")
scan2 <- readOneScanMS2(ftFile = demo_file, 107728)
scan1 <- getRealScanFromList(scan2)
spectra <- slot(scan1, "spectra")
charges <- slot(scan1, "charges")
pep <- "HSQVFSTAEDNQSAVTIHVLQGER"
pep2 <- "[HSQVFSTAEDNQSAVTIHVLQGER]"

## -----------------------------------------------------------------------------
isoCenter <- scan2$isolationWindowCenterMZ
anno <- annotatePSM(
    spectra$MZ, spectra$Prob,
    spectra$Charge,
    pep, 1:2, "C13",
    0.0107, isoCenter, 5.0, TRUE
)
WDPscores <- sapply((0:100) / 100, simplify = TRUE, function(prob) {
    return(scorePSM(spectra$MZ,
        realIntensity = spectra$Prob,
        realCharge = spectra$Charge, parentCharge = charges[1],
        pepSeq = pep2, Atom = "C13", Prob = prob
    ))
})
scores <- sapply((0:100) / 100, simplify = TRUE, function(prob) {
    anno <- annotatePSM(
        spectra$MZ, spectra$Prob,
        spectra$Charge,
        pep, 1:2, "C13",
        prob, isoCenter, 4.0, TRUE
    )
    return(c(XcorrScore = anno$XcorrScore, MVHscore = anno$MVHscore))
})
scores <- t(scores)
scores <- cbind(WDPscores, scores)

## -----------------------------------------------------------------------------
scores_df <- as.data.frame(scores)
scores_df$abundance <- 0:100

# Reshape the data to long format for faceting
scores_long <- scores_df %>%
    pivot_longer(
        cols = c(WDPscores, XcorrScore, MVHscore),
        names_to = "score_type",
        values_to = "score_value"
    )

# Create a publication-quality plot
ggplot(scores_long, aes(x = abundance, y = score_value)) +
    geom_line() +
    facet_wrap(~score_type, scales = "free_y") +
    labs(
        x = expression(paste(~ {}^{
            13
        }, "C %")),
        y = "Score",
        title = expression(paste("Scores Across", ~ {}^{
            13
        }, "C %"))
    ) +
    theme_minimal() +
    theme(
        plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
        axis.title = element_text(size = 14),
        axis.text = element_text(size = 12),
        strip.text = element_text(size = 12, face = "bold"),
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "white", colour = "black")
    )

## -----------------------------------------------------------------------------
demo_file <- system.file("extdata", "X13_4068_2596_8182.FT2", package = "Aerith")
ft2 <- readAllScanMS2(demo_file)
scan2 <- ft2[["2596"]]
scan1 <- getRealScanFromList(scan2)
spectra <- slot(scan1, "spectra")
charges <- slot(scan1, "charges")
pep <- "HYAHVDCPGHADYVK"
pep2 <- "[HYAHVDCPGHADYVK]"
pct <- 0.52

## -----------------------------------------------------------------------------
isoCenter <- scan2$isolationWindowCenterMZ
WDPscores <- sapply((0:100) / 100, simplify = TRUE, function(prob) {
    return(scorePSM(spectra$MZ,
        realIntensity = spectra$Prob,
        realCharge = spectra$Charge, parentCharge = charges[1],
        pepSeq = pep2, Atom = "C13", Prob = prob
    ))
})
scores <- sapply((0:100) / 100, simplify = TRUE, function(prob) {
    anno <- annotatePSM(
        spectra$MZ, spectra$Prob,
        spectra$Charge,
        pep, 1:2, "C13",
        prob, isoCenter, 4.0, TRUE
    )
    return(c(XcorrScore = anno$XcorrScore, MVHscore = anno$MVHscore))
})
scores <- t(scores)
scores <- cbind(WDPscores, scores)

## -----------------------------------------------------------------------------
scores_df <- as.data.frame(scores)
scores_df$abundance <- 0:100

# Reshape the data to long format for faceting
scores_long <- scores_df %>%
    pivot_longer(
        cols = c(WDPscores, XcorrScore, MVHscore),
        names_to = "score_type",
        values_to = "score_value"
    )

# Create a publication-quality plot
ggplot(scores_long, aes(x = abundance, y = score_value)) +
    geom_line() +
    facet_wrap(~score_type, scales = "free_y") +
    labs(
        x = expression(paste(~ {}^{
            13
        }, "C %")),
        y = "Score",
        title = expression(paste("Scores Across", ~ {}^{
            13
        }, "C %"))
    ) +
    theme_minimal() +
    theme(
        plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
        axis.title = element_text(size = 14),
        axis.text = element_text(size = 12),
        strip.text = element_text(size = 12, face = "bold"),
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "white", colour = "black")
    )

## ----session-info-------------------------------------------------------------
sessionInfo()