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

## ----install, eval = FALSE----------------------------------------------------
# if (!requireNamespace("BiocManager", quietly = TRUE))
#     install.packages("BiocManager")
# BiocManager::install("STRprofilerR")
# 
# # You can install the development version from Github
# pak::pak("j-andrews7/STRprofilerR")

## ----read-long----------------------------------------------------------------
library(STRprofilerR)
longFile <- system.file("extdata", "ExampleSTR_long.csv", package = "STRprofilerR")
head(read.csv(longFile, check.names = FALSE))

## ----read-long-2--------------------------------------------------------------
profiles <- readSTRProfiles(longFile, sampleCol = "Sample Name")
profiles

## ----read-both----------------------------------------------------------------
xlsxFile <- system.file("extdata", "ExampleSTR.xlsx", package = "STRprofilerR")

if (requireNamespace("readxl", quietly = TRUE)) {
    both <- readSTRProfiles(c(longFile, xlsxFile), sampleCol = "Sample Name")
    print(both)
}

## ----alleles------------------------------------------------------------------
alleles(profiles)[["marker1"]]

## ----clean--------------------------------------------------------------------
cleanAlleles("10.0,10,13,13.0,14,14 ")

## ----markers------------------------------------------------------------------
harmonizeMarkers(c("Penta D", "Penta_E", "vWA"))

## ----marker-data--------------------------------------------------------------
markerData(profiles)

## ----metadata-----------------------------------------------------------------
dbFile <- system.file("extdata", "main_database.csv", package = "STRprofilerR")
reference <- readSTRProfiles(dbFile)

dim(reference)
head(sampleData(reference), 3)

## ----score--------------------------------------------------------------------
queries <- readSTRProfiles(longFile, sampleCol = "Sample Name")
refs <- readSTRProfiles(
    system.file("extdata", "ExampleSTR_database.csv", package = "STRprofilerR"),
    sampleCol = "Sample Name"
)

scoreProfiles(queries, refs)

## ----score-filter-------------------------------------------------------------
scoreProfiles(queries, refs, minScore = 90)

## ----compare------------------------------------------------------------------
cmp <- compareProfiles(queries, refs)
cmp

## ----summary------------------------------------------------------------------
as.data.frame(summary(cmp))[, c("Sample", "mixed", "topHit", "nextBest")]

## ----compare-self-------------------------------------------------------------
compareProfiles(refs)

## ----mixing-------------------------------------------------------------------
batch <- readSTRProfiles(
    system.file("extdata", "Example_Batch_File.csv", package = "STRprofilerR")
)

flagMixedSamples(batch)

## ----write--------------------------------------------------------------------
outDir <- file.path(tempdir(), "STRprofilerR-vignette")
written <- writeSTRResults(cmp, outDir)

basename(written)

## ----validate-----------------------------------------------------------------
validateClastrMarkers(markers(refs))

## ----clastr, eval = FALSE-----------------------------------------------------
# profile <- c(
#     Amelogenin = "X", CSF1PO = "11,12", D2S1338 = "19,23", D3S1358 = "15,17",
#     D5S818 = "11,12", D7S820 = "10", D8S1179 = "10", D13S317 = "11,12",
#     D16S539 = "11,12", D18S51 = "13", D19S433 = "14", D21S11 = "29,30",
#     FGA = "20,22", PentaD = "11,13", PentaE = "14,16", TH01 = "6,9",
#     TPOX = "8,9", vWA = "17,19"
# )
# 
# hits <- clastrQuery(profile, scoreFilter = 90)
# as.data.frame(hits)[, c("accession", "name", "score", "problem")]

## ----cli-install, eval = FALSE------------------------------------------------
# install.packages("Rapp")
# Rapp::install_pkg_cli_apps("STRprofilerR")

## ----app, eval = FALSE--------------------------------------------------------
# STRprofilerApp()
# 
# # Against your own database instead of the bundled one.
# STRprofilerApp(database = "our_database.csv")

## ----bioc-classes-------------------------------------------------------------
lengths(alleles(batch)[["D5S818"]])

## ----pdx----------------------------------------------------------------------
data(pdxSTRProfiles)
pdxSTRProfiles

pdxCmp <- compareProfiles(batch, pdxSTRProfiles)
as.data.frame(summary(pdxCmp))[, c("Sample", "mixed", "topHit", "nextBest")]

## ----from-df------------------------------------------------------------------
calls <- data.frame(
    Sample = "PDX_042",
    Amelogenin = "X", CSF1PO = "11", D13S317 = "9,14", D16S539 = "12",
    D18S51 = "17", D21S11 = "29,31", D3S1358 = "15", D5S818 = "11",
    D7S820 = "9,11", D8S1179 = "13,15", FGA = "21,22", TH01 = "9.3",
    TPOX = "8", vWA = "14"
)
newModel <- STRProfiles(calls)

scoreProfiles(newModel, pdxSTRProfiles, minScore = 80)

## ----se, eval = requireNamespace("SummarizedExperiment", quietly = TRUE), message = FALSE----
library(SummarizedExperiment)

# A stand-in for an experiment on the three models in the batch.
counts <- matrix(
    rpois(30, lambda = 50),
    nrow = 10,
    dimnames = list(paste0("gene", 1:10), rownames(batch))
)

se <- SummarizedExperiment(
    assays = list(counts = counts),
    colData = DataFrame(
        model = c("J000077451", "J000077960", "J000077451"),
        row.names = rownames(batch)
    )
)

# Each sample's best match in the reference.
pairs <- scores(pdxCmp)
pairs <- pairs[order(pairs$query, -pairs$tanabeScore), ]
best <- pairs[!duplicated(pairs$query), ]
best <- best[match(colnames(se), best$query), ]

se$strMatch <- best$reference
se$strTanabe <- best$tanabeScore
se$strMixed <- summary(pdxCmp)[colnames(se), "mixed"]
se$authenticated <- se$strMatch == se$model & se$strTanabe >= 80 & !se$strMixed
colData(se)

# Keep only authenticated samples for analysis.
se[, se$authenticated]

## ----grow-db------------------------------------------------------------------
inHouse <- c(pdxSTRProfiles, newModel)
inHouse

writeSTRProfiles(inHouse, file.path(tempdir(), "in_house_database.csv"))

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