In milaboratory/mirutil:

knitr::opts_chunk$set(echo = TRUE)

Install package via devtools:

library(ggplot2)
devtools::install_github("milaboratory/mirutil")
library(mirutil); packageVersion("mirutil")
#devtools::unload()
devtools::document()
#devtools::load_all()

Load metadata and samples

metadata <- fread("metadata.txt")
dataset <- metadata %>% read_mixcr_dataset

Compute rearrangement statistics

#library(kableExtra) # for tables

df.stats <- dataset %>%
  compute_rearr_stats

df.stats %>%
  str

Plot segment usage for alpha chain

# packages for plotting
library(ggplot2)
library(forcats) 

plt_segm <- function(m) {
  m %>%
    ggplot(aes(x = fct_reorder(segment, freq.reads), 
             group = sample.id,
             y = freq.reads, color = sample.id)) +
    geom_line() +
    coord_flip() +
    xlab("") + ylab("Segment frequency") +
    scale_color_brewer(palette = "Set1") +
    facet_wrap(~segment.type, scales = "free") +
    theme_bw() +
    theme(legend.position = "bottom",
          axis.text.y = element_text(family = "mono", size = 6),
          strip.background = element_rect(fill = NA, color = NA),
          panel.grid.major.x = element_blank(),
          panel.grid.minor = element_blank())
}

df.stats %>%
  .$segment.usage %>%
  filter(chain == "TRA") %>%
  plt_segm

For beta chain

df.stats %>%
  .$segment.usage %>%
  filter(chain == "TRB") %>%
  plt_segm

Plot insert size distribution

plt_insert_size <- function(m) {
  m %>%
    ggplot(aes(x = insertions,
               color = sample.id)) +
    geom_line(aes(group = sample.id, y = freq.reads)) +
    scale_x_continuous("Insertions", limits = c(0,30)) + ylab("Frequency") +
    scale_color_brewer(palette = "Set1") +
    facet_wrap(~chain+ins.type, scales = "free") +
    theme_bw() +
    theme(legend.position = "bottom",
          strip.background = element_rect(fill = NA, color = NA),
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank())
}

df.stats %>%
  .$insert.size %>%
  plt_insert_size

Plot insert profiles

plt_insert_profile <- function(m) {
  m %>%
    ggplot(aes(x = paste(nt.1, nt.2, sep = "\n"),
               color = sample.id)) +
    geom_line(aes(group = sample.id, y = freq.reads)) +
    xlab("Base #1 / Base #2") + ylab("Frequency") +
    scale_color_brewer(palette = "Set1") +
    facet_wrap(~chain+ins.profile.type, scales = "free") +
    theme_bw() +
    theme(legend.position = "bottom",
          strip.background = element_rect(fill = NA, color = NA),
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank())
}

df.stats %>%
  .$insert.profile %>%
  plt_insert_profile

Plot V/J trimming profiles

plt_deletion_size <- function(m) {
  m %>%
    ggplot(aes(x = deletions,
               group = paste(segment, sample.id),
               color = sample.id,
               y = freq.reads)) +
    geom_line() +
    xlab("Deletion size") + ylab("Frequency") +
    scale_color_brewer(palette = "Set1") +
    facet_wrap(~chain+del.type, scales = "free") +
    theme_bw() +
    theme(legend.position = "bottom",
          strip.background = element_rect(fill = NA, color = NA),
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank())
}

df.stats %>%
  .$deletion.size %>%
  plt_deletion_size

Get distances between samples based on rearrangement statistics

df.dist <- df.stats %>%
  compute_rearr_stat_dist

df.dist %>%
  head

Convert distances to MDS coordinates

df.mds <- df.dist %>% compute_rearr_stat_mds

Plot MDS for segment usage

plt_mds <- function(m) {
  m %>%
    ggplot(aes(x = mds.x, y = mds.y, color = sample.id)) +
    geom_point() +
    xlab("") + ylab("") +
    facet_grid(value.type~chain + type, scales = "free") +
    scale_color_brewer(palette = "Set1") +
    theme_bw() +
    theme(legend.position = "bottom",
          strip.background = element_rect(fill = NA, color = NA),
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank())
}

df.mds %>%
  filter(statistic == "segment.usage") %>%
  plt_mds

Plot MDS for added N-bases (inserts)

df.mds %>%
  filter(statistic == "insert.size") %>%
  plt_mds

Plot MDS for trimmed bases from segment ends

df.mds %>%
  filter(statistic == "deletion.size") %>%
  plt_mds

Plot MDS for insert profile

df.mds %>%
  filter(statistic == "insert.profile") %>%
  plt_mds

milaboratory/mirutil documentation built on May 20, 2019, 1:24 p.m.