一、不求甚解系列
软件下载及配置
conda安装： conda install -c bioconda homer
使用configureHomer.pl完成HOMER软件的配置
# 下载配置文件 wget http://homer.salk.edu/homer/configureHomer.pl # 使用配置文件进行软件配置 perl configureHomer.pl -install # 下载 hg19 人的参考基因组 # 将hg19替换为mm10可下载mm10 小鼠的参考基因组 # 因为课题组人和小鼠的研究居多，这里只以这两者举例 perl configureHomer.pl -install hg19
Motif 分析及peak注释
此处以 MACS2 软件call peak得到的bed文件举例：
HOMER软件接受两种文件格式，这里不展开介绍文件的具体格式，细节可看使用详解系列
findMotifsGenome.pl H3K4Me3.bed hg19 H3K4Me3_motif -len 8,10,12 # peak文件：ChIP-Seq_H3K4Me3_1_homer.bed # 基因组版本：hg19 # 输出路径：ChIP-Seq_H3K4Me3_1_motifDir # motif长度：-len 8,10,12 # motif的软件默认长度为8，10，12 annotatePeaks.pl H3K4Me3.bed hg19 1> ChIP-Seq_H3K4Me3_1_peakAnn.xls 2>H3K4Me3_annLog.txt
结果展示
软件运行后的结果全部放在我们制定的输出文件路径下，将该目录下所有文件下载到本地，打开homerResults.html 文件，即可得到如下结果展示：
二、使用详解
HOMER 一个常用的motif分析软件。它通过比较两个序列集，并使用ZOOPS scoring和超几何分布（或者负二项分布）进行motif的富集分析。它主要用于ChIP-seq和promoter分析，但也可以用于核酸序列的motif分析问题。
HOMER软件可以进行多种类型的motif分析，如 promoter motif analysis ，基因组位置motif分析(ChIP-seq分析中的motif分析)，利用自定义的fasta文件进行motif分析，RNA序列的motif分析（分析CLIP-seq数据中的RNA binding elements）
HOMER进行motif分析时，需要两个数据集：
感兴趣的目标序列，如ChIP-seq分析中的peak文件
背景序列，如ChIP-seq分析中的物种全基因组序列
HOMER包的介绍与安装
HOMER包分为四种：
SOFTWARE： 软件包，里面包括所有的代码以及一些常见的数据文件，如motif矩阵
ORGANISMS：物种特异性数据包，包括一些指定物种的accession conversion数据，gene description数据，GO 分析数据。大部分数据都是以NCBI的gene数据库为基础的
PROMOTERS：启动子数据包，对promoter进行motif分析时所需要的相关数据和promoter序列信息。大部分数据都是基于RefSeq 转录组。
Genomes：基因组数据包，包括一些基因组序列和注释信息
HOMER软件的安装可看不求甚解部分。HOMER软件刚安装之后，是不包含任何序列信息的。我们可以使用 configureHOMER.pl 脚本下载所需数据。 我们可以使用 perl configureHOMER.pl -list 命令得到已有的HOMER包。
有些软件包的名称后面会跟有 -o, -p, -g等参数，这是HOMER为了区分同一物种不同类型数据所添加的。如 human-p。其中 -o 表示organism，-p 表示promoter，-g 表示genome。 这里给出几个软件包安装与删除的例子：
# install 表示安装 perl configureHomer.pl -install mouse # 安装小鼠的promoter数据集 perl configureHomer.pl -install hg19 # 安装人的hg19版本的基因组数据集 # remove 表示删除 perl configureHomer.pl -remove mm10 # 删除小鼠mm10版本的基因组数据
此外，需要注意的是，当我们下载基因组数据时，同时会自动下载物种包的数据（比如，下载hg19数据时，会自动下载human数据）。每次我们下载某个物种的基因组数据或者promoter数据时，它都会自动替我们检查是否下载了对应的物种数据。
motif分析
正如前面所说，HOMER可以进行多种类型的motfi分析。这里我们以最常见的ChIP-seq流程中的peaks motfi分析举例。该分析需要使用 findMotifsGeome.pl 脚本。
该脚本的输入文件支持两种类型的文件，一种是常见的bed文件格式(MACS2软件call peak得到的bed文件即可直接使用)，另一种是HOMER软件指定使用的peak文件格式。 下面详细讲一下这两种文件格式:
bed文件格式：findMotifsGeome.pl 脚本用到的信息为bed文件的前六列，分别是 chr, start, end, peak ID, score(HOMER用不到该信息，但是必须有), strand 。
peak文件格式：该文件格式需要五列信息（使用Tab分隔），分别是 peak ID , chr , start , end ,strand 。
当我们准备分析一个peak/bed文件的motif时，需要运行以下命令：
findMotifsGenome.pl <peak/BED file> -size # [options] # 比如： findMotifsGenome.pl ERpeaks.txt hg18 ER_MotifOutput -size 200 -mask
-mask : 该参数告诉motif分析程序，在得到一个可能的motif之后，在后续的motif分析中是否排除该motif的影响。有点类似于抽样调查中的无放回抽样。
-size: 指定用于motif分析的片段长度，默认为200； -size given 设置片段大小为目标序列长度。越大，motif分析所需要的计算资源越多。
下面介绍一些常用的其他参数：
-p : 设置线程数-S : 所寻找的motif数目，默认为25。25已经不算少了，如果自定义数目，推荐设置更少而不是更多。
-mis : 所允许的错配数，默认为2
-cpg : 在对背景序列和目标序列进行normalization时，使用CpG%矫正而非GC%。
norevopp : 只在peak所在的链搜索motif，不进行反义链搜索。HOMER默认会在两条链上都进行搜索。
-rna : 搜索RNA上的motif（如使用CLIP-seq数据进行分析的时候）
-bg : 指定自定义的背景序列。HOMER默认随机选取基因组序列作为背景。
-h: findMotifsGenome.pl脚本默认使用二项分布进行motif富集分析，这在背景序列多于目标序列时是十分有用的。但是有时我们使用 -bg 参数自定义背景序列时，其数目可能会小于目标序列，此时使用 -h 参数选择超几何分布会更加适合。
-len : motif的长度设置，默认8，10，12，越大越消耗计算资源。
-N: 用于motif分析时所需的序列数目，通常当我们设置 -len过大，内存不够时，会选择减小 -N参数或者 -size 参数。
三、原理
前面我们说到，HOMER软件可以进行多种类型的motif分析，如 promoter motif analysis ，基因组位置motif分析(ChIP-seq分析中的motif分析)，利用自定义的fasta文件进行motif分析，RNA序列的motif分析（分析CLIP-seq数据中的RNA binding elements）。 但是，不管我们如何调用HOMER，对于发现motif的基本步骤都是一致的，下面主要介绍一下各步骤的原理：
原理理解即可，HOMER已经将各个步骤封装起来，不用一步一步分别进行。
预处理
1. 序列提取
若给出的基因组位置信息，则提取出来的是对应的基因组序列；
若给出的是基因accession number，给需要选择适当的promoter区域
2. 背景选择
如果没有自定义背景序列信息（即： -bg 参数），HOMER将自动为用户选取。如果使用的是基因组位置，将从整个基因组序列抓取GC含量一致的序列作为背景序列。若进行的是promoter分析，则将所有的promoter作为背景。
3. GC含量矫正
将目标序列和背景序列对GC含量进行分bin（5%区间），背景序列通过调节权重得到和目标序列同样的GC含量分布。这可以使得HOMER在分析来自CpG岛时的序列时，不会简单的找到那些GC富集的motif。下图是一个ChIP-seq分析中GC含量分bin的例子
4. 自动标准化
除了GC含量之外，目标序列还会存在其他的情况影响分析结果。如外显子上的密码子偏好等生物学现象，由A富集序列优先排列引起的实验偏差。当这些偏差显著时，将会影响HOMER的motif分析。HOMER通过给背景序列分配适当的权重，减少寡核苷酸序列频率（如AA）在背景与目标的差异。
重头预测motif
默认情况下HOMER调用homer2版本，若想使用旧版本，在命令行中添加 -homer1 参数即可
5. 解析输入序列
根据设置的motif长度，将输入序列解析为寡核苷酸，并生成一个寡核苷酸频度表。该表只记录出现的寡核苷酸和其出现的次数，可以增加motif搜索时的效率，但是也会破坏寡核苷酸与其原始序列的关系。
6. 寡核苷酸自动均一化
HOMER除了在第四步中对全局序列进行了自动均一化矫正，还在寡核苷酸频度表中进行了矫正。该方法的主要思想是将那些较短的寡核苷酸与较长的寡核苷酸相平衡，但是该方法由于存在许多与motif同长的寡核苷酸，需要分配数量众多的权重。对于那些极端的序列偏好，该方法不会移除它们。
7. 全局搜索
在构建好寡核苷酸频度表后，HOMER进行全局motif搜索。其基本原理就是若某个motif富集，则其存在的寡核苷酸也同样富集。为了增加灵敏度，HOMER允许搜索motif时的错配，同时跳过多个错配的寡核苷酸节省计算资源。
8. 矩阵优化
9. Mask and Repeat
当最优oligo被优化成motif后，motif 对应的序列从要分析的数据中移除，接下来再分析最优的…..直到设定个数的motif被发现。
搜索已知的motif（homer2版本）
10. 加载moitf数据库
HOMER从以前的数据中加载一个已知的motif列表从而搜索已知的motif还可以通过在命令行中指定的motif( -mknown)或通过编辑homer包中的文件( data/knownTFs/known.motifs )来添加motifs。
11. 扫描所有的motif
HOMER通过扫描所有的序列，确定每个motif的富集度，并计算其显著性。
motif输出的结果分析
12. Motif文件
HOMER的输出有许多后缀名为 motif 的文件，其内容如下：
每个motif信息是一块，均以>开头。>所在行的信息以tab分隔。 motif首行信息解释：
> 一致性序列：如图上的>ASTTCCTCTT
Motif名称：如图上的1-ASTTCCTCTT
比值检测概率的log值：8.059752
P-value得log值：-23791.535714
占位符：如上图得0，不具有任何信息
逗号分隔得富集信息，如：T:17311.0(44.36%),B:2181.5(5.80%),P:1e-10317
T表示带有该motif的目标序列在总的目标序列（target）中得百分比
B表示带有该motif的背景序列在总的背景序列（background）中得百分比
P表示最终富集的p-value
逗号分隔的motif统计信息，如：Tpos:100.7,Tstd:32.6,Bpos:100.1,Bstd:64.6,StrandBias:0.0,Multiplicity:1.13
Tpos：motif在目标序列中的平均位置
Tstd：motif在目标序列中位置的标准差
Bpos：motif在背景序列中的平均位置
Bstd：motif在背景序列中位置的标准差
StrandBias: 链偏好性，在正义链上的motif数与反义链motif数的比值的log
Multiplicity：具有一个或多个结合位点的序列中每个序列的平均出现次数
13 denovo预测的motif
首先对产生的motif进行去冗余，然后将motif转换为向量值从而计算皮尔逊相关系数。HOMER产生的motif结果使用网页展示，该文件名称为
homerResults.html ，而目录homerResults/
存放着所有该网页依赖的文件和图片。生成的motif首先去冗余。以下是一个简单的输出示例：
输出结果按照p-value排序，最后一列是一个链接到motif文件的超链接，可以从这个文件中找到包含此motif的其他序列。在Best
Match/Details 列中，HOMER将会展示与denovo
motif最匹配的已知motif（不可全信，最匹配不一定是最优）。该列还包含有一个名为 More information
的超链接，点击后会出现以下页面：
该页面包含motif的一些基本信息，如链接到motfi文件的超链接，且可以生成motif logo的pdf格式。
接着是与已知motif的match。这部分主要看看denovo motif和已知的motif的相似性，需要检查一下是否合理，如下：
如果点击 Simiar motifs 超链接，将会展示在denovo motif搜寻过程中与该motif相似的其他motifs（这些motif的enrichment value相对较低）。通常我们会点击查看一下是否有一些motifs被不正确的分类，比如几个motif具有几个相同的碱基。具体界面如下：
14. 已知motif的输出
该输出也是以网页形式展现，页面结果根据enrichement进行排序，具体如下：
- 本文固定链接: https://maimengkong.com/kyjc/1331.html
- 转载请注明: : 萌小白 2023年1月13日 于 卖萌控的博客 发表
- 百度已收录
