HiChIP 数据分析: 用HiC-Pro预处理原始数据

用HiC-Pro预处理原始数据

HiC-Pro 是用于处理 Hi-C 数据的 pipeline，它从测序 reads（FASTQ 文件）开始，执行多个步骤，如 alignment、filtering、binning 和 normalization。

最终输出文件是 raw 和 normalized contact matrices，以及关于输入数据质量的附加信息。

hichipper以 HiC-Pro 生成的 aligned 和 filtered read pairs 作为输入，我们只需要执行 alignment 和 filtering 步骤。

在继续 HiC-Pro 之前，还建议对 raw data 进行 standard quality controls，并在必要时进行 trimming，如同在其他任何测序实验中所做的那样。

Input

HiC-Pro 要求原始测序文件被组织在每个样本的子目录中。为此，我们在 fastq/ 文件夹内为每个样本创建一个文件夹，并按如下方式移动相应的 FASTQ 文件：

mkdir fastq/Rad21_Rep1
mv fastq/Rad21_Rep1_*.fastq.gz fastq/Rad21_Rep1/

使用星号 * 可以一次性移动同一样本的两条 mate 文件。

要运行 HiC-Pro，需要提供四种注释文件：

1. Bowtie2 indexes
2. 染色体大小的表
3. 酶切后 restriction fragments 坐标的文件
4. 一个配置文件。

所有这些文件必须与所选参考基因组相同版本相对应，在我们的例子中是 hg19。HiC-Pro 已经提供了部分所需注释文件，例如 hg19 的 chromosome sizes 以及几种 restriction enzymes 的 restriction fragments。

如果在 HiChIP 实验中使用的酶未被包含在内——如我们使用的 DpnII——可以使用位于 HiC-Pro-2.11.1/bin/utils 文件夹中的 digest_genome.py 脚本生成 restriction fragment 文件。digest_genome.py 把限制酶识别的序列作为输入（DpnII 为 ^GATC，其中 ‘^’ 表示切割位点），并将参考基因组的 FASTA 序列作为输入，返回 restriction fragment 坐标文件：

/home/Programs/HiC-Pro-2.11.1/bin/utils/digest_genome.py -r
^GATC -o DpnII_resfrag_hg19.bed genome.fa

生成的 DpnII_resfrag_hg19.bed 文件需要移动到 HiC-Pro 的 /annotation 文件夹中，这样在处理中就会被自动找到：

mv DpnII_resfrag_hg19.bed /home/Programs/HiC Pro_2.11.1/annotation/

Config

HiC-Pro 的配置文件包含了所有必需的设置和指向注释文件的路径。就像 restriction fragment 文件一样，HiC-Pro 的安装文件夹里提供了一个预编译的配置文件，名为 config-hic-pro.txt。

尽管大多数参数可以保留默认值，但其他参数必须根据用户设置进行修改，才能正确运行分析。特别需要修改的有：

• 在 DATA 部分，我们需要指定每条 mate 的 FASTQ 文件是如何命名的，在我们这里为 _1 和 _2。
• 在 ALIGNMENT 部分，必须提供 Bowtie2 indexes 的路径。
• 在 ANNOTATION file 部分，需要指定参考基因组的名称，该名称必须与 Bowtie2 索引文件前缀保持一致。
• 在 DIGESTION 部分，需要指明 restriction fragment 文件的名称（无需完整路径，因为 HiC-Pro 会在其自带的 annotation/ 文件夹中查找）和 ligation site。在我们的例子中，ligation site 为 GATCGATC。

配置文件中更新后的字段如下：

PAIR1_EXT = _1
PAIR2_EXT = _2
BOWTIE2_IDX_PATH = /home/Annotation/Homo_sapiens/UCSC/hg19/
Sequence/Bowtie2Index
REFERENCE_GENOME = genome
GENOME_FRAGMENT = DpnII_resfrag_hg19.bed
LIGATION SITE = GATCGATC

修改后的配置文件可以保存为 config-HiChIP.txt 并放在我们的工作目录中。

Alignment

HiC-Pro 可以在一次运行中处理所有样本，也可以按顺序模式（sequential mode）运行，要求用户通过 -s 选项逐一调用每一步的子集。因为我们希望在每一步都检查结果质量，并且我们不需要跑完整条 pipeline（不会生成 contact matrices），所以我们利用 HiC-Pro 的顺序模式。对于 alignment 步骤（-s mapping），我们使用如下命令：

HiC-Pro=/home/Programs/HiC-Pro-2.11.1/bin/HiC-Pro

$HiC-Pro -c config-HiChIP.txt -i fastq -o HiC_Pro -s mapping -s quality_checks

这些命令会创建名为 HiC_Pro/ 的输出文件夹，其中比对结果位于 bowtie_results/ 子文件夹内。

通过 -s quality_checks，HiC-Pro 会再创建两个子文件夹（hic_results/stats/ 和 hic_results/pic/），用于存放所有比对统计信息并以图形方式展示；所有这些结果文件夹还会进一步按 fastq/ 中发现的每个样本再分出一个子文件夹。

HiC-Pro 采用基于 Bowtie2 的两步比对策略。首先，将两条 mate 各自比对到基因组上的 reads 保存到一个 BAM 文件中。接着，未比对上的 reads 会按配置文件中指定的 ligation site 做 trimming，然后再次用 Bowtie2 对这些被修剪后的 reads 进行比对，结果存入第二个 BAM 文件。最后，两个 BAM 文件合并并保存在 bwt2 文件夹中。由于两条 mate 是独立比对的，统计信息也分别针对每条 mate 进行报告。

如上表所示，绝大多数 reads 都能比对到基因组（平均超过 90%），而平均约 10% 的 reads 是在 trimming 后比对上的。这意味着约 10% 的 reads 是 chimeric，即跨越了连接接头。chimeric reads 的比例取决于样本制备阶段所选片段大小以及 read 长度。如果未比对上的 reads 数量很高（超过 20%），则可能提示文库存在问题（如 read-through adapters 或污染），或者向 HiC-Pro 输入了错误的 ligation site。

鉴于当前的比对结果，我们可以继续进行 filtering。

未完待续，欢迎关注！