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

## ---- message=FALSE, warning=FALSE--------------------------------------------
library(MSstatsPTM)

## ---- eval = FALSE------------------------------------------------------------
#  if (!requireNamespace("BiocManager", quietly = TRUE))
#      install.packages("BiocManager")
#  
#  BiocManager::install("MSstatsPTM")

## ----raw_data_ptm, eval = FALSE-----------------------------------------------
#  
#  annotation <- data.frame('Condition' = c('Control', 'Control', 'Control',
#                                           'Treatment', 'Treatment', 'Treatment'),
#                           'BioReplicate' = c(1,2,3,4,5,6),
#                           'Run' = c('prot_run_1', 'prot_run_2', 'prot_run_3',
#                                     'phos_run_1', 'phos_run_2', 'phos_run_3'),
#                           'Type' = c("Protein", "Protein", "Protein", "PTM",
#                                      "PTM", "PTM"))
#  
#  # Run MSstatsPTM converter with modified and unmodified datasets.
#  raw.input <- ProgenesistoMSstatsPTMFormat(raw_ptm_df, annotation,
#                                            raw_protein_df, fasta_path)

## ----raw_data, eval = FALSE---------------------------------------------------
#  # Add site into ProteinName column
#  raw_ptm_df$ProteinName <- paste(raw_ptm_df$ProteinName,
#                                  raw_ptm_df$Site, sep = "_")
#  
#  # Run MSstats Converters
#  PTM.data <- ProgenesistoMSstatsFormat(raw_ptm_df, annotation)
#  PROTEIN.data <- ProgenesistoMSstatsFormat(raw_protein_df, annotation)
#  
#  # Combine into one list
#  raw.input <- list(PTM = PTM.data,
#                    PROTEIN = PROTEIN.data)

## -----------------------------------------------------------------------------
head(raw.input$PTM)
head(raw.input$PROTEIN)

## ----summarize, message=FALSE, warning=FALSE----------------------------------

MSstatsPTM.summary <- dataSummarizationPTM(raw.input, verbose = FALSE)

head(MSstatsPTM.summary$PTM$ProteinLevelData)
head(MSstatsPTM.summary$PROTEIN$ProteinLevelData)

## ----qcplot, message=FALSE, warning=FALSE-------------------------------------
dataProcessPlotsPTM(MSstatsPTM.summary,
                    type = 'QCPLOT',
                    which.PTM = "allonly",
                    address = FALSE)

## ----profileplot, message=FALSE, warning=FALSE--------------------------------

dataProcessPlotsPTM(MSstatsPTM.summary,
                    type = 'ProfilePlot',
                    which.Protein = "Q9Y6C9",
                    address = FALSE)

## ----model, message=FALSE, warning=FALSE--------------------------------------

# Specify contrast matrix
comparison <- matrix(c(-1,0,1,0),nrow=1)
row.names(comparison) <- "CCCP-Ctrl"
colnames(comparison) <- c("CCCP", "Combo", "Ctrl", "USP30_OE")

MSstatsPTM.model <- groupComparisonPTM(MSstatsPTM.summary, 
                                       data.type = "LabelFree",
                                       contrast.matrix = comparison,
                                       verbose = FALSE)
head(MSstatsPTM.model$PTM.Model)
head(MSstatsPTM.model$PROTEIN.Model)
head(MSstatsPTM.model$ADJUSTED.Model)

## ----volcano, message=FALSE, warning=FALSE------------------------------------
groupComparisonPlotsPTM(data = MSstatsPTM.model,
                        type = "VolcanoPlot",
                        FCcutoff= 2,
                        logBase.pvalue = 2,
                        address=FALSE)

## ----heatmap, message=FALSE, warning=FALSE------------------------------------
groupComparisonPlotsPTM(data = MSstatsPTM.model,
                        type = "Heatmap",
                        which.PTM = 1:30,
                        address=FALSE)