7.ChIP-seq

本章主要介绍 ChIP-seq 相关分析方法，包括 peak calling 和 motif analysis。

染色质免疫共沉淀测序（Chromatin immunoprecipitation sequencing，简称 ChIP-seq）被用于分析蛋白质与 DNA 的交互作用。常见的 ChIP-seq 实验主要研究转录因子和组蛋白修饰如何通过与DNA相互作用来调控基因表达，从而影响表型。ChIP-seq 对于充分了解许多生物过程和疾病状态至关重要。

1) Pipeline

2) Data Structure

2a) getting software & data

1. install software (already available in Docker)

   HOMER

2. data
   我们使用 GSE61210 中的两个样本：

   • Input：GSM1499619（input.bam）
   • IP：GSM1499607（ip.bam）

   如上表所示，我们已经准备好 .bam 文件（位于 Docker 中的 /home/test/chip-seq/input）。用户也可以参考 这里 从 CHIP-seq 数据生成自己的 .bam 文件。

2b) input

2c) output

2c.1) 核心输出文件如下所示：

peak table in peak file

# Column Headers:
#PeakID chr     start   end     strand  Normalized Tag Count    focus ratio   findPeaks Score  Total Tags      Control Tags (normalized to IP Experiment)    Fold Change vs Control   p-value vs Control      Fold Change vs Local    p-value vs Local       Clonal Fold Change
chrIII-1        chrIII  78346   78578   +       69987.1 0.862   5971.000000   5966.0   201.6   29.60   0.00e+00        26.54   0.00e+00        0.50
chrIII-2        chrIII  133     365     +       61226.9 0.775   5364.000000   5227.0   116.3   44.93   0.00e+00        59.92   0.00e+00        0.61
chrI-1  chrI    141663  141895  +       41225.6 0.854   3515.000000     3514.0169.1    20.78   0.00e+00        17.09   0.00e+00        0.50
chrII-1 chrII   165145  165377  +       35334.5 0.845   3015.000000     3018.0171.1    17.64   0.00e+00        14.47   0.00e+00        0.50
chrII-2 chrII   555827  556059  +       34817.1 0.790   2973.000000     2970.0159.6    18.61   0.00e+00        10.12   0.00e+00        0.50

图1. motif table in homerResults.html

2c.2) detailed description of Peak calling output

peak file contains the following columns

2c.3) detailed description of Motif analysis output

Detailed output files of Motif analysis will produce many files, we only explain the main output -- homerResults.html in above. Here we will briefly introduce other files.

3) Running Steps

首先进入到容器（在自己电脑的 Terminal 中运行，详情请参见 这里）：

docker exec -it bioinfo_tsinghua bash

以下步骤均在 /home/test/chip-seq/ 下进行:

cd /home/test/chip-seq/

准备输出目录

mkdir output

3a) Peak Calling

The common call peak software are HOMER and MACS, here we mainly introduce the usage of HOMER. You can get manual for MACS here.

HOMER contains a program called findPeaks that performs all of the peak calling analysis. Before we use findPeaks to call peak, we need to convert our .bam file into tag file by using makeTagDirectory:

makeTagDirectory input/ip    input/ip.part.bam
makeTagDirectory input/input input/input.part.bam

In the end, your output directory will contain several .tags.tsv files, as well as a file named tagInfo.txt. This file contains information about your sequencing run, including the total number of tags considered. This file is used by later peak-calling programs to quickly reference information about the experiment. Then we call peak by using these tag file:

findPeaks input/ip/ -style factor -o output/part.peak -i input/input/

Important parameters

输出文件为 /home/test/chip-seq/output/part.peak, 示例如下

#PeakID chr     start   end     strand  Normalized Tag Count    focus ratio     findPeaks Score Total TagControl Tags (normalized to IP Experiment)       Fold Change vs Control  p-value vs Control      Fold Change vs Local      p-value vs Local        Clonal Fold Change
chrIII-1        chrIII  78346   78578   +       69987.1 0.862   5971.000000     5966.0  201.6   29.60   0.00e+00  26.54   0.00e+00        0.50
chrIII-2        chrIII  133     365     +       61226.9 0.775   5364.000000     5227.0  116.3   44.93   0.00e+00  59.92   0.00e+00        0.61
chrI-1  chrI    141663  141895  +       41225.6 0.854   3515.000000     3514.0  169.1   20.78   0.00e+00 17.09    0.00e+00        0.50
chrII-1 chrII   165145  165377  +       35334.5 0.845   3015.000000     3018.0  171.1   17.64   0.00e+00 14.47    0.00e+00        0.50
chrII-2 chrII   555827  556059  +       34817.1 0.790   2973.000000     2970.0  159.6   18.61   0.00e+00 10.12    0.00e+00        0.50
chrIII-3        chrIII  163527  163759  +       31266.1 0.826   2662.000000     2670.0  186.0   14.35   0.00e+00  14.16   0.00e+00        0.51

3b) Motif Analysis

