## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## ----library_prep, include=FALSE, message=FALSE, warning=FALSE----------------
library(DamageDetective)
library(Seurat)
library(ggplot2)
library(patchwork)

## ----seurat_data, eval=FALSE--------------------------------------------------
# # Retrieve the dataset of interest
# SeuratData::InstallData("pbmc3k")
# data("pbmc3k")
# 
# # Extract the count matrix
# pbmc3k_counts <- pbmc3k[["RNA"]]$counts

## ----retrieve_alignment, echo=FALSE-------------------------------------------
url <- "https://cf.10xgenomics.com/samples/cell-exp/3.0.0/pbmc_1k_v3/pbmc_1k_v3_filtered_feature_bc_matrix.tar.gz"

# Set a temporary directory
temp_dir <- tempdir()
tar_file <- file.path(tempdir(), "pbmc_1k_v3_filtered_feature_bc_matrix.tar.gz")

# Download files
download.file(url, tar_file, mode = "wb", quiet = TRUE)
untar(tar_file, exdir = temp_dir)
extracted_dir <- file.path(temp_dir, "filtered_feature_bc_matrix")
alignment_counts <- Seurat::Read10X(extracted_dir)

## ----demonstrate_alignment_conversion, eval=FALSE-----------------------------
# # Set the file paths relative to location on your device
# matrix_file <- "~/Projects/demonstrations/matrix.mtx.gz"
# barcodes_file <- "~/Projects/demonstrations/barcodes.tsv.gz"
# features_file <- "~/Projects/demonstrations/features.tsv.gz"
# 
# # Construct the sparse matrix
# alignment_counts <- Seurat::ReadMtx(
#   mtx = matrix_file,
#   cells = barcodes_file,
#   features = features_file
# )

## ----sce_data, eval=FALSE-----------------------------------------------------
# # Retrieve multisample dataset
# pbmc_sce <- scRNAseq::fetchDataset("kotliarov-pbmc-2020", "2024-04-18")
# 
# # Extract sample of interest
# metadata <- SummarizedExperiment::colData(pbmc_sce)
# sample_sce <- subset(metadata, sample == "234_d0")
# sample_sce <- rownames(sample_sce)
# 
# # Subset and convert to sparse format
# pbmc_counts <- SummarizedExperiment::assay(pbmc_sce, "counts")
# sample_counts  <- pbmc_counts[, sample_sce]
# sample_counts <- as.matrix(sample_counts)
# sample_counts <- as(sample_counts, "dgCMatrix")

## ----view_penalty-------------------------------------------------------------
penalty_plot <- DamageDetective::simulate_counts(
  count_matrix = alignment_counts,
  ribosome_penalty = 0.01,  
  damage_proportion = 0.05,  
  target_damage = c(0.5, 1), 
  plot_ribosomal_penalty = TRUE,
  seed = 7
)

## ----select_penalty-----------------------------------------------------------
selected_penalty <- DamageDetective::select_penalty(
  count_matrix = alignment_counts,
  max_penalty_trials = 3, # shortened for vignette, default is 10
  seed = 7, 
  verbose = TRUE
)
selected_penalty

## ----eval = FALSE-------------------------------------------------------------
# damage_levels <- list(
#   pANN_50 = c(0.1, 0.5),
#   pANN_100 = c(0.5, 1)
# )

## ----run_detection------------------------------------------------------------
# Run detection
detection_output <- DamageDetective::detect_damage(
  count_matrix = alignment_counts, 
  ribosome_penalty = selected_penalty,
  seed = 7
)

# View output
head(detection_output$output[, -1])
table(detection_output$output$DamageDetective_filter)


## ----session-info, echo=FALSE-------------------------------------------------
sessionInfo()