Hi-C fastq data to valid pair and matrix data format transition.
本文总结了Hi-C原始数据处理流程,以及不同Hi-C数据格式之间的转换。
一、比对
1、工具安装
工具:Hi-CPro 软件: https://github.com/nservant/HiC-Pro
选择最新的安装包下载,按照官网的步骤安装软件及依赖软件。
tar -zxvf HiC-Pro-master.tar.gz
cd HiC-Pro-master
## Edit config-install.txt file if necessary
make configure
make install
2、建立project文件夹(重要)
文件夹一定要按照软件要求的格式创建
tree rawdata
rawdata/
├── CT
│ ├── CT-1_R1.fq.gz
│ ├── CT-1_R2.fq.gz
└── Treat
│ ├── Treat-1_R1.fq.gz
│ ├── Treat-1_R2.fq.gz
3、创建酶切位点文件
如果使用HindIII酶切,则软件中已经有预置,不需要创建。我这里用的是MboI酶切,酶切位点为GATCGATC,需要重新创建。流程如下:
## 需要参考基因组fasta文件
pwd
~/HiC-Pro_3.1.0/annotation ## 需要在annotation文件夹下创建bed
~/HiC-Pro_3.1.0/bin/utils/digest_genome.py -r GATCGATC -o Mus_musculus_MboI.GRCm38.dna.toplevel.bed Mus_musculus.GRCm38.dna.toplevel.fa
4、修改config-hicpro.txt文件
N_CPU = 40
## 比对软件的路径
BOWTIE2_IDX_PATH = /your/bowtie2/index/dir_path
## 参考基因组的路径
REFERENCE_GENOME = Mus_musculus.GRCm38.dna.toplevel
GENOME_SIZE = Mus_musculus.GRCm38.dna.toplevel.chrom.size
## 基因组切割bed文件,及酶切位点碱基序列,这里使用的MboI酶切割
GENOME_FRAGMENT = Mus_musculus_MboI.GRCm38.dna.toplevel.bed
LIGATION_SITE = GATCGATC
## 其他参数可以根据自己的情况修改
5、运行HicPro
runhic=/home/edith/Downloads/HicPro/HiC-Pro_3.1.0/bin/HiC-Pro
input=/home/edith/data/hic/rawdata
output=/home/edith/data/hic/output
conf=/home/edith/data/hic/output/config-hicpro.txt
$runhic -i $input -o $output -c $conf
6、输出文件
hicpro 输出两个文件夹:
bowtie_results/ :
output/bowtie_results/bwt2/ 文件夹下包含比对的最终bamfile文件;
bowtie_results/bwt2_glob/ 和 bowtie_results/bwt2_loc/文件夹可以移;
hic_results/:
hic_results/data 包含 valid pair文件
hic_results/matrix 包含HiC pro格式的稀疏矩阵文件
hic_results/pic 包含输出的质控图片;
二、获得个性化的HiC矩阵
1、HicPro分步运行:
HiCPro的处理流程包括以下6部分:
INPUT DATA TYPE IN STEPWISE MODE
-s mapping .fastq(.gz) files
-s proc_hic .bam files
-s quality_checks .bam files
-s merge_persample .validPairs files
-s build_contact_maps .validPairs files
-s ice_norm .matrix files
如果想 HiCPro 分步运行,只需要命令后面加上-s参数:
$runhic -i ${output}bowtie_results/bwt2/ -o ${output} -c $conf -s proc_hic
$runhic -i ${output}bowtie_results/bwt2/ -o ${output} -c $conf -s quality_checks
$runhic -i ${output}hic_results/data/ -o ${output} -c $conf -s merge_persample
$runhic -i ${output}hic_results/data/ -o ${output} -c $conf -s build_contact_maps
$runhic -i ${output}hic_results/matrix/ -o ${output} -c $conf -s ice_norm
2、HiCPro文件的merge
output/bowtie_results/bwt2/ 文件夹下的 .bam 文件可以修改,如果想将多个HiC文件进行进行merge分析,可以在比对后将bamfile 全部放到一个bwt2/文件夹下,然后运行后续的HiCpro分析内容:
## 将比对文件merge 后文件夹内容如下:
tree bwt2/
├── CT
│ ├── CT-1_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
│ ├── CT-2_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
│ ├── CT-3_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
├── Treat
│ ├── Treat-1_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
│ ├── Treat-2_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
│ ├── Treat-3_Mus_musculus.GRCm38.dna.toplevel.bwt2pairs.bam
## 然后运行hicpro
$runhic -i ${output}bowtie_results/bwt2/ -o ${output} -c $conf -s proc_hic
$runhic -i ${output}bowtie_results/bwt2/ -o ${output} -c $conf -s quality_checks
$runhic -i ${output}hic_results/data/ -o ${output} -c $conf -s merge_persample
$runhic -i ${output}hic_results/data/ -o ${output} -c $conf -s build_contact_maps
$runhic -i ${output}hic_results/matrix/ -o ${output} -c $conf -s ice_norm
3、ValidPair的标准化
利用分步运行,可以在运行完 pro_hic 步骤后获得validpair文件,随后可以修改validpair文件的大小,将HiC数据标准化再进行下游分析。
## 利用shuf 命令随机抽取100000000 个readpair
shuf -n100000000 filename.validPair > filename-100M.validPair
## 然后运行后续的HICPRO 命令
$runhic -i ${output}hic_results/data/ -o ${output} -c $conf -s build_contact_maps
$runhic -i ${output}hic_results/matrix/ -o ${output} -c $conf -s ice_norm
随后,HiCPro的输出文件hic_results/data/文件夹下面的validPair文件可以被转化成.hic文件,在Juicer中可视化;
hic_results/matrix/ 文件夹下面的 _iced.matrix 文件可以使用 HicPro自带的工具转化成 density 矩阵格式,也可以使用其他软件比如 hicexplorer 进行格式转化。下篇笔记整理一下我的 hic数据格式转换的流程。
欢迎关注公,一起交流!
原文:HI-C数据分析 pipeline(一:上游数据预处理