HI-C数据分析 pipeline（二：数据格式转换）

4 个月前

肆灵

由于HIC数据目前很难用一个软件就能将全套流程全部分析完，从数据预处理到后续的TAD计算，loop计算等需要用到众多软件，但是在数据格式方面，各大软件又不一致，所以要进行很多的数据格式转换，本文整理一下HicPro数据结果如何对接其他Hi-C后续分析。

一、hicpro得到的数据

用HicPro跑完HI-C流程后，会得到两个重要的文件：

hic_results/data/ 文件夹下的 allValidPair 文件，该文件可以转化为juicer软件的 .hic 文件进行juicer后续流程分析。

hic_results/matrix/ 文件夹下的sparse互作矩阵，包含两个文件，一个abs.bed文件，一个matrix文件。

二、Sparse matrix to density/cworld matrix

hicpro的输出矩阵很特别，一个abs矩阵包含了基因组切割区域及编号，另外一个matrix文件一共三列，每一行记录一个相互作用。

1、hicpro matrix 数据格式如下：

## abs.bed文件格式
$ head -n 4 FL12_100M_40000_abs.bed
## chr  start  end  id(编号)
chr1       0       40000   1
chr1       40000   80000   2
chr1       80000   120000  3
chr1       120000  160000  4
## matrix文件格式
$ head -n 4 FL12_100M_40000.matrix
## id1  id2  contacts
76      76      286
76      77      138
76      78      62
76      79      71

2、density数据格式行名列名均代表染色质位点，矩阵内容为各位点之间的contact，行数等于列数，可使用hicpro自带的sparseToDense.py 进行转化得到。

## 该命令位于 ~/HiC-Pro_3.1.0/bin/utils/sparseToDense.py
usage: sparseToDense.py [-h] [-b BINS] [-g ORG] [-d] [-i] [-c] [-o OUTPUT] 
                        filename
sparseToDense.py -b merge_abs.bed -o ~/output/merge_density.matrix --perchr
## 使用 --perchr 参数可以输出单个染色体的density 矩阵
## 如果后续矩阵用于计算 DI value，可带上-d参数

3、hicpro 转化成 cworld-dekker 格式，可以用于计算insulation score，inputfile为第2步得到的 perchr density.matrix裸矩阵，即没有加 -d参数进行转换的矩阵。

step1: get cworld header from abs.bed

## step 1: get matrix header for cworld.
## input file : _abs.bed file from hicpro
## cworld header 有两个，行的和列的。
## header 格式如下：
$head -n 4 cworld_40000_x.chr1.header
REV_1|mm10|chr1:0-40000
REV_2|mm10|chr1:40000-80000
REV_3|mm10|chr1:80000-120000
REV_4|mm10|chr1:120000-160000
$head -n 4 cworld_40000_y.chr1.header
FOR_1|mm10|chr1:0-40000
FOR_2|mm10|chr1:40000-80000
FOR_3|mm10|chr1:80000-120000
FOR_4|mm10|chr1:120000-160000
## 直接shi用cat awk sed 将abs.bed转换成相应格式
for i in {`seq 1 19`,X}
cat merge_abs.bed |awk '{if($1=="'${i}'") print $0}' | awk 'BEGIN{OFS=""}{print "REV_",$4,"|mm10|","'${i}'",":",$2,"-",$3}' > cworld_header_x.${i}.header
sed 's/REV/FOR/' cworld_header_x.${i}.header > cworld_header_y.${i}.header
done

step 2 : add header to hicpro density matrix，使用cwold-dekker软件的addMatrixHeaders.pl 命令给裸矩阵加header

## step 2 : add header to hicpro density matrix
## 给第2步中转化好的density矩阵加header，逐个染色体分别计算
cworld-dekker/scripts/perl/addMatrixHeaders.pl -i merge_chr1_matrix.density \
  --xhf cworld_header_x.${i}.header \
  --yhf cworld_header_y.${i}.header \
  -o merge_chr1_cworld

三、HicPro —— h5/juicer/homer/cooler/2d-text 数据转换

1、使用HicPro软件自带的 hicpro2juicebox.sh 将allValidPair转换为.hic文件, hic文件可以使用juicbox进行可视化。

