ChIP-Seq分析之ChIPQC结果含义

刘小泽写于2020.6.20
首先是ChIP-Seq分析的前言介绍部分：
1：了解ChIP-seq的实验流程
 2：继续了解ChIP-Seq
3：关于ChIP-Seq的实验对照与偏差来源
 4：ChIP-Seq的实验设计补充
 5：ChIP-Seq数据库及实战数据介绍
然后开始实战部分：
6：ChIP-Seq计算资源准备与实战数据下载
 7：ChIP-Seq数据质控和过滤
 8：ChIP-Seq数据比对注意事项
 9：ChIP-Seq数据比对实战
 10：使用ChIPQC进行质控

这次将介绍ChIPQC结果的各个指标含义
内容来自：https://hbctraining.github.io/Intro-to-ChIPseq/lessons/06_combine_chipQC_and_metrics.html

首先回顾一下ChIPQC的运行

library(ChIPQC)

## Load sample data
samples <- read.csv('meta/samplesheet_chr12.csv')
View(samples)

## Create ChIPQC object
chipObj <- ChIPQC(samples, annotation="hg19")
# ChIPQC需要调用BiocParallel，windows用户可能报错，可以先运行一句：register(SerialParam())

## Create ChIPQC report
ChIPQCreport(chipObj, reportName="ChIP QC report: Nanog and Pou5f1", reportFolder="ChIPQCreport")

示例结果可以看：https://www.dropbox.com/s/sn8drmjj2tar4xs/ChIPQCreport - full dataset.zip?dl=1

结果展示

首先映入眼帘的是总体报告：

其中有一些指标见过，一些没见过，像是SSD, RiP and RiBL这几列就是ENCODE计划提出的指标。就是评估了信号的分布，在富集区域、整个基因组、已知的artefact regions

总而言之，这些全部的指标可以分为4类：

Read characteristics
Enrichment of reads in peaks
Peak signal strength
Peak profiles

切记：通过这些指标也不意味着实验就是成功的，相反亦然

1 关于测序reads的指标（Read characteristics）

包括了：read depth, read length, duplication rate

如果read depth, read length在样本间差异很大，就要引起注意了

由于之前已经过滤掉了重复reads，所以这里的duplication rate没什么用

2 reads在peaks中的富集情况（Enrichment of reads in peaks）

包括了RiP, SSD, and RiBL

2.1 RiP (Reads in Peaks)

也叫FRiP，表示：the percentage of reads that overlap ‘called peaks’，也就是peaks包含的reads数占reads总数的百分比

可以理解成：信噪比（signal-to-noise）

根据感兴趣蛋白（POI，protein of interest）的不同，RiP值也差异较大：

质量好的转录因子(sharp/narrow peaks)一般得到5%以上的RiP值
质量好的Pol2 (mix of sharp/narrow and dispersed/broad peaks)得到30%以上的RiP
也有质量不错的数据集但RiP很小，小于1%（例如 RNAPIII or a protein that binds few sites）

上面图中看到，Nanog比Pou5f1的RiP值要高，而Pou5f1-rep2更是低的可怜，可以说明的是：Nanog样本富集效果更好

有两张图可以反映：

不过看箱线图发现，虽然Nanog的RiP较高，但这个分布和Pou5f1也相差太远，推测可能与read length 、 depth有关

2.2 SSD（standard deviation of signal）

表示基因组中信号值的标准差，可以反映reads在基因组中覆盖度一致性，越大越离散，就是高的越高，低的越低

我们希望看到：IP样本中这个值较大，说明富集区域信号很强，非富集区域信号较弱，因此它的标准差很大；而control样本最好就是标准差较小，不要有太大的波动

SSD值高虽然说明有的区域信号强，但不一定是ChIP的富集区域，一些blacklist区域也会存在较强的信号

【关于blacklist：】

我们这里的数据显示：Pou5f1比Nanog的SSD值要高，可能说明Pou5f1的富集效果更好，但不能确定，因为还需要确定Pou5f1的SSD高不是由于未知的artifact造成

有一张图可以反映：Coverage histogram

x轴表示每个碱基位置的reads堆积高度（也就是信号值、测序深度）【联想一下IGV中reads堆积在一起，向下延伸，反过来就是一个个的山峰】
y轴表示有多少个位置有这样的堆积高度（log后的）

好的富集结果一般是：有一条尾巴（依然存在很多位点具有较高的测序深度）；而像input样本这种低富集的，主要是包含背景，因此它的y轴很高，同时x轴很低

我们这个数据集中，尤其是Nanog rep2样本，具有更粗壮的尾巴（Heavy tail，意思就是在曲线以下具有更大的空间）。Nanog样本具有更多高深度的位点

综合考量：

Pou5f1的coverage不如Nanog，但SSD高于Nanog。说明Pou5f1存在某一块区域深度较高，但不是整体都高，可能存在blacklist区域

