scenic/scenicplus自定义cistarget数据库

SCENIC/SCENIC+分析转录因子或者增强子时，其中有一个输入文件是motif的调控潜能排序数据库，而软件只配套提供了human、mouse、fly这三个物种，其他物种或者不满意软件提供配套的数据库可以选择自己动手创建。尤其在使用SCENIC+时，软件文档也推荐根据自己的数据创建自有库。

想要手动创建数据库首先需要安装相应的软件，安装过程这里就省略了。然后，需要准备好以下三种输入文件：

1. 包含调控区域的FASTA文件
   序列的名字最后会在分析数据时被用到，所以如果是分析转录因子建议用基因的symbol名，因为转录因子结合位点位于基因启动子区域，这样可以得到直接的TF与靶基因间的调控关系。如果是增强子用染色体区间位置，因为增强子的位置可以离靶基因很远。

>ENSG00000227232.5
ATTTTATTTTATTTTTGTAGAGATGGAATCTCACAATATTGCCCAGGCTGATCTCAAACTCTTGGGCTCAAGCAGTCCTCTCGCGTCAGCCTACCAAAGTGCTAGGATTACAGGCATGAGCCACCACGCCCAGCCAGATTACAAAATACTTTGACATCTTTTCTCTACAGCCCTCAAAAGGA

2. motif矩阵文件
   内容格式类似于MEME数据库里面motif的PFM/PPM与fasta的结合，行是motif从5'->3'每个位置上的碱基，列是每个位置上四种碱基出现的频率或者概率。每个motif单独存放在一个文件里面，文件名后缀以“.cb”结尾，如bergman__Adf1.cb，所有的文件放在一个目录里面：

>bergman__Adf1
0       0       100     0
0       100     0       0
0       0       0       100
0       50      50      0
0       100     0       0
0       50      0       50
0       0       50      50
0       100     0       0
0       50      0       50
0       0       100     0
0       50      0       50
33.33333333     33.33333333     0       33.33333333

3. motif列表文件
   普通的列表文件，可以控制哪些motif最终用于生成数据库，这里可以直接motif的文件名。程序会读取每一行然后到存放motif的目录里去查找文件并读取使用：

predrem__nrMotif0.cb
predrem__nrMotif1000.cb
predrem__nrMotif1001.cb
predrem__nrMotif1002.cb
predrem__nrMotif1003.cb
predrem__nrMotif1004.cb
predrem__nrMotif1005.cb
predrem__nrMotif1006.cb
predrem__nrMotif1007.cb
predrem__nrMotif1008.cb

准备好以上三种输入，就可以构建自己的数据库了：

python create_cistarget_motif_databases.py \
       -f tss_2kb.fa \
       -M motif_cb \
       -m motif_list.txt \
       -o hg38 \
       -t 20

由于小序列易受到局部背景干扰的影响，通过扩展小序列在两侧引入背景序列，可以使得小序列的评分效果更好。所以，官方文档创建SCENIC+数据库的流程就增加了这一步，生成的调控区域的FASTA文件就包含了背景序列。不过，文档也说明了添加此填充是完全可选的，因为这不会对分析结果产生太大影响。如果引入了背景序列则构建数据库时需要-b X参数。

bash create_fasta_with_padded_bg_from_bed.sh \
        GRCh38.primary_assembly.genome.fa \
        hg38_chrlen.txt \
        consensus_regions.bed \
        hg38_1kb_bg_padding.fa \
        1000 \
        yes