生信小工具:Orthofinder使用教程

2021-03-02 生信菜鳥團

最近Orthofinder2開始陸續更新，文章已經投發到biorxiv上，在網上搜了一圈，關於Orthofinder的使用的中文教程幾乎是空白的，這周就藉此機會和大家一起來學習一下這款軟體。

Orthofinder介紹

OrthoFinder是比較基因組學中的實用的，運行快速，準確的全面的工具。它的主要功能是，找到了正交群和直系同源物，推斷出所有正交群的根基因樹，並識別那些基因樹中的所有基因重複事件。它還為所分析的物種推斷出有根的物種樹，並將基因重複事件從基因樹比對到物種樹的分支中。

另外，OrthoFinder還為比較基因組分析提供全面的統計數據。相比OrthoMCL， OrthoFinder使用簡單安裝都很簡單，安裝只需要用過conda就行，運行只需一組FASTA格式的蛋白質序列文件（每個物種一個）。

基礎知識介紹

為了讓大家更好的理解使用Orthofinder，在正式介紹其使用的流程前，先來回顧一下與之相關的生物學知識。

直系同源物和旁系同源物

Orthologs(直系同源物)是在兩個物種的最後共同祖先（LCA）中來自單個基因的一對基因。直系同源物是同源性基因，是物種形成事件的結果。Paralogs（旁系同源物)是同源基因，是重複事件的結果。下圖就可以看到，不同物種間的alpha-chain gene互為Orthologs(直系同源物)。這時候可以引用一個新名詞orthogroup （正交群）就用來形容自一組物種的LCA中的單個基因的基因組（在圖中就是alpha chain gene）。然後同一物種間alpha 和beta chain gene互為Paralogs(旁系同源物)。最後所有這些關係都可以由OrthoFinder來識別。

為什麼要研究正交群？

正交群中的所有基因都來自單個祖先基因。因此，正交群中的所有基因都有類似的序列和功能。由於基因重複和丟失在進化中經常發生，一對一的直系同源物很少見，通過分析orhtogroup所有直系同源的情況（一對一，多對一，多對多），我們可以分析數據的所有情況。

直系同源物是由系統發育定義的。它不能通過胺基酸含量，密碼子偏好，GC含量或其他序列相似性測量定義。在沒有系統發育的情況下使用這些分數來定義，直向同源物的方法只能提供猜測。從序列相似性提供粗略的類比猜測直系同源，類似於從氣味中猜測顏色。因此，真正識別直向同源物的唯一方法是通過分析系統發育樹。確保直系分配正確的唯一方法是對所考慮的物種的最後共同祖先的單個基因下降的所有基因進行系統發育重建。這組基因是正交群。因此，定義直系同源的唯一方法是分析正交群。

安裝

廢話說了一大篇，是時候開始安裝軟體了。在這裡推薦大家使用conda來安裝，簡單明了，不用擔心其他Dependencies的安裝

conda install orthofinder

接著準備好測試數據，為了展示的方便，我將原始測試數據重新命名：

##重命名

cd /work1/ExampleData

mv Mycoplasma_agalactiae.faa Maga.faa

mv Mycoplasma_gallisepticum.faa Mgal.faa

mv Mycoplasma_genitalium.faa Mgen.faa

mv Mycoplasma_hyopneumoniae.faa Mhyo.faa

##創建新的數據目錄

mkdir Mycoplasma/

cp *faa Mycoplasma/

正式運行Orthofinder，相當簡單的操作, -f輸入目錄，裡面包含你需要運行的蛋白質fasta文件， -t 所用到的CPU數目。基本的用法就如下，更多的可以去manual中查看。

orthofinder -f Mycoplasma/ -t 2

生成的結果會存儲於Results_XXX文件中，現在簡單看看裡面有啥。正交群以兩種不同的格式返回，一種是制表符分隔的表格格式（* .csv），另一種格式與orthoMCL（* .txt）的輸出格式相同。還應該有一個名為WorkingDirectory的目錄，其中包含運算過程的中間文件，例如blast結果。

cd Mycoplasma/Results_Nov10/

ll -thrl

