利用ROSE鉴定超级增强子

1.安装ROSE

代码格式（下载并解压）：

wget https://bitbucket.org/young_computation/rose/get/1a9bb86b5464.zip
unzip 1a9bb86b5464.zip

2.数据的预处理 先把sam文件转为bam:

samtools view -bS SRR12344059.sam > SRR12344059.bam

再对bam文件进行sort （排序）:

samtools sort -@ 8 SRR12344059.bam -o SRR12344059.sorted.bam

为bam文件建立索引:

samtools index -@ 8 SRR12344059.sorted.bam

3.使用ROSE 先把call peak后的narrowpeak改成bed文件，直接在Xftp中重命名。 注意要在python 2中使用，创建一个环境（py27）即可

conda create -n py27

记得激活环境py27 

source activate py27

需要注意一定要到软件的安装目录去运行，因为会在运行目录查找annotaton这个文件夹下的物种注释文件。 还要注意R是否可以调用，如下就代表调用没问题 ######我遇到的问题：

可能是R版本有问题，解决办法：

conda update r-base -y

##########当一切准备好了就可以用ROSE鉴定超级增强子了 要在ROSE_main.py所在的位置进行以下代码的运行:

python ROSE_main.py -i /data/hj/m_embryo_H3K27ac/raw_data/trim/bidui/peak_call/narrow_peak/SRR12344059.bam_peaks.narrowPeak.bed -r /data/hj/m_embryo_H3K27ac/raw_data/trim/bidui/SRR12344059.sorted.bam -o /data/hj/m_embryo_H3K27ac/raw_data/trim/bidui/peak_call/narrow_peak/super_enhancer -g MM10 -c /home/bijf/GBM_NSC_Astrocytes/NSC_H3k27ac/fastq_data/NSC_Input_sort.bam

ROSE的具体参数：

python ROSE_main.py 
-g HG18 
-i HG18_MM1S_MED1.gff 
-r MM1S_MED1.hg18.bwt.sorted.bam 
-c MM1S_WCE.hg18.bwt.sorted.bam 
-o out_dir 
-s 12500 
-t 2500

-g指定参考基因组版本，用于检索对应的物种注释文件； -i输入gff文件，一般为使用MACS鉴定得到的Med1富集区域； -r排序后的bam文件，同时需为bam添加index； -c指定Input样本的bam文件； -o代表输出的文件； 可选参数： -s指定合并增强子的距离，默认为12.5kb, 小于该距离的两个增强子会合并为一个区间； -t指定距离TSS的距离，如果一个peak与某个转录起始位点的距离小于指定的距离，则有可能是一个启动子，这种潜在的启动子会被过滤掉； -c CONTROL_BAM，control样本的bam文件
4.生成结果

结果解读： 1.gff 文件夹 是：该文件为输入gtf文件的副本， 其中STITCHED.gff该文件为通过在STITCHING_DISTANCE将INPUT_CONSTITUENT_GFF拼接在一起创建的gff文件；文件列数如下：chrom, name, [blank], start, end, [blank], [blank], strand, [blank], [blank], name其中 name 字段的命名方式为：拼接起来的区域数+最左端区域ID。 2.mappedGFF文件夹 是：MAPPED.gff每个bam文件通过bamToGFF的输出文件，并包含以下列：（成分ID，测试区域，平均读取密度（单位为每百万位元每百万映射的单位读数密度）） 其中12KB_STITCHED_SRR9719746.sorted.bam_MAPPED.gff每个bam文件通过bamToGFF的输出文件，该文件中对增强子区域进行了拼接，包含以下列：
（拼接增强子ID，测试区域，平均读取密度（单位为百万映射每单位拼接增强子数）） 3.STITCHED_ENHANCER_REGION_MAP.txt 是：bamToGFF计算后得到的拼接增强子密度文件，包含以下列：（拼接增强子ID，染色体，拼接增强子起始位置，拼接增强子末端位置，拼接数，BAM信号等级，BAM信号） 4.Enhancers_withSuper.bed 是：可以加载到UCSC浏览器中可视化的增强子bed文件 5.Rplots.pdf 是：所有增强子散点图 6.AllEnhancers.table.txt 是：增强子列表，包含每个增强子的排名和是否为超级增强子 7.SuperEnhancers.table.txt 是：超级增强子的排名，是第8个文件的子集。包含以下列：（拼接增强剂ID，染色体，拼接增强子起始位点，拼接增强子末端，拼接数，缝合在一起的成分的大小，RANKING_BAM的信号，RANKING_BAM的等级，超增强子的二进制（1）与典型（0）） 5.用ROSE做注释enhancer 需要ROSE_geneMapper.py

Usage: ROSE_geneMapper.py [options] -g [GENOME] -i [INPUT_ENHANCER_FILE]

# 参数解释
-i INPUT    输入ROSE_mian.py生成的enhancer table文件
-g GENOME   输入genome信息（MM9,MM8,HG18,HG19等）  
-o OUT      输出路径
# 可选参数
-l GENELIST 要过滤的基因列表
-w WINDOW   搜索基因距离，默认值为50,000bp
-f format   如果使用此参数，将保持原输入文件格式输出

具体使用：

python ROSE_geneMapper.py -i /data/hj/rose/young_computation-rose-feb35cb1d955/super_enhancer/SRR9719746_AllEnhancers.table.txt 
-g MM10 -o /data/hj/rose/young_computation-rose-feb35cb1d955/super_enhancer/enhancer 

结果：  结果解读： 1.ENHANCER_TO_GENE.txt 包含：enhancer重叠基因、附近基因以及最近的基因列表 2.GENE_TO_ENHANCER.txt 包含：以每个基因为列名和其相关的增强子位置信息列表 6.在R做GO分析

#调用R包
library(clusterProfiler)
library(DOSE)
library(org.Mm.eg.db)
#enhancer GO富集分析
#首先读取文件，第一行作为列名
enhancer.relate.gene<-read.table(file="C:/Users/hj990/Desktop/SRR9719746_AllEnhancers_GENE_TO_ENHANCER.txt",header=T)
#txt里的最后一列注明是否为super enhancer，是则为1；如果不是则空白。取为1的所有行，即挑选出所有的super enhancer。
gene<-enhancer.relate.gene[which(enhancer.relate.gene$IS_SUPER=='1'),]
#查看super enhancer的数据集
head(gene)
#将第一列提出来，因为第一列是基因名
gene_name<- gene[,1]
#查看基因名
head(gene_name)
#调用R包
library("stringr")
library("clusterProfiler")
#GO的BP分析
ego_bp <- enrichGO(gene_name,OrgDb      = org.Mm.eg.db,keyType    = 'SYMBOL',ont        = "BP",pAdjustMethod = "BH", pvalueCutoff = 0.001,)
#泡泡图
dotplot(ego_bp, font.size=10)

#GO图
plotGOgraph(ego_bp)

#bar图
barplot(ego_bp)