增强子（Enhancer）位于结构基因附近，是一类非编码DNA顺式作用元件，在真核生物的发育过程中通过结合转录因子、辅因子以及染色质复合物作用于启动子，可以激活或增强基因的转录。简单说：增强子是能够增加启动子活性从而增加基因转录频率的DNA序列。
超级增强子 (Super enhancer, SE) 是一类具有超强转录激活特性的顺式调控元件，富集高密度的关键转录因子、辅因子和增强子表观修饰标记，与普通增强子 (Typical enhancer, TE) 相比，超级增强子区域跨度范围通常可达8~20Kb，远高于普通增强子的跨度范围，超级增强子具有强大的调控功能。
目前对增强子的鉴定，主要采用染色质免疫共沉淀技术（ChIP-seq）针对活性增强子相关联的因子或组蛋白修饰进行检测，如转录因子、转录辅激活因子(如Mediator、p300)、组蛋白修饰H3K27ac 和H3K4me1等。在分析方法上，首先对所得增强子进行连接，主要依据为在基因组范围内，单个增强子间距离如在12.5 kb 范围内，则合并为一个增强子，即缝合增强子(Stitched enhancer)。缝合增强子和其余的单个增强子按照ChIP-seq所测信号水平的强度排序，绘制获得一张曲线图，该曲线上斜率为1 的切线的切点所得的信号值为区分超强增强子和普通增强子之间的阈值，高于该值为超级增强子，其余的则称为普通增强子。我们使用ROSE软件进行增强子鉴定。
1、宽峰检测
与同项目窄峰脚本相同
macs2 callpeak -t ../mapping/NC_IP.uni.dedup.bam -c ../mapping/NC_In.uni.dedup.bam \ -n NC \ -g $effective_genome_size \ -f AUTO \ --nomodel \ --extsize $fragment_size \ --broad \ --outdir . \ 2> NC_broad_running_status
2、enhancer鉴定
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/geneMapper.sh .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_bamToGFF.py .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_bamToGFF_turbo.py .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_callSuper.R .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_geneMapper.py .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_main.py .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_main_turbo.py .
cp /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_utils.py .
ln -s /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/annotation/ .
awk '{{if($1~"chr")print $1"\\t"$4"\\tBY_MACS\\t"$2"\\t"$3"\\t.\\t.\\t.\\t"$4;else print "chr"$1"\\t"$4"\\tBY_MACS\\t"$2"\\t"$3"\\t.\\t.\\t.\\t"$4}}' \\
{workdir}/peak/{experiment}_peaks.broadPeak \\
> {workdir}/superEnhancer/{experiment}_peaks.broadPeak.gff \\
python ROSE_main.py \\
-g {rose_anno} \\
-i {workdir}/superEnhancer/{experiment}_peaks.broadPeak.gff \\
-r {workdir}/mapping/{IP}.uni.dedup.bam \\
{if_input_is_here} \\
-s 12500 \\
-t 2000 \\
-o {workdir}/superEnhancer/{experiment}/ \\
Rscript /PUBLIC/source/RNA/IPseq/script/ROSE_superEnhancer/ROSE_callSuper.R \\
{workdir}/superEnhancer/{experiment}/ \\
{workdir}/superEnhancer/{experiment}/*_MAP.txt {experiment}_chip \\
NONE
python ROSE_geneMapper.py \\
-i {workdir}/superEnhancer/{experiment}/{experiment}_chip_AllEnhancers.table.txt \\
-g {rose_anno} \\
-o {workdir}/superEnhancer/{experiment} \\
使用脚本如下：
路径：/TJPROJ6/RNA_SH/script_dir/superenhancer
python superenhancer.py --experiment_name con_1,LA_1,con_2,LA_2 --experiment_design C_109_CHIP_1:C_109_INPUT_1,LA_109_CHIP_1:LA_109_INPUT_1,C_109_CHIP_2:C_109_INPUT_2,LA_109_CHIP_2:LA_109_INPUT_2 --mapdir /TJPROJ6/RNA_SH/script_dir/superenhancer/test/mapping/ --workdir /TJPROJ6/RNA_SH/script_dir/superenhancer/test --fa /TJPROJ1/RNA/reference_data/Animal/Homo_sapiens/GRCh37_Ensemble_75/IP/Homo_sapiens.GRCh37.75.dna_sm.primary_assembly.fa --rose_anno hg38
#---- 包括以下文件 (1) *_Plot_point.png 横坐标表示enhancer的排序,纵坐标表示enhancer区域对应的信号强度 (2) *AllEnhancer.table.txt 所有enhancer区域的结果汇总表，包括typical-enhancer和super-enhancer。 (3) *_SuperEnhancers.table.txt super-enhancer的结果汇总表。 (4) *_ENHANCER_TO_GENE.txt enhancer的注释结果 #----表头信息 (1) *AllEnhancer.table.txt 和 *_SuperEnhancers.table.txt REGION_ID: enhancer的ID，由两部分组成，14_H27ac_1_peak_1_lociStitched：14表示该enhancer由14个broad peak区域组成；H27ac_1_peak_1_lociStitched表示这14个组成peak的其中1个的peakID。 CHROM : enhancer所在的染色体。 START : enhancer的起点。 STOP : enhancer的终点。 NUM_LOCI : stitch的区段个数，表示这个enhancer区域是由多少个peak区段所组成的。 CONSTITUENT_SIZE : 组成该enhancer区域的peak区段长度的加和。 *.bam : enhancer区域对应的信号值（分别是IP和Input样本中的rpm值）。 isSuper : 是否是super-enhancer，是=1，否=0。 enhancerRank : enhancer的排序（按照IP_RPM-Input_RPM的值进行排序）。 (2) *_ENHANCER_TO_GENE.txt:enhancer的注释结果 OVERLAP_GENES: 与增强子有overlap的基因 PROXIMAL_GENES： 与增强子临近的基因，指转录起始位点TSS与超增强子边界的距离在50kb以内的基因（不包括overlap的基因）。 CLOSEST_GENE： 与增强子最为临近的基因，指转录起始位点TSS与超增强子边界距离1000kb以内的基因当中基因起始位置距离超增强子中心最近的基因。
问题：
1. 增强子分析的时候，是怎么选定“区域”的呢？
答： 区域的选择，首先增强子区域都是基于peak区域来选择的。计算peak区域上的信号值。
2. 超级增强子是按照信号区分的？信号值怎么选的呢？
答：有些增强子是紧密地成簇出现的，为了精确的识别到这种紧密成簇的增强子，认为两两增强子之间的距离如果在12.5kb范围之内，则把他们合并起来，当成一个stitched enhancer。将所有增强子区域按照ChIP-seq信号水平排序，画一张曲线图（结果中已反馈），横坐标是增强子的排序，纵坐标是每个增强子对应的信号水平，在曲线上斜率为1的切点的位置所对应的信号值是区分超级增强子与普通增强子的阈值。大于这个阈值，就是超级增强子，否则是普通增强子。
3. 怎么比较两个增强子是不是相同的呢？应该比较什么值呢？如下图一样，红色和蓝色有相同和不同的增强子。
答：通过第2个问题中所回答的那样，识别到了超级增强子和普通增强子之后，每一个增强子都是有特定的基因组区间位置的。如果两个增强子的基因组位置有交集或者有一定比例的交集（比如50%的区域是overlap的），则认为这两个增强子是相同的。
4 表头REGION_ID； CONSTITUENT_SIZE是什么意思？有用吗？PC9_CHIP_30h.uni.bam是peak的信号值？
答： REGION_ID：是每一个增强子的ID。
CONSTITUENT_SIZE：如果某个增强子只有1个区域组成，那么CONSTITUENT_SIZE就是该区域的长度； 如果增强子是stitched enhancer，那么CONSTITUENT_SIZE是组成该stitched enhancer的区域的长度的加和。
PC9_CHIP_30h.uni.bam：是每个enhancer区域对应的PC9_CHIP_30h.uni.bam数据中的信号值。
5、增强子选择的时候，从TSS前后2000bp的位置开始计算，因为有启动子区干扰，我看结果文件似乎不是这样的……
答：因为有启动子区域（TSS上下游2000bp）的干扰，如果某个enhancer完全包含在启动子区域的话，那么这个enhancer是不会被用来做stitching的。 整个enhancer的识别方法是从“Super-Enhancers in the Control of Cell Identity and Disease. Cell”这篇文献中得到的，脚本链接为：https://bitbucket.org/young_computation/rose。对于ROSE_main.py这个脚本的参数“-t”的意思是排除TSS周围一定长度范围内的区域，这个参数我们给定的是2000，也就是过滤掉那些在TSS周围2000bp区域的enhancer。
6 接合12.5Kb以内的信号，算一个增强子。
答：如果两个增强子之间的距离如果在12.5kb范围之内，则合并成一个stitched enhancer。
7 超级增强子的长度在25kb~35kb，超过是不是就应该不要了，增强子的中位值不会超过5000，一般也就800或者2000吧，超级增强子一般20000吧，我按CONSTITUENT_SIZE算了一下中位值，太高了，超级增强子长度中位值20多万了，增强子也过万，麻烦您们再确定下，我是不是搞错了，先只检查PC9的结果好了 。
答：对于超级增强子的长度是多少，这个目前我们不清楚这个长度的明确反范围。 这个项目的售后应该是2016年完成的，目前已经删掉了售后路径和结果文件。我这边找到了另一个项目（NH160987）的PC9_CHIP_0h的enhancer结果，对这个结果的中的super enhancer和普通enhancer的CONSTITUENT_SIZE进行了统计，结果如下：super-enhancer长度的中位值是42.3kb。普通enhancer的长度的中位值是2.2kb。
## super-enhancer
Min. 1st Qu. Median Mean 3rd Qu. Max. 14350 33670 42350 54900 63240 479800
## 普通enhancer
Min. 1st Qu. Median Mean 3rd Qu. Max. 263 1163 2262 3975 4671 42080
8 PROXIMAL_GENES和CLOSEST_GENE分别与超增强子相距多远的距离？
答：PROXIMAL_GENES是指转录起始位点TSS与超增强子边界的距离在50kb以内的基因（不包括overlap的基因）；CLOSEST_GENE是指转录起始位点TSS与超增强子边界距离1000kb以内的基因当中基因起始位置距离超增强子中心最近的基因。
#bin name chrom strand txStart txEnd cdsStart cdsEnd exonCount exonStarts exonEnds score name2 cdsStartStat cdsEndStat exonFrames 0 NR_075077 chr1 - 67092175 67134971 67134971 67134971 10 67092175,67096251,67103237,67111576,67113613,67115351,67125751,67127165,67131141,67134929, 67093604,67096321,67103382,67111644,67113756,67115464,67125909,67127257,67131227,67134971, 0 C1orf141 unk unk -1,-1,-1,-1,-1,-1,-1,-1,-1,-1, 0 NM_001276352 chr1 - 67092175 67134971 67093579 67127240 9 67092175,67096251,67103237,67111576,67115351,67125751,67127165,67131141,67134929, 67093604,67096321,67103382,67111644,67115464,67125909,67127257,67131227,67134971, 0 C1orf141 cmpl cmpl 2,1,0,1,2,0,0,-1,-1, 0 NM_001276351 chr1 - 67092175 67134971 67093004 67127240 8 67092175,67095234,67096251,67115351,67125751,67127165,67131141,67134929, 67093604,67095421,67096321,67115464,67125909,67127257,67131227,67134971, 0 C1orf141 cmpl cmpl 0,2,1,2,0,0,-1,-1, 0 NM_000299 chr1 + 201283451 201332993 201283702 201328836 15 201283451,201293941,201313165,201316552,201317571,201318617,201319815,201320266,201321977,201323012,201324427,201324940,201325753,201328761,201330073, 201283904,201294045,201313560,201316697,201317779,201318795,201319878,201320381,201322133,201323189,201324581,201325127,201325838,201328868,201332993, 0 PKP1 cmpl cmpl 0,1,0,2,0,1,2,2,0,0,0,1,2,0,-1,
使用gff3ToGenePred处理gff/或者使用gtfToGenePred处理gtf gff3ToGenePred /public/database/annotation/GRCm38/GRCm38.gff GRCm38_refseq.ucsc /TJPROJ2/GB/PUBLIC/software/GB_TR/mRNA/miniconda3/envs/prepare_data/bin/gtfToGenePred
完成后将 ref_gene.txt 文件拷贝到/biodata/02.software/ROSE/annotation 下。且改为 基因组_refseq.ucsc 的名称
并在ROSE_main.py 里的字典加入加入的基因组的ID