美文网首页CHIP分析
Homer

Homer

作者: shine_9457 | 来源:发表于2021-01-29 00:58 被阅读0次

    使用 conda 安装 Homer

    默认安装 homer 最新版
    conda install homer -y # -y 表示直接安装，不询问确认与否

    查看 HOMER 的数据库文件

    which homer # 返回 homer 所在路径，我这里是~/miniconda3/envs/chipseq/bin/homer cd打不开

    Linux 查找软件安装路径的五种方法
    • find -name homer # 最推荐
    • locate homer # 相当于find -name ,速度更快。不搜索具体目录，而是搜索一个数据库（/var/lib/locatedb）
    • whereis homer#只能用于程序名的搜索，而且只搜索二进制文件（参数-b）、man说明文件（参数-m）和源代码文件（参数-s）
    • which homer #在PATH变量指定的路径中，搜索某个系统命令的位置，并且返回第一个搜索结果
    • type -p homer# 相当于which homer
      参考：HOMER安装和使用攻略，如何获取Motif? 
    我这里**在主目录**使用了  `find -name homer`
结果:
./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer
./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/cpp/backup/homer
./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/cpp/homer
./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/bin/homer
./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/bin/homer
./miniconda3/envs/chipseq/bin/homer
./miniconda3/envs/chipseq/share/homer
./miniconda3/envs/chipseq/share/homer/bin/homer
./miniconda3/envs/chipseq/share/homer/cpp/backup/homer
./miniconda3/envs/chipseq/share/homer/cpp/homer

    1 cd./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer
    2 perl configureHomer.pl -install hg38/mm10 # 下载参考数据库
    3 cd ./data/genomes # 如果下载成功，可以看到有一个 hg38/mm10 的文件夹
    4 du -sh hg38/mm10 # 查看文件夹大小，为5.7G/。
    zhu:下载很慢
    我是使用conda (conda 使用基础命令安装的homer，但是我不知道版本于是 查一下怎么看软件版本
    1 conda --version #查看conda的版本
    我现在的版本：==》conda 4.9.2
    2 conda list homer #查看软件版本

    Name                    Version                   Build  Channel
homer                     4.11        https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/bioconda

    HOMER isn't part of bioconda (partially due to its unorthodox code base/organization, although there is a HOMER package that may or may not work), but regardless, bioconda can help install the required packages and other NGS software even if HOMER itself must be installed manually.
    homer官网说明

    关于homer

    • HOMER 基于 C++ 和 Perl 语言
    • 用于 motif 查找和 数据分析的工具 需要两个序列作为参数：
    • （1）参考序列：hg19、mm10 等基因组序列、promoter 序列、自定义的 FASTA 序列
    • （2）所要分析的序列：DNA 或 RNA 序列
      HOMER 适用于在大规模数据中寻找 DNA 或 RNA 序列的 motif
    关于Motif以及实战

    motif是比较有特征的短序列，会多次出现的，一般认为它的生物学意义重大，做完CHIP-seq分析之后，一般都会寻找motif 。查找有两种，一种是de novo的，要求的输入文件的fasta序列，一般是根据peak的区域的坐标提取好序列 。另一种是依赖于数据库的搜寻匹配，会将现有的ChIP-seq数据进行整合，提供更全面，更准确的motif数据库。
    两个工具：
    （1）Homer
    （2）R包bioconductor
    (3) MEME 或

    参考文献A survey of motif finding Web tools for detecting binding site motifs in ChIP-Seq data

    这里主要介绍使用HOMER寻找motif

    1 提供包含基因组坐标的文件 比如 peak文件或BED文件
    （1）分析peak文件中富集的motif：
    findMotifsGenome.pl <peak/BED file> <genome> <output directory> -size

    简单例子：
    findMotifsGenome.pl ERpeaks.txt hg18 ER_MotifOutput/ -size 200 -mask
    -mask 使用repeated-mask序列
    -size 设置motif长度

    homerMotifs.motifs <#> : these are the output files from the de novo motif finding, separated by motif length, and represent separate runs of the algorithm.
homerMotifs.all.motifs : Simply the concatenated file composed of all the homerMotifs.motifs<#> files.
motifFindingParameters.txt :         记录执行findMotifsGenome.pl的命令
knownResults.txt :                   记录motif的统计数据，text file(open in EXCEL).
seq.autonorm.tsv :                    autonormalization statistics for lower-order oligo normalization.
homerResults.html :                   de novo motif finding的格式化输出

    实践

    1 下载数据库，并锁定位置
    find -name mm10

    ./miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/data/genomes/mm10

    2 先找到自己处理好的bed文件位置

    /home/yangjy/project/epi/align/mergeBam

    2.* 发现在之前处理过后有两套_summits.bed文件，决定批量移动一下，只处理一种
    使用mv命令：
    2.    1 查看当前目录

    /home/yangjy/project/epi/align/mergeBam/

    2.*2 查看要移动到的目录

    /home/yangjy/project/epi/peaks

    先cd到文件当前所在目录
    使用mv *.log_summits.bed /home/yangjy/project/epi/peaks

    3.0 所以处理好的bed文件所在目录为：

    /home/yangjy/project/epi/peaks

    3 来到处理目录下
    cd ~/project/epi/motif
    3.1 写批量脚本，先echo一下

    for id in /home/yangjy/project/epi/peaks/*.bed;
do
echo $id
done

    成功的话：出现：

    /home/yangjy/project/epi/peaks/Control.log_summits.bed
/home/yangjy/project/epi/peaks/H2Aub1.log_summits.bed
/home/yangjy/project/epi/peaks/H3K36me3.log_summits.bed
/home/yangjy/project/epi/peaks/mm10.tss.bed
/home/yangjy/project/epi/peaks/Ring1B.log_summits.bed
/home/yangjy/project/epi/peaks/RNAPII_8WG16.log_summits.bed
/home/yangjy/project/epi/peaks/RNAPII_S2P.log_summits.bed
/home/yangjy/project/epi/peaks/RNAPII_S5P.log_summits.bed
/home/yangjy/project/epi/peaks/RNAPII_S7P.log_summits.bed

    好的，下一步：
    加上找文件名的代码

    echo $id
file=$(basename $id )
sample=${file%%.*} 
echo $sample

    于是：==》
    在

    for id in /home/yangjy/project/epi/peaks/*.bed;
do
echo $id
file=$(basename $id )
sample=${file%%.*} 
echo $sample 
done

    有：

    /home/yangjy/project/epi/peaks/Control.log_summits.bed
Control
/home/yangjy/project/epi/peaks/H2Aub1.log_summits.bed
H2Aub1
/home/yangjy/project/epi/peaks/H3K36me3.log_summits.bed
H3K36me3
/home/yangjy/project/epi/peaks/mm10.tss.bed
mm10
/home/yangjy/project/epi/peaks/Ring1B.log_summits.bed
Ring1B
/home/yangjy/project/epi/peaks/RNAPII_8WG16.log_summits.bed
RNAPII_8WG16
/home/yangjy/project/epi/peaks/RNAPII_S2P.log_summits.bed
RNAPII_S2P
/home/yangjy/project/epi/peaks/RNAPII_S5P.log_summits.bed
RNAPII_S5P
/home/yangjy/project/epi/peaks/RNAPII_S7P.log_summits.bed
RNAPII_S7P

    继续 对bed文件进行处理
    首先我们要取重要的信息，第4,1,2,3列的信息，使用awk,print命令：
    awk '{print $4"\t"$1"\t"$2"\t"$3"\t+"}' $id >homer_peaks.tmp
    使用less命令打开任意一个文档看一下
    less /home/yangjy/project/epi/peaks/H2Aub1.log_summits.bed
    看到：

    chr1    4492667 4492668 H2Aub1.log_peak_1       4.02077
chr1    4493157 4493158 H2Aub1.log_peak_2       8.94022
chr1    4493467 4493468 H2Aub1.log_peak_3       4.02077
chr1    4496553 4496554 H2Aub1.log_peak_4       2.72721
chr1    4497102 4497103 H2Aub1.log_peak_5       4.02077
chr1    4497409 4497410 H2Aub1.log_peak_6       4.02077
chr1    5916554 5916555 H2Aub1.log_peak_7       5.07861
chr1    5917032 5917033 H2Aub1.log_peak_8       5.96592
chr1    9546061 9546062 H2Aub1.log_peak_9       1.92244

    第4,1,2,3列的信息应该分别表示,eg：H2Aub1.log_peak_1 chr1 4492667 4492668
    验证一下：

    for id in /home/yangjy/project/epi/peaks/*_summits.bed;
do
echo $id
file=$(basename $id )
sample=${file%%.*} 
echo $sample  
awk '{print $4"\t"$1"\t"$2"\t"$3"\t+"}' $id >homer_peaks.tmp
done

    查看命令：head homer_peaks.tmp
    可以看到Homer Peak/Positions file有4列：
    第一列: 染色体
    第二列: 起始位置
    第三列: 终止位置
    第四列: 链的方向(+/- or 0/1, where 0="+", 1="-")

    RNAPII_S7P.log_peak_1   chr1    4085494 4085495 +
RNAPII_S7P.log_peak_2   chr1    4559682 4559683 +
RNAPII_S7P.log_peak_3   chr1    4559912 4559913 +
RNAPII_S7P.log_peak_4   chr1    4560134 4560135 +
RNAPII_S7P.log_peak_5   chr1    4560750 4560751 +
RNAPII_S7P.log_peak_6   chr1    4769884 4769885 +
RNAPII_S7P.log_peak_7   chr1    4770731 4770732 +
RNAPII_S7P.log_peak_8   chr1    4771139 4771140 +
RNAPII_S7P.log_peak_9   chr1    4771620 4771621 +
RNAPII_S7P.log_peak_10  chr1    4772002 4772003 +

    为什么要这样，因为Homer软件的说明书要求这样的格式；才能给做注解；

    接下来开始做注释：

    为什么给了mm10就可以自动做注释；
    下载数据之前先查看数据列表，看Homer已经有哪些数据包：

    语法公式：perl /path-to-homer/configureHomer.pl -list
    具体事例：perl /home/yangjy/miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/configureHomer.pl -list
    等同于：perl /home/yangjy/miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/.//configureHomer.pl -list
    主要是给系统一个configureHomer.pl路径
    save to :/home/yangjy/miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer//update.txt
    启用安装：

    同前面描述，不介意再来一遍：（Homer数据包安装使用-install参数）

    perl /path-to-homer/configureHomer.pl -install mouse #下载老鼠启动子数据
    perl /path-to-homer/configureHomer.pl -install mm10 #下载 mm10小鼠参考基因组
    perl /path-to-homer/configureHomer.pl -install hg19 #下载 hg19人的参考基因组

    实例：perl /home/yangjy/miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/configureHomer.pl -install mm10
    很慢，可以挂载后台：
    nohupperl /home/yangjy/miniconda3/pkgs/homer-4.11-pl526hc9558a2_3/share/homer/configureHomer.pl -install mm10&
    综上：

    for id in /home/yangjy/project/epi/peaks/*_summits.bed;
do
echo $id
file=$(basename $id )
sample=${file%%.*} 
echo $sample  
awk '{print $4"\t"$1"\t"$2"\t"$3"\t+"}' $id >homer_peaks.tmp  
findMotifsGenome.pl homer_peaks.tmp mm10 ${sample}_motifDir -len 8,10,12
annotatePeaks.pl    homer_peaks.tmp mm10  1>${sample}.peakAnn.xls 2>${sample}.annLog.txt 
done

    ==》输出文件
    其实control的bed为0 没有必要做注释；
    下载看看文件如何：pwd---/home/yangjy/project/epi/motif
    cd motifDir;会看到以下文件

    homerMotifs.all.motifs
    homerMotifs.motifs10
    homerMotifs.motifs12
    homerMotifs.motifs8
    homerResults
    homerResults.html
    knownResults
    knownResults.html
    knownResults.txt
    motifFindingParameters.txt
    seq.autonorm.tsv

    打开后缀为.html的文件，可以看到


    注释.png
    单独分析 寻找富集Motifs

    findMotifsGenome.pl命令用于在基因组区域中寻找富集Motifs
    命令使用格式为：
    findMotifsGenome.pl <Homer Peak/Positions file> <genome> <output directory> -size # [options]

    findMotifsGenome.pl homer_peaks.tmp mm10 ${sample}_motifDir -len 8,10,12
annotatePeaks.pl    homer_peaks.tmp mm10  1>${sample}.peakAnn.xls 2>${sample}.annLog.txt
    #参数解释
    • -输入文件：awk处理好的Homer Peak/Positions file
    • -参考基因组：这里是mm10
    • -输出文件：给一个路径和输出文件的名字
    • -len：motif大小设置，默认8,10,12；越大需要的计算资源越多

    参考：

    补充

    通过meme去来找motif，需要bed格式的peaks的坐标来获取fasta序列。
    MEME，链接：http://meme-suite.org/

    相关文章

      网友评论

        本文标题：Homer

        本文链接：https://www.haomeiwen.com/subject/foyftltx.html

        延伸阅读

        深度阅读

        • 您也可以注册成为美文阅读网的作者，发表您的原创作品、分享您的心情！

        栏目导航

        热点阅读

        CHIP分析

        关于我们|服务条款|联系我们|Homer|投稿指南|网站地图|RSS订阅|排版工具|手机版

        提供经典美文摘抄,优美散文欣赏,现代诗歌精选,短篇小说,心情随笔,表白情书范文,故事会在线阅读欣赏

        Copyright © 2014-2023 Haomeiwen.com All Rights Reserved. 好美文阅读网 版权所有

        备案信息：桂公网安备 45052102000051号 · 桂ICP备13007215号-3

        本站所收录作品、热点评论等信息部分来源互联网，目的只是为了系统归纳学习和传递资讯

        所有作品版权归原创作者所有，与本站立场无关，如不慎侵犯了你的权益，请联系我们告知，我们将做删除处理！