ChIP-Seq全部分析流程（转）

http://www.bio-info-trainee.com/2773.html

 

ChIP-seq 实验和数据获得：

•   将蛋白交联到DNA上；
•   通过超声波剪切DNA链；
•   加上附上抗体的磁珠用于免疫沉淀靶蛋白；（抗体很重要）
•   接触蛋白交联；纯化DNA
•   送去测序 

ChIP-seq分析相关软件：

•   sratoolkit、fastQC、bowtie2、samtools、macs2、htseq-count、bedtools

ChIP-seq分析流程：

•  质量控制（FastQC）
•  序列比对（bowtie2/bwa）
•  peak calling（MACS2）
•  peak注释（ChIP-seeker）
• 

ChIP-seq分析案例流程：

   1、安装相关软件（sratoolkit、fastQC、bowtie2、samtools、macs2、htseq-count、bedtools）

   2、下载相关数据

•   样本数据
•   参考基因组的bowtie2索引数据和注释文件数据
•   参考基因组数据

    3、序列比对

     将得到的fastq文件用bowtie2比对小鼠参考基因组上

bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620204.fastq | samtools sort -O bam -o analysis/alignment/ring1B.bam
bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620205.fastq | samtools sort -O bam -o analysis/alignment/cbx7.bam
bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620206.fastq | samtools sort -O bam -o analysis/alignment/suz12.bam
bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620207.fastq | samtools sort -O bam -o analysis/alignment/RYBP.bam
bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620208.fastq | samtools sort -O bam -o analysis/alignment/IgGold.bam
bowtie2 -p 6 -3 5 --local -x reference/mm10 -U SRR620209.fastq | samtools sort -O bam -o analysis/alignment/IgG.bam

   结果：ls -lh

     

计算比对率 （samtools flagstat）

这里我是通过samtools flagstat 计算的比对率，结果如下：

ring1B : 60.49%
cdx7: 87.52%
suz12: 67.01%
RYBP: 84.17%
IgGold: 57.80%
IgG: 82.80%   

用IGV查看

从官网下载IGV, 解压即可使用，linux下 igv.sh打开IGV界面，Windows 下点击igv.bat

●首先载入参考基因组，可以载入自己下载好的参考基因组，也可选择IGV中含有的参考基因组，ref_Genome 必须是fasta格式 

● 载入比对的文件，比对的文件必须先经过sort 和 index, 才可加载。

samtools sort a.bam a.sort 
samtools index a.sort.bam

igv.sh

● 比对可视化结果：IgG是对照组，其他组有的峰在对照组没有，即为peak

4、 peak calling （MACS2）

macs2 callpeak -c IgGold.bam -t suz12.bam -q 0.05 -f BAM -g mm -n suz12 2suz12.macs2.log
macs2 callpeak -c IgGold.bam -t ring1B.bam -q 0.05 -f BAM -g mm -n ring1B 2/ring1B.macs2.log
macs2 callpeak -c IgGold.bam -t cbx7.bam -q 0.05 -f BAM -g mm -n cbx7 2cbx7.macs2.log
macs2 callpeak -c IgGold.bam -t RYBP.bam -q 0.01 -f BAM -g mm -n RYBP 2RYBP.macs2.log

    5、结果注释和可视化

结果的注释用的是Y叔 的 Chipseeker包。

ChIPseeker的功能分为三类: ● 注释：提取peak附近最近的基因， 注释peak所在区域 ● 比较：估计ChIP peak数据集中重叠部分的显著性；整合GEO数据集，以便于将当前结果和已知结果比较 ● 可视化： peak的覆盖情况；TSS区域结合的peak的平均表达谱和热图；基因组注释；TSS距离；peak和基因的重叠。

下载chipseeker 有关包

• download packages

source ("https://bioconductor.org/biocLite.R")
biocLite("ChIPseeker")
biocLite("org.Mm.eg.db")
biocLite("TxDb.Mmusculus.UCSC.mm10.knownGene")
biocLite("clusterProfiler")
biocLite("ReactomePA")
biocLite("DOSE")

● loading packages

library("ChIPseeker")
library("org.Mm.eg.db")
library("TxDb.Mmusculus.UCSC.mm10.knownGene")
txdb <- TxDb.Mmusculus.UCSC.mm10.knownGene
library("clusterProfiler")

读入bed文件

ring1B <- readPeakFile("F:/Chip-seq_exercise/ring1B_peaks.narrowPeak")

Chip peaks coverage plot

查看peak在全基因组的位置

covplot(ring1B,chrs=c("chr17", "chr18"))   #specific chr

ring1B