hic2juice=HiC-Pro_3.0.0/bin/utils/hicpro2juicebox.sh
genomesize=/path/to/Mus_musculus.GRCm38.dna.toplevel.chrom.size
juicertool=/home/edith/juicer/scripts/juicer-1.6/CPU/common/juicer_tools.jar
$hic2juice -i /path/to/sample.allValidPair -g $genomesize -j $juicertool -o /path/to/output

2、使用HiCexplorer软件的 hicConvertFormat 将.hic/ hicpro matrix文件转换成其他格式。

usage: hicConvertFormat --matrices MATRICES [MATRICES ...] --outFileName
                        OUTFILENAME [OUTFILENAME ...] --inputFormat
                        {h5,cool,hic,homer,hicpro,2D-text} --outputFormat
                        {cool,h5,homer,ginteractions,mcool,hicpro}
                        [--correction_name CORRECTION_NAME]
                        [--correction_division] [--store_applied_correction]
                        [--chromosome CHROMOSOME] [--enforce_integer]
                        [--load_raw_values]
                        [--resolutions RESOLUTIONS [RESOLUTIONS ...]] [--help]
                        [--chromosomeSizes txt file] [--version]
                        [--bedFileHicpro BEDFILEHICPRO [BEDFILEHICPRO ...]]

注意使用.hic文件进行格式转换时，需要指定-r参数，因为 .hic 文件中包含所有分辨率的矩阵，如果不指定 -r，转换成cool格式时默认转成.mcool文件，然而hicpro的默认输出分辨率包括40kb, cooler似乎不支持40kb分辨率，不指定 -r 后续可能会出现报错。

例:

hicConvertFormat -m merge.allValidPairs.hic --inputFormat hic --outputFormat cool -r 100000 -o merge.100k.cool
## 如果input文件是hicpro，则需要额外指定 abs bed文件
## 挨个转换一次
for format in {cool, h5, homer, ginteractions}
hicConvertFormat -m merge_iced.matrix -bf merge_abs.bed --inputFormat hicpro --outputFormat $format -o outdir/${format}/merge_100k_matrix.${format} -r 100000
done

四、HIC文本数据格式总结

1、hicpro格式，2部分，上面亦有总结：

## abs.bed文件格式
$ head -n 2 FL12_100M_40000_abs.bed
## chr  start  end  id(编号)
chr1       0       40000   1
chr1       40000   80000   2
## matrix文件格式
$ head -n 2 FL12_100M_40000.matrix
## id1  id2  contacts
76      76      286
76      77      138

2、density格式裸矩阵，即没有header的矩阵，行列相等的大型矩阵，类似以下4X4的矩阵：

3、Dixon et.al. 计算DI的矩阵格式：

chr15   0       40000      0.0     0.0     0.0     0.0
chr15   40000      80000      0.0     0.0     0.0     0.0
chr15   80000      120000      0.0     0.0     0.0     0.0
chr15   120000     160000      0.0     0.0     0.0     0.0

4、cworld-dekker矩阵，为裸矩阵加了header以后的格式

## cworld::dekker
## Dekker Lab
## Contact:     Bryan R. Lajoie
## https://github.com/blajoie
## Version:     0.41.1
## Date:        03:05:17 pm, 06/14/2022
## Host:        cas592.pi.sjtu.edu.cn
## Tool:        addMatrixHeaders.pl
## inputMatrix = '/lustre/home/acct-medhdl/medhdl/data/hic/output-100M/density/FL12_100M/iced/FL12_100M_500000_matrix.density'
## output = '/lustre/home/acct-medhdl/medhdl/data/hic/output-100M/density/FL12_100M/iced/FL12_100M_500000'
## verbose = '0'
## xHeaderFile = '/lustre/home/acct-medhdl/medhdl/data/hic/header_file/cworld_500000_x.header'
## yHeaderFile = '/lustre/home/acct-medhdl/medhdl/data/hic/header_file/cworld_500000_y.header'
        REV_1|mm10|chr1:0-500000        REV_2|mm10|chr1:500000-1000000        REV_3|mm10|chr1:1000000-1500000        REV_4|mm10|chr1:1500000-2000000
FOR_1|mm10|chr1:0-500000        1        0        0        0
FOR_1|mm10|chr1:500000-1000000        0        1        0        0