knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

library(utilitarian)

Read in SAM files as tibble

toy <- read_sam(system.file("extdata", "toy.sam", package = "utilitarian"))

toy

Create Tidy CIGAR Data

Results in a long table of tidy data showing start and end positions of the operations along with their English translations.

Due to the facetting of CIGAR visualizations on Reference names (RNAME) it is wise to group them together if the organism is what matters. For example in the case of the human genome vs sars coronavirus grouping all the human chromosomes together is useful.

sam %>% 
  mutate(RNAME = case_when(RNAME == "MN908947.3" ~ "hcov",
                           RNAME == "*" ~ "unmapped",
                           TRUE ~ "human")) %>% 
  tidy_cigar()
toy.tc <- tidy_cigar(toy)

toy.tc

Plot a CIGAR

Plotting a cigar from a single read in the tidied cigar table is done using plot_cigar() and specifying the read. Plots are facetted based on the references and the mapped position of the read. Returns a plot.

plot_cigar(toy.tc, "r002")

Plot all CIGARs

A wrapper around calling plot_cigar() on all the reads in the SAM file. Will warn and stop if more than 20 elements to be created, override this with the lg=TRUE argument.

p <- plot_all_cigars(toy.tc)

toy[toy$QNAME == "r001", "CIGAR"]
p[["r001"]]
toy[toy$QNAME == "x3", "CIGAR"]
p[["x3"]]

Postscript

SAM headers contain useful information but for the purposes of cigar visualization they are ignored. Having access to them can be helpful so there is an included function to pull headers (as I couldn't figure out a way to do so with readr).

Show Reference Lengths

h <- read_sam_headers(system.file("extdata", "toy.sam", package = "utilitarian"))

h$References

Show Command line params and software details

h$Mapping


TheZetner/utilitarian documentation built on Aug. 13, 2022, 12:31 p.m.