1. 程式人生 > 其它 >access 匯入txt 找不到可安裝的isam_Homer軟體的介紹最全面而詳細的找motif教程

access 匯入txt 找不到可安裝的isam_Homer軟體的介紹最全面而詳細的找motif教程

技術標籤:access 匯入txt 找不到可安裝的isammac上用perl xml檔案找不到

7e8f69b3-401f-eb11-8da9-e4434bdf6706.other

HOMER是一套用於Motif查詢和二代資料分析的工具。HOMER中的工具是使用Perl和C++編寫的,是Linuxcommandlinebased。HOMER這個軟體是一個大雜燴,能解決幾乎所有的高通量測序資料的分析。(這裡,我們僅僅是介紹Motif查詢這個功能)

7f8f69b3-401f-eb11-8da9-e4434bdf6706.png1.HOMER軟體安裝和配置

homer軟體的安裝

使用conda安裝,一句話搞定:conda install -c bioconda homer

使用configureHomer.pl完成Homer軟體的配置

##首先下載configureHomer.pl

wgethttp://homer.salk.edu/homer/configureHomer.pl

##使用configureHomer.pl配置Homer

perlconfigureHomer.pl-install

configureHomer.pl配置Homer的語法是:(下一節講具體應用)

perl/path-to-homer/configureHomer.pl[options]

我們看下這些options都有哪些:

-list 列出可用的包 -install 安裝homer需要用到的資料包 -version 安裝homer包時,可以指定包版本 -remove 移除包 -update

更新所有包到最新版本 -reinstall 強制重新安裝所有已經安裝過的包 -all 安裝所有包 -getFacts (add humor to HOMER - to remove delete contents of data/misc/)-check檢查第三方軟體:samtools,DESeq2,edgeR-make 重新配置和編譯可執行檔案 -sun SunOS系統,使用gmake 和 gtar代替make 和 tar -keepScript 不更新configureHomer.pl-url安裝時,使用的資源地址,預設:http://homer.ucsd.edu/homer/Hubs&BigWigsettings(withreadexistingsettingsfromconfig.txtifupgrading): -bigWigUrl
(Setting for makeBigWigs.pl) -bigWigDir (Setting for makeBigWigs.pl) -hubsUrl (Setting for makeMultiWigHub.pl) -hubsDir (Setting for makeMultiWigHub.pl)Configurationfiles:下載update.txt,更新config.txt

7f8f69b3-401f-eb11-8da9-e4434bdf6706.png2. HOMER包介紹和安裝

homer包分為4種:

  • SOFTWARE:homer工具包,包含一些常用資料。

  • ORGANISMS:物種特異性的資料,包含Geneaccession,genedescriptions,andGOanalysis資訊。大多數資料來自於NCBIGenedatabase。下載promoter或genome資料時,會自動下載對應Organism包。

  • PROMOTERS:Promoter序列資訊和Promoter富集分析的檔案。大多數資料來自RefSeqtranscript。

  • GENOMES:基因組序列及其註釋資訊。

這裡需要特別注意下:雖然conda安裝好了Homer,但並沒有包含參考序列或註釋資料。所以需要使用configureHomer.pl下載各種資料包。

下載資料之前先檢視資料列表,看Homer已經有哪些資料包:

語法公式:perl/path-to-homer/configureHomer.pl-list
具體事例:perl/home/jhuang/miniconda2/share/homer-4.9.1-6/configureHomer.pl-list
等同於:perl/home/jhuang/miniconda2/share/homer-4.9.1-6/.//configureHomer.pl-list##主要是給系統一個configureHomer.pl路徑

如上述介紹,Homer資料包安裝使用-install引數

perl/path-to-homer/configureHomer.pl-installmouse#下載老鼠啟動子資料
perl/path-to-homer/configureHomer.pl-installmm10#下載mm10小鼠參考基因組
perl/path-to-homer/configureHomer.pl-installhg19#下載hg19人的參考基因組
具體例項:
nohupperl/home/jhuang/miniconda2/share/homer-4.9.1-6/.//configureHomer.pl-installhg19&#太慢了所以掛後臺了

7f8f69b3-401f-eb11-8da9-e4434bdf6706.png3.利用ChIP-Seq結果在基因組區域中尋找富集的Motifs

