陈作舟张俊芳
(上海海洋大学水产与生命学院 上海 201306)
一个Chip-seq的生物信息分析流程
陈作舟*张俊芳
(上海海洋大学水产与生命学院 上海 201306)
本文简单的介绍了一个Chip-seq 的生物信息分析流程的设计和实现。
生物信息 Chip-seq 基因组
随着二代测序技术的持续优化,测序成本大幅度降低,生命科学的各个领域发展了针对二代测序的各种应用,其中有一类就是通过免疫沉淀富集基因组的片段,使研究者得以研究基因组的不同状态之间的差异,例如不同组织基因表达的差异是如何实现的,癌症和正常细胞的基因组状态有何差别,细胞的基因组是如何适应环境温度变化的。Chip-seq是通过染色质免疫共沉淀技术(ChIP)富集目的蛋白结合的DNA片段,继而对富集得到的DNA片段进行高通量测序。目前已经有很多工具能够对这类数据进行分析,各有优缺点,并没有形成统一的模式,例如有的仅仅提供在线分析,有的仅仅提供某些核心环节的分析,为此,我们建立了一个Chip-seq的分析流程,该流程整合了若干生物信息工具以及若干R语言包,现简述如下,以供生命科学的研究人员参考。
2.1 测序数据与基因组比对
假定我们已经得到斑马鱼的两组转录组因子的Chip-seq测序数据,A和B,以及它们的未进行免疫沉淀测序的对照结果序列文件(control)inputAB。 A和B可以预先通过FastQC或FastX等工具来控制数据质量。然后我们将A和B对斑马鱼的基因组进行比对,比对工具有Bowtie/Bowtie2, BWA和STAR等。这里我们用Bowtie2举例说明:
bowtie2 -p 4 -x dr.genome -U A.fatsq S A.sam
bowtie2 -p 4 -x dr.genome -U B.fatsq S B.sam
bowtie2 -p 4 -x dr.genome -U inputAB.fatsq S inputAB.sam
其中-p 代表使用的CPU核心数量, -x代表对应的基因组, -U代表输入的fastq序列文件, -S代表输出的比对结果,该结果为SAM格式。SAM格式可通过Samtools和Bedtools等工具转化成Bed格式。
2.2 得到免疫沉淀的峰文件(Peak Calling)
有多个工具可以执行Peak Calling工作,比较常用的有MACS/MACS2、SICER等,在这里我们以MACS举例说明。
macs14 -t A.bed -c inputAB.bed -f BED -g dr -n A --keep-dup=1
macs14 -t B.bed -c inputAB.bed -f BED -g dr -n A --keep-dup=1
其中-t 代表前一步得到的基因组比对(alignment)文件,-c代表免疫沉淀的control文件,-f BED表示输入文件的格式为BED,-g代表基因组的类型,这里用的是斑马鱼,-n代表输出文件名,--keep-dup代表重复的测序计算的次数。
2.3 合成峰(Peak Merging)
Bedtools merge -i AB_peaks >AB_peaks.merged
利用Bedtools软件包中的merge功能,将A和B的两组峰合成一组。以合成的峰组作为一个公共的可比较的对象(Reference)来进行后续分析。
2.4 分析1:覆盖情况分析及文氏图
Bedtools coverage -a AB_peaks.merged -b A_peaks.bed>A_coverage
Bedtools coverage -a AB_peaks.merged -b B_peaks.bed>B_coverage
利用Bedtools软件包中的coverage功能,计算A和B对Reference的覆盖情况。可以将得到的结果利用R语言自带的或者第三方的文氏图相关软件包进行作图,例如“VennDiagram”,“Vennerable”等。
该分析适用于以下几种情形:
(1)两个或多个具有潜在相关性的DNA结合蛋白,例如转录因子,我们需要研究它们的相关性情况,为相互作用提供证据(此为Chip-seq分析)。(2)一个DNA结合蛋白在细胞不同状态的结合情况(此为Chip-seq分析)(3)以上两种情况的结合(此为Chip-seq分析)(4)以上的DNA结合蛋白更换为组蛋白修饰(此为Chip-seq分析)
2.5 分析2:差异分析
以上覆盖情况分析的着重点在于不同库(library)的峰在的基因组位置上的异同,从而为寻找出基因组的不同状态之间的生物学差异提供线索,而差异分析目的是进一步得到不同库的公共的峰的相对表达量的差异。先利用第4步得到的合成峰作为参考区域,然后计算参考区域的Reads覆盖情况,均一化(Normalize)以后进行统计分析,一般如果没有重复(Replica)的话,使用Fishers' Exact test或Chi-square test,如果有重复,则使用专用的R语言软件包如edgeR等。
2.6 分析3:差异分析后的Gene Ontology(GO)富集分析
我们利用GO.db这一R语言模块检索GO的上下级关系。富集检测利用Fishers' Exact test。
2.7 分析4:通过多Peak的联合分析发现整体上调或下调的通路
先将斑马鱼的基因注释到KEGG通路上,然后利用Wilcoxin Rank test检测出差异分布的通路。
2.8 分析5:峰的注释及分析
利用UCSC的斑马鱼基因组对峰进行注释,一般根据峰和基因组元件(如启动子区域)的重叠情况进行注释。然后得到峰在不同基因组元件上的分布情况。
2.9 分析6:免疫沉淀峰(Chip-seq)与基因表达谱(RNA-seq)的联合分析
免疫沉淀峰和基因表达谱的联合分析有多种类型,我们主要对不同组织的峰的组合所对应的基因表达情况进行分析,例如两种转录因子组合所在区域的基因表达是否提高。
[1]Park PJ. Chip-seq: advantages and challenges of a maturing technology. Nat Rev Genet.2009(10):669-680.
[2]Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010(10):106.
Q753
A
1674-2060(2016)03-0333-01
感谢:本文受上海市青年教师资助计划(项目编号:ZZHY13001,陈作舟)、上海市人才发展资金(项目编号:201457,张俊芳)和上海市上海高校高峰高原学科建设计划资助。
张俊芳,(1976—),女,山西太原人,博士,教授,基因组学与表观遗传学。
陈作舟(1979—),男,浙江杭州人,博士,高级工程师,主要从事生物信息分析工作。