那么是否真的存在blacklist区域呢？还有再看一个指标：

2.3 RiBL (Reads overlapping in Blacklisted Regions)

也就是与已知blacklist有交集的reads占比。这个值越低越好

黑名单区域一般也是唯一比对，因此常规的去重复操作对它无效。这些区域一般是：着丝粒，端粒和卫星重复序列

黑名单区域的危害是：confound peak callers and fragment length estimation，因此需要追踪并去除比对到这些区域的reads

我们的数据中，RiBL的比例看上去还比较合理，并没有出奇的高。因此高SSD可能是因为存在更多容易破碎的开放染色质区域，或者存在hyper-ChIPable区域，与很多不相干的蛋白也能产生富集，导致假阳性

当然，如果在peak calling之前去掉了黑名单区域，就没必须分析RiBL了

3 peak信号强度（Peak signal strength）

主要包括：FragLength and RelCC（又称Relative strand cross-correlation coefficient or RSC）

一般，RelCC在所有的ChIP样本中大于1，表示具有较高的信噪比；FragLength也应该与文库制备过程中设定的片段长度接近

3.1 Strand cross-correlation

一个高质量的ChIP实验，会在POI附近形成非常显著的reads富集，会在正负链发现双峰分布

Cross-Correlation scores的计算：Pearson’s linear correlation between coverage for each complementary base. These Pearson correlation values are computed for every peak for each chromosome and values are multiplied by a scaling factor and then summed across all chromosomes，就是先在正负链生成两个向量，表示某个碱基位点的reads数量，然后求这两个向量的相关性，并逐渐沿着shift size移动，最后得到一个相关性表

最后这个cross-correlation值算好，就会画在y轴上，x轴就表示shift size

一般这个cross-correlation plot会产生两个峰：

a peak of enrichment corresponding to the predominant fragment length (highest correlation value)
a peak corresponding to the read length (“phantom” peak)

我们这里的数据中，Nanog 和 Pou5f1都能看到两个峰：

看Y轴：Nanog中最大相关性的值比Pou5f1要高，表示Nanog中信号的量更大【峰两侧reads数越多，信号越强，计算的相关性值也就越高】
再看对应的x轴，当cross-correlation 值最大时，对应的 shift size可以评估片段长度

3.2 RelCC（又称RSC）

它就是根据cross-correlation的最大、最小值计算的

RSC值低可能是由于ChIP的质量差、测序reads质量差导致错配多、测序深度不够【其实可以理解为：RSC值低=》就是相关性计算的值低=》正负链没有足够的reads =》也就是上述原因】

另外，数据集的结合位点太少（比如小于200）也会导致低的RSC【这个也很好理解，结合位点少，更别提位点正负链富集的reads数量了】。结合位点少的原因可能是生物因素（比如某一个因子在某一个特定组织中就这么几个位点）

Cross-Correlation Plots的例子

强信号：

下面这个例子是人类细胞的CTCF 转录因子(zinc-finger transcription factor)。使用一个好的抗体，转录因子一般会富集45,000 - 60,000个peaks。红线表示真正的peak，蓝色线表示read length

弱信号：
抗体不是特别有效，得到的峰也比较分散，在185-200bp间存在真的峰，另一个蓝色则是read length。对于弱信号的数据，read-length peak将占据主导地位

没有信号：
表示实验失败或者input样本，基本看不到fragment length这个峰
也就是在特定的结合位点附近，没有富集到reads

4 peak的类型

4.1 Relative Enrichment of Genomic Intervals (REGI)

将peaks与基因组注释结合起来，看看reads主要富集在哪些区域

我们的数据中，“Promoters500” and “All5UTRs”的富集程度最高，也符合预期（Nanog和Pou5f1作为转录因子应该结合在这块区域）

4.2 peak的形状（Peak Profile）

这个形状根据抗体的类型存在差异：transcription factor, histone mark, or other DNA-binding protein such as a polymerase

总结

样本内多个重复比较：按说应该是类似的变化趋势，如果出现较大差异，要留意
样本间比较：比如Pou5f1的SSD和RelCC值高，好像说明富集效果好，但coverage plots、cross-correlation plots又比Nanog低
低质量数据的来源可能是：
- Strength/efficiency and specificity of the immunoprecipitation，ChIP的质量实际上就是抗体的特异性与富集程度体现，如果对预定目标反应性差或使用了非特异性抗体，引起与其他DNA相关蛋白的交叉反应，也会效果不好
- Fragmentation/digestion：超声处理的方式可能导致不同的片段大小分布。因此如果input没有和IP样本一起超声处理，那么最好不用

欢迎关注我们的公众号～_～　　
我们是两个农转生信的小硕，打造生信星球，想让它成为一个不拽术语、通俗易懂的生信知识平台。需要帮助或提出意见请后台留言或发送邮件到jieandze1314@gmail.com

Welcome to our bioinfoplanet!