首先對bed檔案進行處理

用MACS2軟體進行peak calling,也就是找比對結果(bam檔案)的peaks,MACS2找到的peaks會存放在生成的*.bed檔案裡。homer軟體找motif需要讀取我們這裡得到的bed格式的peaks檔案。homer軟體不僅可以找到motif,還可以註釋peaks。

homer軟體在讀取bed檔案時,需要提取對應的列作為輸入檔案。我們要對MACS找到的peaks記錄檔案,還需提取對應的列給HOMER作為輸入檔案,提取操作為:

awk'{print$4" "$1" "$2" "$3" +"}'sample_peaks.bed>sample_homer.bed

不理解沒關係,我們拿我跑的chipseq資料做示範,看awk對bed檔案做了什麼。
使用命令less -S ChIP-Seq_H3K4Me3_1_trimmed_rmdup_summits.bed看下awk處理之前bed檔案長什麼樣子:

$less-SChIP-Seq_H3K4Me3_1_trimmed_rmdup_summits.bed
chr12920229203ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_1260.18610
chr1714432714433ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_2165.19276
chr1762743762744ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_3195.10117
chr1839761839762ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_436.37314
chr1859839859840ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_5122.88363
chr1876411876412ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_620.77421
chr1877742877743ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_777.55285
chr1894490894491ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_8256.37842
chr1896069896070ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_9183.38295

使用awk對bed檔案進行處理,處理完的檔案命名為ChIP-Seq_H3K4Me3_1_homer.bed

awk'{print$4" "$1" "$2" "$3" +"}'ChIP-Seq_H3K4Me3_1_trimmed_rmdup_summits.bed>ChIP-Seq_H3K4Me3_1_homer.bed
$less-SChIP-Seq_H3K4Me3_1_homer.bed
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_1chr12920229203+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_2chr1714432714433+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_3chr1762743762744+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_4chr1839761839762+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_5chr1859839859840+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_6chr1876411876412+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_7chr1877742877743+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_8chr1894490894491+
ChIP-Seq_H3K4Me3_1_trimmed_rmdup_peak_9chr1896069896070+

把bed檔案變成這個樣子,用於Homer軟體的讀取,我們稱這樣的檔案叫做 Homer Peak/Positions file
我們可以看到Homer Peak/Positions file有4列:

  • 第一列:染色體

  • 第二列:起始位置

  • 第三列:終止位置

  • 第四列:鏈的方向(+/-or0/1,where0="+",1="-")

尋找富集Motifs

findMotifsGenome.pl命令用於在基因組區域中尋找富集Motifs
命令使用格式:

findMotifsGenome.pl<HomerPeak/Positionsfile><genome><outputdirectory>-size#[options]

使用例項:尋找我們剛剛生成的Homer Peak/Positions file檔案:ChIP-Seq_H3K4Me3_1_homer.bedmotif使用命令:

findMotifsGenome.plChIP-Seq_H3K4Me3_1_homer.bedhg19ChIP-Seq_H3K4Me3_1_motifDir-len8,10,12
#引數解釋-輸入檔案:awk處理好的HomerPeak/Positionsfile-參考基因組:這裡是hg19-輸出檔案:給一個路徑和輸出檔案的名字-len:motif大小設定,預設8,10,12;越大需要的計算資源越多

上述命令(找motif)每一樣品需要執行30-40分鐘後,得到資料夾ChIP-Seq_H3K4Me3_1_motifDir,資料夾裡的內容有:

-rw-rw-r--1jhuangjhuang30KApr1121:56homerMotifs.all.motifs-rw-rw-r--1jhuangjhuang9.9KApr1121:19homerMotifs.motifs10-rw-rw-r--1jhuangjhuang12KApr1121:56homerMotifs.motifs12-rw-rw-r--1jhuangjhuang8.6KApr1121:10homerMotifs.motifs8drwxrwxr-x2jhuangjhuang12KApr1121:59homerResults-rw-rw-r--1jhuangjhuang144KApr1121:59homerResults.htmldrwxrwxr-x2jhuangjhuang4.0KApr1121:08knownResults-rw-rw-r--1jhuangjhuang226KApr1121:08knownResults.html-rw-rw-r--1jhuangjhuang45KApr1121:08knownResults.txt-rw-rw-r--1jhuangjhuang81Apr1120:56motifFindingParameters.txt-rw-rw-r--1jhuangjhuang1.9KApr1120:57seq.autonorm.tsv

