Chip-seq 流程

sra-tools 下載安裝

#加壓包
tar zxf sra-tools.tar.gz

#將來(lái)需要調(diào)用的命令都在bin文件夾中，更改環(huán)境變量(全局)，>>的意思是追加不是覆蓋欢顷，以下路徑需要自己主動(dòng)變更捉蚤，因?yàn)槭俏易约旱南到y(tǒng)布持，你只要找到你存放sratoolkit的路徑替換掉/home......ubuntu64/之間的內(nèi)容即可题暖。
export PATH=$PATH:/home/decen/software/sratoolkit.2.9.2-ubuntu64/bin >> ~/.bashrc
source ~/.bashrc
#調(diào)用一下命令試試看胧卤，prefetch SRR6819004,立即關(guān)閉終端，文件太大
#但是我嘗試幾次好像都是臨時(shí)生效惜纸，關(guān)閉終端還是要重新配置變量
export PATH=$PATH:/share/home/chenli-lyo/soft/sratoolkit.2.9.6-1-ubuntu64/bin >> ~/.bashrc
export PATH="/share/home/chenli-lyo/soft/sratoolkit.2.9.6-1-ubuntu64/bin:$PATH"
export PATH="share/home/chenli-lyo/soft/FastQC/fastqc$PATH"
#方法二：第二個(gè)方法更有效
vi ~/.bashrc  #用vi/vim編輯器修改bashrc文件
i  #開(kāi)始插入內(nèi)容堪簿，移動(dòng)光標(biāo)到最底部哪审，不管前面有任何內(nèi)容，都要重啟一行錄入以下虑瀑，并且你要找到你存放sratoolkit的路徑替換掉/home......ubuntu64/之間的內(nèi)容即可滴须。
export PATH=$PATH:/home/username/local/app/sratoolkit/bin
ESC, :wq  #退出vi編輯器并保存文件
source ~/.bashrc  #讓配置生效
#關(guān)閉linux扔水，重啟動(dòng)
--------------------- 

fastqc 下載安裝

#在下列網(wǎng)站上下載對(duì)應(yīng)的fastqc壓縮包：http://www.bioinformatics.babraham.ac.uk/projects/download.html#fastqc
#解壓縮：
unzip fastqc.zip
#解壓縮后生成了一個(gè)名為：fastqc的文件
cd FastQC #進(jìn)入fastqc文件魔市，可以看到里面有一個(gè)fastqc執(zhí)行文件，
chmod 755 fastqc(將fastqc設(shè)置為可執(zhí)行程序)
 ./fastqc #即可運(yùn)行fastqc程序待德。如想要在任何目錄下運(yùn)行fastqc程序枫夺，則需將fastqc程序的路徑添加至環(huán)境變量.bashrc中即可较坛。
echo 'export PATH=~/soft/FastQC:$PATH'  >>~/.bashrc
source ~/.bashrc
#檢測(cè)是否安裝成功
fastqc -h 
#開(kāi)始質(zhì)控
fastqc SRR4034951.fastq.gz 

lsf集群提交作業(yè)

#使用bsub提交命令丑勤，例如提交shell 腳本
bsub bash -i xxx.sh

fastqc質(zhì)控結(jié)果解讀

Per tile sequence quality

Position specific failures of flowcells

介紹

當(dāng)Per tile sequence quality顯示fail或者warning确封，表明測(cè)序的lane或某個(gè)run中出現(xiàn)出現(xiàn)了部分故障爪喘，從而影響一些特定的區(qū)域和循環(huán)，進(jìn)而使測(cè)序數(shù)據(jù)的質(zhì)量下降稠诲。另外侦鹏，如果read的3’端的質(zhì)量是好的略水，就意味著存在瞬時(shí)質(zhì)量損失(Transient quality loss)的區(qū)域難以被剪切處理。

質(zhì)控過(guò)濾

#過(guò)濾
java -jar /share/home/chenli-lyo/soft/Trimmomatic-0.39/trimmomatic-0.39.jar SE -phred33 SRR3085650.fastq.gz TrimmedSRR3085650.fastq.gz ILLUMINACLIP:TruSeq3-SE:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:50

samtools的安裝

參考http://www.reibang.com/p/6b7a442d293f

#解壓samtools
tar -xjf *.tar.bz2
#進(jìn)入samtools文件夾
cd samtools-1.9
#編譯渊涝，安裝
make
make prefix=/share/home/chenli-lyo/soft/samtools-1.9 install
#加入環(huán)境變量
echo 'export PATH=/share/home/chenli-lyo/soft/samtools-1.9/bin:$PATH' >>~/.bashrc
source ~/.bashrc

比對(duì)

• HISAT2/STAR/Tophat: suitable for transcriptome based or RNA-seq alignment (splice-aware aligner)

