rock
is an R package that (hopefully) helps with the day to day
bioinformatics life at UMCCR (UniMelb Centre for Cancer Research).
You can do the following:
Create Perl circos plots using structural variant calls from Manta, and copy number variant calls from CNVkit or PURPLE.
Create CNV profiles in horizontal facets for multiple samples or callers (piano plots). Can also zoom into specific chromosomes, and include an ideogram when specifying a single chromosome.
Generate bedgraph files for viewing the copy number segments in IGV as a bar plot.
Generate IGV files for viewing SNP values in IGV as a scatter plot.
You can install the development version of rock
from
GitHub with:
# install.packages("devtools") # if not pre-installed
devtools::install_github("pdiakumis/rock") # master version
devtools::install_github("pdiakumis/rock@v1.2.3") # release v1.2.3
devtools::install_github("pdiakumis/rock@abcd") # commit abcd
require(rock)
There is a conda package version at https://anaconda.org/pdiakumis/r-rock which is updated regularly.
You need to create a conda environment, and then install with:
conda install -c pdiakumis r-rock
We can generate circos plots using the original circos software package, written in Perl.
circos
needs to be installed in your PATH
for it to work straight
from within rock
. I’ve had a terrible time trying to install circos
either from source or from conda. The only way that I’ve found to
consistently work is through Docker.
I have a docker image available on
DockerHub, which you can
pull (docker pull pdiakumis/circos:0.69-6
) and then run as shown
below.
Start by preparing the Manta and CNVkit calls. The required input files
will be written to
outdir
:
manta <- system.file("extdata", "HCC2218_manta.vcf", package = "pebbles")
cnvkit <- system.file("extdata", "HCC2218_cnvkit-call.cns", package = "pebbles")
outdir <- "man/figures/perl_circos"
rock::circos_prep(outdir = outdir, manta = manta, cnv = cnvkit)
#> Warning in dir.create(outdir, recursive = TRUE): 'man/figures/perl_circos'
#> already exists
#> Exporting Manta and/or CNV circos files to 'man/figures/perl_circos'.
#> Registered S3 method overwritten by 'R.oo':
#> method from
#> throw.default R.methodsS3
#> Copying circos templates to 'man/figures/perl_circos'. 'template' is cnvsv.
Now comes the fun part of running the circos
command. Depending on if
you’ve managed to install a working version of circos
or if you want
to use the Docker version, you have the following options:
circos
in your R session, just run the following R
command. If you’re an RStudio user, you can make sure it recognises
the user’s PATH by opening the RStudio app via the terminal, or
perhaps following the suggestions here:
https://stackoverflow.com/questions/31121645plot_circos(outdir = outdir, name = "my_fabulous_plot")
circos -nosvg -conf <outdir>/circos_simple.conf -outputdir <outdir> -outputfile my_fabulous_plot_circos.png
outdir
:docker container run --rm -v $(pwd):/data pdiakumis/circos:0.69-6 -conf /data/circos.conf -outputdir /data
knitr::include_graphics("man/figures/perl_circos/foo_circos_cnvkit_manta.png")
We can generate ‘piano’ plots to compare CNV calls from multiple callers or samples.
Start by preparing the SV and CNV calls.
manta <- system.file("extdata", "HCC2218_manta.vcf", package = "pebbles")
cnvkit <- system.file("extdata", "HCC2218_cnvkit-call.cns", package = "pebbles")
facets <- system.file("extdata", "HCC2218_facets_cncf.tsv", package = "pebbles")
titan <- system.file("extdata", "HCC2218_titan.segs.tsv", package = "pebbles")
purple <- system.file("extdata", "HCC2218_purple.cnv.tsv", package = "pebbles")
truth <- system.file("extdata", "HCC2218_truthset_cnv_bcbio.tsv", package = "pebbles")
sv_manta <- prep_manta_vcf(manta)
cn_facets <- prep_facets_seg(facets)
cn_cnvkit <- prep_cnvkit_seg(cnvkit)
cn_purple <- prep_purple_seg(purple)
cn_truth <- prep_truth_seg(truth)
cn_titan <- prep_titan_seg(titan) # titan needs -1 for this case
cn_titan$cnv$tot_cn <- cn_titan$cnv$tot_cn - 1
cnv_list <- list(truth = cn_truth, cnvkit = cn_cnvkit, facets = cn_facets, purple = cn_purple, titan = cn_titan)
plot_piano(cnv_list = cnv_list)
plot_piano(cnv_list = cnv_list, chromosomes = c("1", "7", "8"), hide_x_lab = FALSE)
plot_piano(cnv_list = cnv_list, chromosomes = c("1", "7", "8"),
seg.col = c("orange", "lightblue", "blue", "pink"), hide_x_lab = FALSE)
require(patchwork)
#> Loading required package: patchwork
plot_ideogram(chrom = "13") +
plot_piano(cnv_list = cnv_list,
chromosomes = "13", hide_x_lab = FALSE) +
plot_layout(ncol = 1, heights = c(1, 15))
cn_fname <- system.file("extdata", "HCC2218_purple.cnv.tsv", package = "pebbles")
cnv <- read_cnv(cn_fname)
cnv2igv(cnv, out_file = "~/Desktop/tmp/cnv_segs4igv.bedgraph", track_name = "cnv_segs2")
knitr::include_graphics("man/figures/README-cnv2igv_output.png")
bed <- system.file("extdata", "HCC2218_baf.tsv", package = "pebbles")
bedval2igv(bed, out_file = "~/Desktop/tmp/baf1.igv", track_name = "baf", col = "purple")
# example for COLO829 whole-genome BAFs
knitr::include_graphics("man/figures/README-bedval2igv_output.png")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.