其中會生成一個詳細版的網頁結果,生成這些檔案的具體作用下節講述

858f69b3-401f-eb11-8da9-e4434bdf6706.other

常用引數:
-bg:自定義背景序列;
-size:用於motif尋找得片段大小,預設200bp;-sizegiven設定片段大小為目標序列長度;越大需要得計算資源越多;
-len:motif大小設定,預設8,10,12;越大需要得計算資源越多;
-S:結果輸出多少motifs,預設25;
-mis:motif錯配鹼基數,預設2bp;
-norevopp:不進行反義鏈搜尋motif;
-nomotif:關閉重投預測motif;
-rna:輸出RNAmotif,使用RNAmotif資料庫;
-h:使用超幾何檢驗代替二項式分佈;
-N:用於motif尋找得背景序列數目,default=max(50k,2xinput);耗記憶體引數

使用 annotatePeaks.pl 對peaks進行註釋

命令使用格式:

annotatePeaks.pl<HomerPeak/Positionsfile><genome>1>output.peakAnn.xls2>output.annLog.txt

使用例項:註釋ChIP-Seq_H3K4Me3_1_homer.bedpeaks使用命令:

annotatePeaks.plChIP-Seq_H3K4Me3_1_homer.bedhg191>ChIP-Seq_H3K4Me3_1_peakAnn.xls2>$ChIP-Seq_H3K4Me3_1_annLog.txt

7f8f69b3-401f-eb11-8da9-e4434bdf6706.png4.HOMER Motif 分析基本步驟和結果分析

Homer主要被用於ChIP-Seq分析,但是核酸序列motif尋找問題都可以嘗試使用。

注意:這裡大部分是在介紹Homer軟體背後是怎麼工作的,不是說我們在使用Homer軟體找motif時的步驟,我們使用的時候直接呼叫Homer軟體的內部命令就好

4.1 預處理

提取序列
提供的資料是基因組位置資訊,就需要提取對應的DNA資訊;提供基因號時,需要選擇啟動子區域。對應著就是我們前面用awk處理bed檔案,最後得到要求的那四列。背景選擇
未指定背景序列時,HOMER會自動選擇。(上面chipseq處理的時候就沒有指定背景序列)
對基因組某些區域進行分析時,從基因組隨機選擇GC含量一致的序列作為背景序列。
對啟動子進行分析時,除用來分析外的所有啟動子將被作為背景。
自定義背景使用引數-bgGC 標準化
目標序列(對應著上面的就是Homer Peak/Positions file)和背景序列會基於GC含量按5%作為bin檢視GC含量的分佈。背景序列會得到權值,從而使得其GC含量分佈與目標序列一致。
ChIP-Seq實驗得到序列GC含量。自動標準化
需要分析的序列除了GC含量會帶來誤差,其他的生物學現象,外顯子中密碼子偏好性或測序實驗中偏好性都會影響分析。對於足夠強的偏差,HOMER會自動追蹤目標序列和背景中顯著差異的特徵序列,並通過調整背景序列的權重來平衡輸入資料和背景中短寡聚核酸序列不平衡。
短寡聚核酸序列長度可以通過引數-nlen指定。

4.2 重頭預測Motifs

預設情況下,HOMER呼叫homer2進行motif分析;通過引數"-homer1"可以指定老版本工具。

將輸入序列解析為寡聚核苷酸序列
將輸入序列按照motif長度期望值解析為寡聚核苷酸序列,以及建立Oligo資料表。Oligo資料表中記錄著每條oligo在目標序列和背景中被發現的次數。Oligo 自動標準化全域性搜尋階段
Oligo表格資訊構建好之後,HOMER對富集的Oligo進行全域性搜尋。如果一個Motif是富集的,那麼屬於這個Motif的Oligo也應該會富集。首先,HOMER會搜尋可能富集的Oligo。HOMER允許錯配。
使用引數-mis 調節允許的錯配數目Mask and Repeat
當最優oligo被優化成motif後,motif對應的序列從要分析的資料中移除,接下來再分析最優的…..直到25(預設值,"-S")個motifs被發現。
比如,我們這裡處理chipseq時的情況