• Bowtie/Bowtie2/bwa: suitable for genome based alignment (ChIP-seq/WGS)

  使用Bowtie2 進(jìn)行比對(duì)

  #在bowtie2的文件夾里跨释，輸入：
  gmake NO_TBB=1 
  
  echo 'PATH=/share/home/chenli-lyo/soft/bowtie2-2.3.5.1:$PATH' >>~/.bashrc
  source ~/.bashrc
  #bowtie2安裝完成
  #需要注意的是：
  #這條命令把bowtie2 生成的sam文件通過(guò)管道|傳遞到samtools胸私，將sam轉(zhuǎn)換為bam文件岁疼，省去中間sam文件的空間占用
  #genome_index 指的是用于bowtie2的索引文件（如下圖），而不是參考基因組本身缆娃，構(gòu)建過(guò)程參考后文疙驾。
  #genome_index 需要指定路徑及其共用文件名函荣，比如我的索引文件放在/data/ref/bowtie2/mm10目錄下傻挂，但是需要輸入的參數(shù)為/data/ref/bowtie2/mm10/mm10。最后一個(gè)mm10指的是共用文件名金拒。
  #參數(shù)說(shuō)明
  #-q: 輸入文件為fastq
  #--phred33: 測(cè)序堿基的質(zhì)量體系兽肤，現(xiàn)在基本都是33
  #-p: 線程數(shù)
  #--no-unal：不保留未必對(duì)上的記錄
  #-x：索引前綴
  #-S：sam格式輸出
  bowtie2 -p 16 -3 5 --local -x /share/home/chenli-lyo/soft/mm9/mm9 -U SRR4034952.fastq.gz | samtools sort -O bam -o /share/home/chenli-lyo/ChipEpop/SRR4034952.bam
  #參考https://blog.csdn.net/u011262253/article/details/79833969
  

Call Peaks

用MACS2 call peak

需要Python 2.7

#安裝Python 2.7
#參考 https://blog.csdn.net/qq_23113053/article/details/61203557

#安裝numpy
python2 pip install numpyXXX
#numpy下載要注意版本號(hào)x86
#安裝之前需前置numpy幢码，numpy要前置Django
#解壓Django python2.Django-1.11之前的
python2.7 -m pip install --user Django-1.11.22-py2.py3-none-any.whl
 #同樣安裝pillow
python2.7 -m pip install --user Pillow-6.1.0-cp27-cp27mu-manylinux1_x86_64.whl
#安裝numpy
python2 -m pip install --user numpy-1.16.4-cp27-cp27mu-manylinux1_x86_64.whl

#MACS2安裝
#參考http://www.mamicode.com/info-detail-1658533.html
cd 進(jìn)入macs2 文件夾
#沒(méi)有sudo權(quán)限，所以用prefix安裝辕坝，參考https://blog.csdn.net/yuan_lo/article/details/48289317
##參考此博客，解決所有問(wèn)題https://blog.csdn.net/Zephyr_Hu/article/details/81836347


cd soft/MACS2-2.1.2/
export PYTHONPATH=$PYTHONPATH:/share/home/chenli-lyo/soft/softbin
export PYTHONPATH=$PYTHONPATH:/share/home/chenli-lyo/soft/softbin/lib64/python2.7/site-packages/
python2.7 setup.py install --prefix=/share/home/chenli-lyo/soft/softbin
echo 'PATH=/shareshare/home/chenli-lyo/soft/softbin/masc2:$PATH' >>~/.bashrc
echo 'PATH=/shareshare/home/chenli-lyo/soft/softbin/bin/masc2:$PATH' >>~/.bashrc
source ~/.bashrc

#需離線安裝python-devel
#下載python-devel-2.7.5-76.el7.x86_64.rpm
#解壓python-devel-2.7.5-76.el7.x86_64.rpm
rpm2cpio python-devel-2.7.5-76.el7.x86_64.rpm | cpio -idvm
vim  ~/.bashrc
export PATH=$PATH:$share/home/chenli-lyo/soft/usr/bin/
：wq
source ~/.bashrc
echo 'PATH=/share/home/chenli-lyo/soft/bowtie2-2.3.5.1:$PATH' >>~/.bashrc
source ~/.bashrc

callpeak用法

# 常規(guī)的peak calling
macs2 callpeak -t ChIP.bam -c Control.bam -f BAM -g hs -n test -B -q 0.01
# 較寬的peak calling
macs2 callpeak -t ChIP.bam -c Control.bam --broad -g hs --broad-cutoff 0.1

--- 這次用的代碼

macs2 callpeak -t SRR4034951.bam -c SRR3085650.bam -f BAM -g hs -n test -B -q 0.01

macs2 callpeak -t /share/home/chenli-lyo/ChipEpop/SRR4034951.bam? -c SRR3085650.bam? BAM -g hs -n test -B -q 0.01
macs2 callpeak -t /share/home/chenli-lyo/ChipEpop/SRR4034951.bam? -c SRR3085650.bam? -g hs -n test -B -q 0.01 #BAM 系統(tǒng)無(wú)法識(shí)別，故刪除，不知其影響

注釋 Rstudio中用ChIPseeker

下載ChIPseeker顾画，配置工作環(huán)境

#下載ChIPseeker
source ("https://bioconductor.org/biocLite.R")
biocLite("ChIPseeker")
# 下載人的基因和lincRNA的TxDb對(duì)象
biocLite("GenomicFeatures")
biocLite("GenomeInfoDb")
biocLite("GenomicRanges")
biocLite("org.Mm.eg.db")
biocLite("TxDb.Mmusculus.UCSC.mm9.knownGene")
biocLite("clusterProfiler")
biocLite("ReactomePA")
biocLite("DOSE")
#loading packages
library("ChIPseeker")
library("GenomicFeatures")
library("GenomeInfoDb")
library("GenomicRanges")
library("org.Mm.eg.db")
library("TxDb.Mmusculus.UCSC.mm9.knownGene")
txdb <- TxDb.Mmusculus.UCSC.mm9.knownGene
library("clusterProfiler")

#讀取文件
usp7 <- readPeakFile("./usp7callpeak/test_summits.bed")