inst/doc/tbspVignette.R

In TBSignatureProfiler: Profile RA-Seq Data Using TB Pathway Signatures

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

## ----setup--------------------------------------------------------------------
suppressPackageStartupMessages({
  library(TBSignatureProfiler)
  library(SummarizedExperiment)
})

## ----Run shiny app, eval = FALSE----------------------------------------------
#  TBSPapp()

## ----loading data-------------------------------------------------------------
## HIV/TB gene expression data, included in the package
hivtb_data <- TB_hiv

### Note that we have 25,369 genes, 33 samples, and 1 assay of counts
dim(hivtb_data)

# We start with only one assay
assays(hivtb_data)

## ----add assays---------------------------------------------------------------
## Make a log counts, CPM and log CPM assay
hivtb_data <- mkAssay(hivtb_data, log = TRUE, counts_to_CPM = TRUE)

### Check to see that we now have 4 assays
assays(hivtb_data)

## ----run data-----------------------------------------------------------------
## List all signatures in the profiler
data("TBsignatures")
names(TBsignatures)

## We can use all of these signatures for further analysis
siglist_hivtb <- names(TBsignatures)

## ----run profiler-------------------------------------------------------------
## Run the TBSignatureProfiler to score the signatures in the data
out <- capture.output(ssgsea_result <- runTBsigProfiler(input = hivtb_data,
                                                 useAssay = "log_cpm",
                                                 signatures = TBsignatures,
                                                 algorithm = "ssGSEA",
                                                 combineSigAndAlgorithm = TRUE,
                                                 parallel.sz = 1))

## ----show scores, message = FALSE---------------------------------------------
## New colData entries from the Profiler
sigs <- c("Anderson_42", "Anderson_OD_51", "Berry_393")
ssgsea_print_results <- as.data.frame(colData(ssgsea_result))[, c("Disease", sigs)]
ssgsea_print_results[,2:4] <- round(ssgsea_print_results[, 2:4], 4)

DT::datatable(ssgsea_print_results)

## ----all sigs heatmap, message = FALSE, fig.height = 9------------------------

# Colors for gradient
colors <- RColorBrewer::brewer.pal(6, "Spectral")
col.me <- circlize::colorRamp2(seq(from = -2, to = 2, 
                                   length.out = 6), colors)

signatureHeatmap(ssgsea_result, name = "Heatmap of Signatures, 
                 ssGSEA Algorithm", 
                 signatureColNames = names(TBsignatures),
                 annotationColNames = "Disease",
                 scale = TRUE,
                 showColumnNames = TRUE,
                 choose_color = col.me)

## ----Boxplots all, message = FALSE, results = 'hide', fig.height = 9, fig.width=15----

signatureBoxplot(inputData = ssgsea_result, 
                 name = "Boxplots of Signatures, ssGSEA", 
                 signatureColNames = names(TBsignatures),
                 annotationColName = "Disease", rotateLabels = FALSE)

## ----compareAlgs ex, warning = FALSE, message = FALSE-------------------------
# Heatmap
compareAlgs(hivtb_data, annotationColName = "Disease",
            scale = TRUE,
            algorithm = c("GSVA", "ssGSEA",
                          "singscore", "PLAGE", "Zscore"),
            useAssay = "log_counts",
            signatures = TBsignatures[1],
            choose_color = col.me, show.pb = FALSE,
            parallel.sz = 1)

# Boxplot
compareAlgs(hivtb_data, annotationColName = "Disease",
            scale = TRUE,
            algorithm = c("GSVA", "ssGSEA",
                          "singscore", "PLAGE", "Zscore"),
            useAssay = "log_counts",
            signatures = TBsignatures[1],
            choose_color = col.me, show.pb = FALSE,
            parallel.sz = 1, output = "boxplot")

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