Harvard FAS Informatics出品的ATAC-seq测序指南

Posted

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了Harvard FAS Informatics出品的ATAC-seq测序指南相关的知识,希望对你有一定的参考价值。

参考技术A github链接: harvardinformatics/ATAC-seq

参考文献: ATAC-seq: A Method for Assaying Chromatin Accessibility Genome-Wide

推荐使用bowtie2进行序列比对(alignment)。常用的参数有以下几个:

默认输出格式问SAM文件,包含了每个序列的比对信息。SAM文件可以压缩为BAM文件,并使用SAMtools进行排序。最佳的分析流程为:

ATAC-seq数据一般含有相当比例的线粒体DNA,详见 该讨论 。可以 使用CRISPR技术去除线粒体基因组污染 。也可以参考最近发表的 Omni-ATAC方法 用去污剂去除线粒体DNA,这种方法对于多数研究人员可能更易于操作,但是不要按照他们的分析流程来分析数据。

不管使用何种实验方案,得到的测序数据中多多少少都会有线粒体DNA。因为线粒体基因组中没有ATAC-seq感兴趣的峰,这些数据会使后续分析变得复杂,因此,推荐在下一步分析之前去除线粒体DNA序列。可以通过下面两种方法中任一种实现:

PCR重复序列是完全一致的reads。

二代测序得到的短序列中,有一些序列会匹配到基因组的多个区域。一些学者在分析数据时会通过samtools view中的-q参数将这些非特异性的序列去除。对于这种非特异性的序列,bowtie2或bwa只会报告一个比对的结果,同时会给一个低匹配质量的(mapping quality, MAPQ)评分,匹配质量定义为: -10 * log10Prmapping position is wrong。可以通过下面命令来去除MAPQ<10的序列。

Model-based Analysis of ChIP-Seq (MACS2)是一款用于检测基因富集区域的软件。尽管是为ChIP-seq设计的,但也适用于ATAC-seq测序及其他有小峰的基因组富集分析。MACS2软件的主程序是callpeak。
可以将前期处理好的比对序列文件作为MACS2的输入文件。然而,一定要记住:比对上的序列只是ATAC技术产生的DNA片段的一部分。因此,必须考虑如何用MACS2来解读这些比对序列。

对于双端测序书来说,比对序列主要分两种类型:成对匹配的和单一匹配的序列。当使用MACS2分析数据时,应事先决定哪种比对序列可用于后续分析。主要有以下三个选择:

输出为BED文件,可以用 MACS2 -f BEDPE 命令处理。

除了前文中描述的一些参数,MACS2还有一些其他参数可以设置。下面是一些可以考虑的参数:

标准的 macs2 callpeak 程序有3个输出文件。如果有 -n NAME 参数,输出文件分布为:NAME_peaks.xls,NAME_peaks.narrowPeak,和NAME_summits.bed。最有用的文件是NAME_peaks.narrowPeak,它是一个bed文档,包含了所有峰的基因组坐标和不同的统计值(倍数改变、p值和q值等)。

一旦完成一组样本的MACS2分析之后,可根据实验设计进行一些后续分析。

可以将peak文件在基因组中可视化。模式生物的ATAC-seq数据,peak文件(NAME_peaks.narrowPeak)可以直接上传到 UCSC genome browser 中进行查看。如果peak文件没有表头的话,可以将下面的内容添加至首行: track type=narrowPeak 。
另外,可以使用IGV进行可视化。peak文件可以直接通过File -->Load from File选项上传。如果要对BAM文件可视化,需要用samtools对BAM文件进行排序并构建索引。

可以使用BEDTools来比较一组peak文件的异同。例如: bedtools intersect 可以比较两个peak文件中相同的峰;寻找两个峰文件中差异的峰,如实验组和对照组,可以通过 bedtools subtract 命令实现。

ChIPseeker 最早是用于注释ChIP-seq数据的,但是也适用于ATAC-seq数据的注释。ChIPseeker可以对基因的位置和其他特征进行注释,详细使用方法见 ChIPseeker的使用指南 。

HOMER可用于寻找motif。可以将peak文件作为HOMER文件的输入,来检测已知的motif和新的motif。

Andrews S. (2010). FastQC: a quality control tool for high throughput sequence data. Available online at: http://www.bioinformatics.babraham.ac.uk/projects/fastqc