HOMER contains a program called findMotifsGenome.pl that can find enriched motifs in ChIP-Seq peaks

findMotifsGenome.pl output/part.peak sacCer2 output/part.motif.output -len 8

Important parameters

• Region Size (-size <#>, -size <#>,<#>, -size given, default: 200) The size of the region used for motif finding is important. If analyzing ChIP-Seq peaks from a transcription factor, Chuck would recommend 50 bp for establishing the primary motif bound by a given transcription factor and 200 bp for finding both primary and "co-enriched" motifs for a transcription factor. When looking at histone marked regions, 500-1000 bp is probably a good idea。
• Motif length (-len <#> or -len <#>,<#>,..., default 8,10,12) In general, it's best to try out enrichment with shorter lengths (i.e. less than 15) before trying longer lengths.
• Number of motifs to find (-S <#>, default 25) Specifies the number of motifs of each length to find. 25 is already quite a bit.

最重要的输出文件为 /home/test/chip-seq/output/part.motif.output/homerResults.html, 示例输出参见 这里

4) Tips/Utilities

4a) Preparation .bam from ChIP-seq data

If you want to get these two .bam files by yourself, you can follows these steps.

1. download data
   1. The fastq data for yeast ChIP-seq was downloaded from GSE61210 .
   2. Input data was downloaded from GSM1499619;
   3. IP data was downloaded from GSM1499607.

2. build yeast bowtie index

   Yeast sacCer2 genome data was downloaded from UCSC http://hgdownload.soe.ucsc.edu/goldenPath/sacCer2/bigZips/chromFa.tar.gz.

   Index was build with commad:

   tar -xvf chromfa.tar.gz
   cat *.fa >yeast.allchrom.fa
   mkdir bowtie_index_yeast
   bowtie-build yeast.allchrom.fa bowtie_index_yeast/sacCer2
   

3. mapping

   bowtie -p 4  -m 1  -v 3  --best --strata bowtie_index_yeast/sacCer2 \
       -q input/ip.fastq -S input/ip.sam
   

4. sampling

   As the .sam file is too big for tutorial example, so we selected parts of them as example file.

   samtools sort input/ip.sam >input/ip.sorted.bam
   samtools index input/ip.sorted.bam
   samtools view input/ip.sorted.bam chrI chrII chrIII -b >input/ip.part.bam
   

4b) peak file header

The top portion of the peak file will contain parameters and various analysis information

图2. peak file 元信息

Two fields need special notification:

• genome size represents the total effective number of mappable bases in the genome;
• Approximate IP effeciency describes the fraction of tags found in peaks versus genomic background. This provides an estimate of how well the ChIP worked. Certain antibodies like H3K4me3, ERa, or PU.1 will yield very high IP efficiencies (<20%), while most rand in the 1-20% range. Once this number dips below 1% it's a good sign the ChIP didn't work very well and should probably be optimized.

4c) peak calling using MACS

MACS can also be used to call peaks (the program is already installed in Docker).

macs2 callpeak -t input/ip.part.bam -c input/input.part.bam --outdir output/macs_peak \
    --name=yeast_macs_p05 --format=BAM --gsize=1.2e7 --tsize=50 --pvalue=1e-5

The main output is /home/test/chip-seq/output/macs_peak/yeast_macs_p05_peaks.xls, which is a tabular file containing information about called peaks. You can open it in excel and sort/filter using excel functions. Information include:

• chromosome name
• start position of peak
• end position of peak
• length of peak region
• absolute peak summit position
• pileup height at peak summit, -log10(pvalue) for the peak summit (e.g. pvalue =1e-10, then this value should be 10)
• fold enrichment for this peak summit against random Poisson distribution with local lambda, -log10(qvalue) at peak summit

4d) reference

HOMER call peak: http://homer.ucsd.edu/homer/ngs/peaks.html
HOMER find motif: http://homer.ucsd.edu/homer/ngs/peakMotifs.html MACS is introduced in "Identifying ChIP-seq enrichment using MACS"

5) Homework and more

1. 解释 findPeaks 和 findMotifsGenome.pl 主要参数的含义。

2. 提交文件: 提交Snf1蛋白在DNA上结合的 peak的位置, fold change, p-value等信息 和 motif 及p-value信息 (peak 和 motif最好放到一个word/pdf文件中提交）。

   • Snf1蛋白是染色体重塑复合体中的重要组成部分之一，为了探究Snf1蛋白会特异性识别哪些DNA序列，Workman小组在酵母中做了该蛋白的ChIP-seq实验。
   • 其中Snf1蛋白IP实验数据为 ip.chrom_part.bam，Input背景数据为input.chrom_part.bam。请同学们根据本节课所学内容寻找Snf1蛋白在DNA上结合的 peak 和 motif。
   • 作业的数据：bam 文件在 Docker 的 /home/test/chip-seq/homework/ 文件夹中。
   • 筛选标准 1. peak:Fold Change (vs Control) >=8，且 p-value (vs Control) < $10^{-8}$ ; 2. mofit: p-value < $10^{-10}$ ;