As requested by some folks at the AIRR meeting last week I am happy to provide the R code for the VDJ plots I have presented in my talk. A fancier way of making these plots (such as Fig. S2 in our J Immunol paper: http://www.jimmunol.org/content/jimmunol/194/1/252.full.pdf?with-ds=yes) needs some longer code which I can also provide if people are interested.
Load libraries
x <- c(“tidyverse”, “alakazam”, “ggthemes”)
sapply(x, require, character.only = T)
rm(x)
Download example files from below URL:
Load .tab data into R
out <- readChangeoDb("~/Downloads/Changeo_Example/S43_db-pass_parse-select.tab")
Modify V/D/J_CALL
out$V_GENE <- as.character(substr(lapply(strsplit(as.character(out$V_CALL), “\"), “[”, 1), 8, 20))
out$D_GENE <- as.character(substr(lapply(strsplit(as.character(out$D_CALL), "\”), “[”, 1), 8, 20))
out$J_GENE <- as.character(substr(lapply(strsplit(as.character(out$J_CALL), “\*”), “[”, 1), 8, 20))
Add Isotype information from CPRIMER and make a factor
out$ISOTYPE <- out$CPRIMER %>% gsub("-PCR", “”, .) %>% factor(levels = c(“IgD”, “IgM”, “IgA”, “IgG”, “IgE”))
VDJ plot
a <- out %>% group_by(ISOTYPE, MID, V_GENE, D_GENE, J_GENE) %>% summarise(Frequency = sum(DUPCOUNT))
a[0, ] %>%
ggplot(aes(V_GENE, D_GENE, group = J_GENE, color = J_GENE, size = Frequency)) +
theme_few(8) + theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5), strip.text = element_text(size = 12, face = “bold”)) +
geom_point(data = a %>% filter(J_GENE == “IGHJ1”), position = position_nudge(x = -.2, y = .1), alpha = .5, shape = 16) +
geom_point(data = a %>% filter(J_GENE == “IGHJ2”), position = position_nudge(x = 0, y = .1), alpha = .5, shape = 16) +
geom_point(data = a %>% filter(J_GENE == “IGHJ3”), position = position_nudge(x = .2, y = .1), alpha = .5, shape = 16) +
geom_point(data = a %>% filter(J_GENE == “IGHJ4”), position = position_nudge(x = -.2, y = -.1), alpha = .5, shape = 16) +
geom_point(data = a %>% filter(J_GENE == “IGHJ5”), position = position_nudge(x = 0, y = -.1), alpha = .5, shape = 16) +
geom_point(data = a %>% filter(J_GENE == “IGHJ6”), position = position_nudge(x = .2, y = -.1), alpha = .5, shape = 16) +
facet_grid(MID ~ ISOTYPE)