898f69b3-401f-eb11-8da9-e4434bdf6706.other

計算已知Motifs是否富集

3.5.1 匯入Motif庫

為了搜尋輸入資料中已知Motifs ,HOMER 可以輸入已知Motifs 資料,可以使用HOMER 預設的 ("data/knownTFs/known.motifs"),也可以是自己構建("-mknown")。

3.5.2 篩選每一個Motif

對於每個motif,HOMER計算丰度(包含motif的序列/backgroundsequences),ZOOPS(zerooroneoccurencepersequence)計數以及使用超幾何檢驗或二項式計算顯著性。

4.3 Motif 的結果分析

我們chipseq,找motif的輸出檔案如下圖:

8b8f69b3-401f-eb11-8da9-e4434bdf6706.other

這些輸出結果包括:

4.3.1 Motif Files

*.motif

檔案包含motifs的資訊

"*.motif"檔案格式:

>ATATTGCG1-ATATTGCG10.055970-947.0213230T:434.0(1.88%),B:22.3(0.08%),P:1
0.9970.0010.0010.001
0.0010.0010.0010.997
0.9970.0010.0010.001
0.0010.3630.0010.635
0.0010.0010.0010.997
0.0010.0010.9970.001
0.0010.9970.0010.001
0.0010.0010.9970.001
>AAAATCGG2-AAAATCGG10.633913-833.8916830T:458.0(1.99%),B:34.6(0.13%),P:1
0.9970.0010.0010.001
0.9970.0010.0010.001
0.9970.0010.0010.001
0.9970.0010.0010.001
0.3510.0010.0010.647
0.0010.9970.0010.001
0.0010.0010.9970.001
0.0010.0010.9970.001

一個motif的資訊分為一塊。motif資訊首行是motif各種統計資訊;其他行對應各個A/C/G/T的佔比。

motif資訊首行解析:

  1. ">"+序列(可能是空白)example:>ASTTCCTCTT

  2. Motif名字example:1-ASTTCCTCTTorNFkB

  3. 檢測閾值對數值example:8.059752

  4. 富集P-value對數值example:-23791.535714

  5. 0 用於老版本格式的佔位符

  6. T:17311.0(44.36%),B:2181.5(5.80%),P:1e-10317

  7. T:#(%)-包含motif的目標資料序列數除以目標資料序列總數

  8. B:#(%)-包含motif的背景序列數除以背景序列總數

  9. P:#-富集p-value

  10. 用逗號分隔Motif統計量,example:Tpos:100.7,Tstd:32.6,Bpos:100.1,Bstd:64.6,StrandBias:0.0,Multiplicity:1.13

  11. Tpos:motif在目標序列中的平均位置(0=startofsequences)

  12. Tstd:motif在目標序列中位置的標準偏差

  13. Bpos:motif在背景序列中的平均位置(0=startofsequences)

  14. Bstd:motif在背景序列中位置的標準偏差

  15. StrandBias:logratioof+strandoccurrencesto-strandoccurrences.

  16. Multiplicity:具有一個或多個結合位點的序列中每個序列的平均出現次數

4.3.2重頭預測(de novo)的 motif
首先會對motif進行去冗餘,將每個motif的概率矩陣轉換為向量,求motif之間的Pearson相關性。HTML結果:

908f69b3-401f-eb11-8da9-e4434bdf6706.other

表格中,BestMatch/Details項中:
MoreInformation:與預測的motif相似的的已知motifs
SimilarMotifsFound:與預測的motif相似的的其它預測motifs

4.3.3 已知 motif 的富集情況(knownResults.html)

918f69b3-401f-eb11-8da9-e4434bdf6706.other

參考資料:
1.http://homer.ucsd.edu/homer/ngs/peakMotifs.html
2.https://www.jianshu.com/p/1384173c353b
3.https://www.jianshu.com/p/4af579b2e532
4.https://www.jianshu.com/p/2bf020cc8f90
5.https://www.jianshu.com/p/1602c2621c2b
6.https://www.jianshu.com/p/1602c2621c2b

由學徒獨立完成,O(∩_∩)O謝謝