Buenrostro JD, Giresi PG, Zaba LC, Chang HY, Greenleaf WJ. Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nat Methods. 2013 Dec;10(12):1213-8.

Buenrostro JD, Wu B, Chang HY, Greenleaf WJ. ATAC-seq: A Method for Assaying Chromatin Accessibility Genome-Wide. Curr Protoc Mol Biol. 2015 Jan 5;109:21.29.1-9.

Corces MR, Trevino AE, Hamilton EG, Greenside PG, Sinnott-Armstrong NA, Vesuna S, Satpathy AT, Rubin AJ, Montine KS, Wu B, Kathiria A, Cho SW, Mumbach MR, Carter AC, Kasowski M, Orloff LA, Risca VI, Kundaje A, Khavari PA, Montine TJ, Greenleaf WJ, Chang HY. An improved ATAC-seq protocol reduces background and enables interrogation of frozen tissues. Nat Methods. 2017 Oct;14(10):959-962.

Heinz S, Benner C, Spann N, Bertolino E, Lin YC, Laslo P, Cheng JX, Murre C, Singh H, Glass CK. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Mol Cell. 2010 May 28;38(4):576-89.

Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012 Mar 4;9(4):357-9.

Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R; 1000 Genome Project Data Processing Subgroup. The Sequence Alignment/Map format and SAMtools. Bioinformatics. 2009 Aug 15;25(16):2078-9.

Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet.journal. 2011;17:10-2.

Montefiori L, Hernandez L, Zhang Z, Gilad Y, Ober C, Crawford G, Nobrega M, Jo Sakabe N. Reducing mitochondrial reads in ATAC-seq using CRISPR/Cas9. Sci Rep. 2017 May 26;7(1):2451.

Quinlan AR. BEDTools: The Swiss-Army Tool for Genome Feature Analysis. Curr Protoc Bioinformatics. 2014 Sep 8;47:11.12.1-34.

Yu G, Wang LG, He QY. ChIPseeker: an R/Bioconductor package for ChIP peak annotation, comparison and visualization. Bioinformatics. 2015 Jul 15;31(14):2382-3.

MGI数据库挖掘 | MGI-Mouse Genome Informatics | InWeb database

做生物信息,遗传发育,分析数据的时候总是要narrow down分析范围,高通量数据尤其是基因表达,在庞大的confounder面前,缩小分析范围是必须的,否则你会一直在混沌中游荡。

 

看一篇文章:2018 - Identification of genes associated with Hirschsprung disease, based on whole-genome sequence analysis, and potential effects on enteric nervous system development

 

Known Genes of ENS Development and Their Interactome
Genes with mutations reported to cause colonic aganglionosis in mutant mice according to the Mouse Genomics database were considered to be known ENS genes. ENS interactome was defined by genes encoding proteins that show protein–protein interaction with known ENS genes (see Supplementary Methods).

 

对于某个特定的疾病,我们总想找到与它直接相关的基因,怎么找呢?

MGI的数据库里就有这种信息,怎么才算直接相关呢?最具逻辑的就是突变数据,如果某个基因上的突变导致了某个疾病,那这个因果关联就肯定存在!!!

 

怎么找呢?看官方教程:

How do I find mutations that cause a specific combination of phenotypes?

 

这个基因list可能包含的基因太少,那如何扩大呢?

另一个数据库就是PPI了,找出与直接基因有蛋白互作的其他的基因,也就是以上方法中的interactome。

文章中庸的InWeb - A scored human protein–protein interaction network to catalyze genomic interpretation

发在NM上了,整合了PPI数据,但是好像很久没更新了,也没仔细去看它是怎么设计的算法。

目前我是直接用的STRING数据库。

 

以上是关于Harvard FAS Informatics出品的ATAC-seq测序指南的主要内容,如果未能解决你的问题,请参考以下文章

CF1151FSonya and Informatics

CF1151F Sonya and Informatics

MGI数据库挖掘 | MGI-Mouse Genome Informatics | InWeb database

CF1151F Sonya and Informatics (计数dp+矩阵优化)

Codeforces Round #775 (Div. 2,based on Moscow Open Olympiad in Informatics) - D. Integral Array - 题解

Codeforces Round #626 (Div. 2, based on Moscow Open Olympiad in Informatics)