-rw-r 1 hhu hhu  12K Nov 10 16:50 Orthogroups.GeneCount.csv

-rw-r 1 hhu hhu  59K Nov 10 16:50 Orthogroups.csv

-rw-r 1 hhu hhu  89K Nov 10 16:50 Orthogroups.txt

-rw-r 1 hhu hhu  110 Nov 10 16:50 Orthogroups_SpeciesOverlaps.csv

-rw-r 1 hhu hhu  33K Nov 10 16:50 Orthogroups_UnassignedGenes.csv

drwxr-x--- 6 hhu hhu 4.0K Nov 10 16:50 Orthologues_Nov10

-rw-r 1 hhu hhu 2.5K Nov 10 16:50 SingleCopyOrthogroups.txt

-rw-r 1 hhu hhu 1.5K Nov 10 16:50 Statistics_Overall.csv

-rw-r 1 hhu hhu 2.5K Nov 10 16:50 Statistics_PerSpecies.csv

drwxr-x--- 2 hhu hhu 4.0K Nov 10 16:50 WorkingDirectory

使用R簡單分析Orthofinder的結果

安裝加載以下R包

install.packages("micropan")

install.packages("phangorn")

library(micropan)

library(phangorn)

讀取Orthofinder中正交群的文件。結果文件將正交群存儲為制表符分隔表，其中每一行是正交群，其名為OG，後跟7位數字（例如OG0000123）。表中的每列對應一個物種，並且與fasta文件具有相同的名稱。表中的單元格包含序列ID。由於正交組可以包含來自單個物種的幾個基因（即，旁系同源物），因此表中的單個cell可以包含由逗號分隔的幾個序列ID。下面是讀取結果文件的例子：

# Define orthofinder results directory (may vary)

resultsDir <- "orthofinder/Mycoplasma/Results_Oct22"

grpsTbl <- read.table( file.path(resultsDir,"OrthologousGroups.csv"),

                      header=T, sep = "\t", row.names = 1, stringsAsFactors=F)

# Split the comma-separated paralogs

grps <- strsplit( as.matrix(grpsTbl),", ")

# grps is now a list of character vectors..

# Add dimensions to grps:

dim(grps) <- dim(grpsTbl)

# Add row and column names:

dimnames(grps) <- dimnames(grpsTbl)

獲取特定正交群的信息：

