How to make chromosome color maps for bed ranges

Question

I have genomic .bed file data of 4 different types; type A,B,C,D. These are some examples-

Type A:

1   101380000   101710000   A
1   110085000   110320000   A

Type B:

1   100930000   101335000   B
1   10430000    10560000    B

Type C:

1   101775000   101825000   C
1   103570000   104070000   C

Type D:

1   153325000   153475000   D
1   156355000   156415000   D

This data is for chr1, but I such data for all the chromosomes. I want to draw a chromosome map (possibly using the Gviz package?) which can assign different colors to the 4 types of genomic ranges and make me an idiogram like plot for the chromosomes I wish to see, so I know where my ranges physically are, on the chromosome.

Can this be done using the Gviz package in R? If not, what is the best way and code to do this?

Answer 1

You can use karyoploteR to plot your regions on an ideogram quite easily.

Disclosure: I'm one of the authors of the tool.

For this example I'll start creating 4 GRanges with random regions, but with your data you should load your bed into R and split it into four GRanges.

library(karyoploteR)


typeA <- createRandomRegions(nregions = 100, length.mean = 5e6, length.sd = 3e6)
typeB <- createRandomRegions(nregions = 20, length.mean = 15e6, length.sd = 10e6)
typeC <- createRandomRegions(nregions = 1000, length.mean = 5e5, length.sd = 2e5)
typeD <- createRandomRegions(nregions = 100, length.mean = 2e6, length.sd = 5e6)

colors <- c(A="#FFD700", B="#00BFFF", C="#FF7F50", D="#90EE90")

And we can plot them. We have to start creating an empty ideogram with plotKaryotype (we are working with human here, but you can use any genome you want) and then plot the regions in the GRanges with kpPlotRegions. As a bonus, if two regions overlap it will plot them one above the other.

kp <- plotKaryotype(genome = "hg19")

kpPlotRegions(kp, data=typeA, r0=0, r1=0.2, col=colors["A"])
kpPlotRegions(kp, data=typeB, r0=0.25, r1=0.45, col=colors["B"])
kpPlotRegions(kp, data=typeC, r0=0.5, r1=0.7, col=colors["C"])
kpPlotRegions(kp, data=typeD, r0=0.75, r1=0.95, col=colors["D"])
legend(x = "right", legend = names(colors), fill = colors, border=NA)

And you'll get something similar to this

In this example we have fairly large regions and not a lot of them. As @heathobrien pointed out, plotting the actual regions in the whole genome plot might not be the best visualization option. You can either use zoom to create different plots of smaller parts of the genome or change your visualization to show the density of regions along the genome with kpPlotDensity.

Plotting the density and tweaking a bit the visualization (and changing the example to have 10 times more regions) with this codei

typeA <- createRandomRegions(nregions = 1000, length.mean = 5e6, length.sd = 3e6)
typeB <- createRandomRegions(nregions = 200, length.mean = 15e6, length.sd = 10e6)
typeC <- createRandomRegions(nregions = 10000, length.mean = 5e5, length.sd = 2e5)
typeD <- createRandomRegions(nregions = 1000, length.mean = 2e6, length.sd = 5e6)

kp <- plotKaryotype(genome = "hg19", plot.type = 4, 
                    ideogram.plotter = NULL, labels.plotter = NULL, main="Density")
kpAddCytobandsAsLine(kp)
kpAddChromosomeNames(kp, srt=45)

kpPlotDensity(kp, data=typeA, r0=0, r1=0.2, col=colors["A"], window.size     = 10e6)
kpAddLabels(kp, labels = "Type A", r0=0, r1=0.2)
kpPlotDensity(kp, data=typeB, r0=0.25, r1=0.45, col=colors["B"], window.size = 10e6)
kpAddLabels(kp, labels = "Type B", r0=0.25, r1=0.45)
kpPlotDensity(kp, data=typeC, r0=0.5, r1=0.7, col=colors["C"], window.size = 10e6)
kpAddLabels(kp, labels = "Type C", r0=0.5, r1=0.7)
kpPlotDensity(kp, data=typeD, r0=0.75, r1=0.95, col=colors["D"], window.size = 10e6)
kpAddLabels(kp, labels = "Type D", r0=0.75, r1=0.95)

it would look something like this

The plots could be tweaked a lot changing various parameters and you could add other data on top of it if needed. You can find more information on the Bioconductor landing page or at the karyoploteR's Tutorial and Examples page.