Skip to main content

Hi-C Visualisation Example

You may find that you just want to visualise some already processed/published Hi-C data. This should be uploaded somewhere like GEO or ENCODE and you can download the processed files and visualise them. You can do this on the UCSC genome browser and in R using plotgardender.

Depending on how they’ve analysed it you may have different data formats. If it was analysed using HiC-Pro you will have .matrix files which are tab delimited text files with the format:

bin1    bin2    count
2   23  146.817177
2   75  60.727011
2   810 11.043244
2   945 9.987969
2   2908    10.181468
2   5536    12.709525

with a complementary .bed file which contains the regions that each of the bins corresponds to:

chr start   end bin
chr1    0   10000   1
chr1    10000   20000   2
chr1    20000   30000   3
chr1    30000   40000   4
chr1    40000   50000   5
chr1    50000   60000   6
chr1    60000   70000   7
chr1    70000   80000   8

The matrix is how many contacts are counted between each of the pairs of bins. You can visualise this as a heatmap, with the intensity of the colour representing the number of contacts. Plotgardener has a function to do this called plot_hic_matrix.

The files are very large and so you may want to subset to a chromosome or region of interest first. I’ve done this with awk to subset to chromosome 4.

# Get bin IDs for target chromosome (bin ID is column 4)
awk '$1 == "chr4" {print $4}' GSE102197_bin_10k.bed > chr4_bin_ids.txt

# Load all bin IDs from chr4_bin_ids.txt into an array called "bins"
# Then go through the matrix file and print only the rows where BOTH interacting bins are chr4 bins
awk 'NR==FNR { bins[$1]=1; next } $1 in bins && $2 in bins' chr4_bin_ids.txt \
    GSE102197_adult_10k_iced.matrix > chr4_matrix.txt

This will give you a smaller matrix file with only the interactions between bins on chromosome 4. You can then visualise this with plotgardener or in the UCSC genome browser.

library(plotgardener)
library("TxDb.Hsapiens.UCSC.hg19.knownGene")
library("AnnotationHub")

##### HIC plotting with plotgardener #####

##### Load in the matrix and bin files #####
Iced.matrix= read.delim("chr4_matrix.txt", sep = "\t", 
                        header=TRUE, col.names = c("Bin1", "Bin2", "Interaction_freq"), skip=0)

bin_index= read.delim("GSE102197_bin_10k.bed", sep = "\t", 
                      header=TRUE, col.names = c("Chr", "Start", "End", "Bin"), skip=0)

### Subset to the region of interest which is chr4:30000000-38978218
bin_index = bin_index[bin_index$Chr == 'chr4',]
bin_index = bin_index[bin_index$Start >= 30000000 & bin_index$End <= 38978218,]

### Add the genomic coordinates to the matrix file by merging with the bin file. This will allow us to visualise the interactions in the context of the genome. We will merge twice, once for each of the interacting bins.
bin_index_short = cbind.data.frame(bin_index$Start, bin_index$Bin)

colnames(bin_index_short) = c("Start","Bin")

data_1 = na.omit(merge(Iced.matrix, bin_index_short, by.x = "Bin1", by.y = "Bin", all.x=TRUE, all.y=TRUE))
colnames(data_1) = c("Bin1", "Bin2", "Interaction_freq", "Start1")

data_2 = na.omit(merge(data_1, bin_index_short, by.x = "Bin2", by.y = "Bin", all.x=TRUE, all.y=TRUE))
colnames(data_2) = c("Bin2", "Bin1", "Interaction_freq", "Start1", "Start2")

HiC_data = na.omit(cbind.data.frame(data_2$Start1, data_2$Start2, data_2$Interaction_freq))
colnames(HiC_data) = c("chr4", "chr4.1", "Interaction_freq")

### Now we have a matrix with the genomic coordinates of the interacting bins and the interaction frequency. We can visualise this with plotgardener.
{
  pageCreate(width = 8, height = 6, default.units = "inches")
  
  ## Set genomic and dimension parameters in a `params` object
  params_c <- pgParams(chrom = "chr4", chromstart = 37000000, chromend = 38978218, 
                       assembly = "hg19", x = 0.5, width = 7, default.units = "inches")
  
  hicPlot = plotHicTriangle(params = params_c, data = HiC_data,  resolution = 10000, zrange = c(0, 100),
                            palette = colorRampPalette(brewer.pal(n = 9, "Blues")),
                            y = 1, height = 3, just =c("left", "top"))
  
  annoHeatmapLegend(plot = hicPlot, x = 8, y = 2.2, width = 0.13, height = 1.2, just = c("right", "top"))
  
  plotText(label = "Erythroid Tiled-C", fontcolor = "black", fontsize =10,fontface = "bold",rot = 45,
           x = 1.75, y = 2.2, just = c("center"), default.units = "inches")
  
  ## Plot genes
  genes_gm <- plotGenes(params = params_c, stroke = 0.05, fontsize = 8,
                        geneHighlights = data.frame("gene" = c("Hba-a1","Hba-a2"),"color" = c("red", "red")),
                        geneBackground = "#081D58", strandLabels = TRUE, 
                        y = 4.2, height = 0.4)
  
  ## Annotate genome label
  annoGenomeLabel(plot = genes_gm, params = params_c, scale = "bp", fontsize = 10, y = 0.75)
  
  #Remove guide
  pageGuideHide()

}

Visualising in UCSC Genome Browser

From the HiC-Pro outputs, the .matrix and .bed files can be used to create a .hic file which can be visualised in the UCSC genome browser. This is done using a tool from the HiCPro suite under HiC-Pro/bin/utils/hicpro2juicebox.sh. You will need to provide the path to the juicer tools jar file which helpfully is already installed on the HPC cluster. You will also need a chromosome sizes file which can be sourced from the databank.

HICPRO_PATH/bin/utils/hicpro2juicebox.sh \
    -i hic_results/data/sample/sample_allValidPairs \
    -g chrom.sizes \
    -j /package/juicebox/2.10.01/lib/juicer_tools.2.10.01.jar \
    -o output_dir

The resulting .hic file can then be uploaded to the UCSC genome browser as a custom track. For example loading it as a custom track in the UCSC genome browser would look like this:

track type=hic name="My Hi-C" bigDataUrl=http://datashare/address/path/to/data.hic