R:http://www.bioconductor.org/packages/release/bioc/html/exomePeak.html
安装:# 看来作者更新的挺好
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("exomePeak", version = "3.8")
安装有问题……
用bioclite安装好了:
source("http://bioconductor.org/biocLite.R")
biocLite("exomePeak")
使用:
library('exomePeak')
GENE_ANNO_GTF="gencode.vm19.GRCm38.all.ano.gtf"
f1="R18060088-QMJ-Library-4-OHT-IP1_combined_R.bam"
f2="R18060088-QMJ-Library-4-OHT-IP2_combined_R.bam"
IP_BAM=c(f1,f2)
f3="OHT1.sorted.bam"
f4="OHT2.sorted.bam"
INPUT_BAM=c(f3,f4)
f5="R18060088-QMJ-Library-DMSO-IP1_combined_R.bam"
f6="R18060088-QMJ-Library-DMSO-IP2_combined_R.bam"
TREATED_IP_BAM=c(f5,f6)
f7="DMSO1.sorted.bam"
f8="DMSO2.sorted.bam"
TREATED_INPUT_BAM=c(f7,f8)
result = exomepeak(GENE_ANNO_GTF=GENE_ANNO_GTF, IP_BAM=IP_BAM, INPUT_BAM=INPUT_BAM,TREATED_IP_BAM=TREATED_IP_BAM, TREATED_INPUT_BAM=TREATED_INPUT_BAM)
UCSC_TABLE_NAME = "?" #如下:
gtf给的例子:
如果只需要找peaks:
result_OHT = exomepeak(GENOME = "mm10",UCSC_TABLE_NAME = "ensGene", IP_BAM=c(f1,f2), INPUT_BAM=c(f3,f4),EXPERIMENT_NAME = "exomePeak_OHT")
result_DMSO = exomepeak(GENOME = "mm10",UCSC_TABLE_NAME = "ensGene", IP_BAM=c(f5,f6), INPUT_BAM=c(f7,f8),EXPERIMENT_NAME = "exomePeak_DMSO")
结果:
# 画分布图:
library(m6Amonster)
peak <- read.table("dmso.peak.bed",sep = "\t", stringsAsFactors = F)
plotMetaGene(peak,gtf = "mm10.gtf")