grps["OG0000018", ]```

可以清晰的看到，在對應的fa文件中，在該特定的正交群OG000018下所包含的fasta序列，其結果如下：

## $Maga.faa

## character(0)

##

## $Mgal.faa

## [1] "Mgal_AAP56913.2" "Mgal_AAP56907.2" "Mgal_AAP56908.2"

##

## $Mgen.faa

## [1] "Mgen_AAC71531.2" "Mgen_AAC71529.1" "Mgen_AAC71563.2"

##

## $Mhyo.faa

## character(0)

接著讓我們來看開正交群的數目，為了獲得特定正交群的序列，需要在每個正交群中構建每個物種的基因數量的矩陣：

grpSize <- apply(grps,2,sapply,length)

我們可以在正交群OG0000018中檢查每個物種中的旁系同源物的數量：

grpSize["OG0000018", ]

## Maga.faa Mgal.faa Mgen.faa Mhyo.faa

##        0        3        3        0

我們還可以繪製每個正交群中存在的基因數量：

plot(table(apply(grpSize>0,1,sum)),ylab="no. orthogroups",xlab="no. genomes represented")

另外一寫有趣的事我們可以查看，每個物種中只有一個直系同源的正交群（1：1：1：1正交群）：

sum(apply(grpSize==1,1,all)) # number of 1:1:1:1 groups

## [1] 254

正好只存在於兩個物種（1：1：0：0正交群）的正交群：

sum(apply(grpSize==1,1,sum)==2 & apply(grpSize > 1,1,sum)==0)

## [1] 155

讓我們找出這個數字背後所對應的物種：

# get the 1:1:0:0 groups

is1100 <- apply(grpSize==1,1,sum)==2 & apply(grpSize > 1,1,sum)==0

spcs <- sub("\\.faa","",colnames(grps)) # get species names

# Count the combination of species-pairs in the 1:1:0:0 groups

table(apply(grpSize[is1100, ], 1, function(grpSizeRow){

 paste(spcs[ grpSizeRow==1 ], collapse="&")

}))

##

## Maga&Mgal Maga&Mgen Maga&Mhyo Mgal&Mgen Mgal&Mhyo Mgen&Mhyo

##        15         2        57        70         8         3

基本的用法和分析就講解到這裡，當然這只是最基本的東西。其結果涉及到的下遊分析還有很多很多，例如構建直系同源基因樹等等的知識，大家可以自行進行探討。

星期三是我的專題日哦，如果你喜歡我的文章，請點一下文末的小贊給與我更多創作的動力，如果你想了解更多，歡迎點擊閱讀原文，查看更多我的筆記~

▼ 如果你生信基本技能已經入門，需要提高自己，請關注下面的生信技能樹，看我們是如何完善生信技能，成為一個生信全棧工程師。

▼ 如果你是初學者，請關注下面的生信菜鳥團，了解生信基礎名詞，概念，紮實的打好基礎，爭取早日入門。

相關焦點

  • 從零開始學生信-orthofinder的安裝
    orthofinder的使用很簡單，可以先行查看help信息。【注意】：orthofinder由python2編寫，請在安裝或使用前將python3注釋掉。【下載】：conda install orthofinder # 用conda 2 進行安裝，一次安裝不成功者請重試（附帶軟體很多），提示更新conda時，請更新。
  • 使用OrthoFinder進行直系同源基因分析
    談論到直系同源基因分析的時候，大部分教程都是介紹OrthoMCL，這是2003年發表的一個工具，目前的引用次數已經達到了3000多，但這個軟體似乎在
  • 一個很實用分析orthofinder結果的工具
    但是Orthofinder不能對每個group進行功能分析和數據統計分析，這裡筆者開發兩個工具，可以進行簡單的功能和數據分析。code 下載連結：https://gitlab.bioinformatics.unibe.ch/troder/orthofinder_tools/-/tree/masterhttps://github.com/davidemms/OrthoFinder/issues/4511.統計相同基因： orthogroup_to_gene_name.py
  • Phylogenomics-DensiTree繪製詳細教程
    數據：某屬物種的全基因組數據，基於胺基酸序列，先使用Orthofinder進行了同源基因鑑定，然後提取單拷貝直系同源基因，使用MAFFT進行多序列比對，Trimal提取保守區，過濾掉低於50aa序列，然後使用IQtree進行單基因樹構建，最後使用DensiTree軟體進行DensiTree繪製。# 1.
  • 生信小工具:bedtools的使用(2)
    不用擔心，bedtools早就想到這一點了，bedtools complement就是為此而生的。例如，讓我們找到基因組的所有非外顯子（即內含子或基因間）區域。注意，要做到這一點，你需要一個「基因組」文件，告訴bedtools文件中每個染色體的長度。
  • 生信小工具之:Bedtk
    這次推文繼續和大家分享近期一些實用的生信小工具。今天給大家帶來的是Heng Li大神，最新寫的一款小工具Bedtk。
  • 100個在線生信小工具
    ，生信分析基本上已經成為了一個繞不開的過程，面對高通量測序的大量數據，我們可能需要在Linux系統中使用專門的生信分析工具完成，這些工具通常學習成本較高，對於大部分研究人員來說很難快速的獨立完成。但是在研究過程中，我們除了需要對大數據進行分析之外，面對少量的數據進行一些統計、識別、比對、功能分析、可視化也是需求量很大的工作內容，這是一些功能單一、操作簡單的在線分析工具就成為了更優的選擇。本篇文章將介紹100個在線的生信分析小工具，涵蓋了幾乎所有可能遇到的分析內容。
  • 生信技能樹論壇-生信基礎版塊介紹-計算機基礎
    本周我們將為大家帶來論壇-生信基礎版塊的介紹。讓我們使用論壇的搜索，看看這個版塊有些什麼吧，這裡對本版塊的帖子簡單分了類，所以可以試試在本版搜索以下關鍵詞（如圖）...看看有沒有你想要的(#^.^#)生信技能樹論壇-生信基礎版塊介紹-計算機基礎與Perl相關與Shell相關與Python相關與JAVA相關與Matlab相關與R語言相關與Javascript
  • 生信軟體安裝神器
    生信草堂將會與更多的優秀微信公眾號合作，把-優秀的微信推文呈現給大家，希望可以幫助讀者更多的了解生信技術，培養和提高讀者的生信分析能力！號外，號外，號外你想和生信分析大神做好朋友麼？你想讓自己的生信分析走上快車道麼？那就趕快加入我們的生信交流微信群吧！
  • Conda | 輕鬆安裝生信工具
    最近有很多朋友問我生信工具安裝的問題，對於初學者來說，工具安裝是一個非常頭疼的問題。
  • 轉錄組生信入門教程--每周一更(1)
    從本周開始，我們公眾號開始轉載生信媛徐洲更的一組推文，推出一組轉錄組生信分析入門教程。
  • 100個在線生信小工具-收藏貼
    寫在前面在與生物相關的研究中，生信分析基本上已經成為了一個繞不開的過程，面對高通量測序的大量數據，我們可能需要在
  • 生信技巧第5課-生信人必須安裝的軟體
    https://www.rstudio.com/products/rstudio/download/Xshell https://www.netsarang.com/products/xsh_overview.htmlwinscp https://winscp.net/eng/docs/lang:chsWindows10系統可以直接啟動ubuntu，自行搜索教程如何開啟
  • 生信分析連載-差異基因篩選必備知識點
    GEO2R主要是對系列數據（series）進行分析，但不是所有系列數據都能用GEO2R工具進行分析，比如測序數據就不能使用GEO2R，對於這類數據，工具欄中的「Analyze with GEO2R」按鈕不會顯示。GEO2R可以在不編程的情況下實現簡單的差異表達分析...
  • 可重複的生信分析系列二:Conda的介紹
    Conda可以說是版本控制和生信工具安裝的一大神器。相信大家對它了解肯定不少，但是又該怎麼樣利用它，進行可重複的分析呢？今天繼續講第二部分 Conda的介紹。本節教程將會使用到docker，去安裝minconda的鏡像。
  • 生信學習第一步!免費領取一套docker虛擬機!
    生信分析環境搭建是生信學習路上的一隻攔路虎！
  • 生信小技巧:並行運行的秘密
    生信學習入門後，我們都掌握了如何在後臺運行命令。一般我們可以通過 forloop，將一系列的命令逐一運行。
  • 生信小技巧:使用linux中的screen命令去同時管理你多個項目
    在這些情況下，你都可以使用很溜的screen讓你的工作變得更加輕鬆！另外screen有一個很大的好處就是，他的運行時掛在伺服器或者電腦的後臺，所以當你在學校辦公室工作使用screen運行分析到一半，然後你關閉terminal的窗口後回家後，可以通過screen的功能繼續你的工作和分析。
  • 生信小工具專題:BBTools/BBMap Suite 的使用 (2)
    接著上一次內容繼續介紹，BBTools、BBMap中的一些實用的小工具。BBMap Read Merger合併雙末端（PE）reads,在預期中這些reads有重疊的位置（例如16S擴增子測序）。通過重疊檢測將配對的reads合併為單個read。如果具有足夠的覆蓋率，還可以通過kmer擴展合併非重疊的reads。
  • 生信分析人員如何系統入門R(2019更新版)
    五年前作為一個初出茅廬的菜鳥生信工程師苦於沒有專業交流社群，遂自建了生信菜鳥團QQ群和博客，一點一滴積累了數萬人氣，進而和若干圈內好友組建了生信技能樹聯盟，三年前的直播生物信息學編程活動細節還歷歷在目，QQ群微信群記帳錄製視頻忙的不亦樂乎，因此產生了程式語言系統入門系列教程，如下：現在回過頭來看，很多教程已然過時，當然並不是說的知識點過時，其實linux基本上幾十年都沒有怎麼變動過基礎知識的