0.簡介 

生信分析中經常要根據指定條件查詢相似序列，比如構建多個樣品間的非冗餘基因集、分析樣品間的相似程度等等，cd-hit這款軟體就可以用較短的時間解決此類問題，可以對單個數據集進行去冗餘，包括DNA/RNA序列和蛋白序列，也可以對兩個資料集進行比較。其工作原理可概述為：將所有序列按照引數設定進行聚類，並將每一組聚類中的最長序列作為代表序列進行輸出，同時給出每組聚類下的每個序列名可供相似度分析使用。下面我們來簡單介紹一下它的使用方法。

1.下載與安裝

網址：http://cd-hit.org ；http://www.bioinformatics.org/cd-hit/ ；https://github.com/weizhongli/cdhit/archive/V4.6.2.tar.gz；

這是一個在linux系統下使用的工作，我們可以給自己的電腦裝一個雙系統或者在windows下使用linux的虛擬機器。然後我們可以執行下面的命令進行解壓（注意我們要將路徑先切換到安裝包所在的檔案目錄下，或者在執行命令時使用完整路徑）。

gzip -d cdhit-4.2.tar.gz

然後進入到解壓後的資料夾(我解壓後的資料夾為cdhit-4.2，同樣要注意我們的檔案路徑問題，如果上面使用的是完整路徑，最好這裡也使用完整路徑，比如我使用完整路徑是‘cd /home/zpf/cdhit-4.2’)

cd cdhit-4.2

最後編譯一下就可以了，執行make

make

然後我們就可以使用這個工具了。

2.輸入檔案格式

Cd-hit的輸入檔案僅有一個fasta格式檔案 ，一般來說cd-hit是將幾個樣品的基因或蛋白序列進行聚類，所以需要將這些樣品的序列彙總到一起作為輸入檔案，可在linux系統下通過cat命令實現：

cat a.fasta b.fasta c.fasta > all.fasta

 其中a.fasta，b.fasta，c.fasta為fasta格式的三個樣品基因或蛋白序列，all.fasta為彙總後的序列，在分析中作為cd-hit的輸入序列。值得注意的是，在三個樣品序列中不能有序列名相同的序列，否則會出現錯誤。因此，一般在分析時會在各樣品序列名前新增樣品名，這樣即可避免重複。序列名是fasta檔案中以“>”開頭的行空格之前的內容，如下圖中藍色線圈出部分。

3.輸出檔案

Cd-hit有兩個輸出檔案：一個是隻含有所有代表序列（即去冗餘後的序列）的fasta檔案，其格式參看圖1；另一個是以.clstr結尾的聚類資訊檔案，其格式如圖2。

 以“>”開頭的是一個聚類組。每組下面按序號排列，如上圖中Cluster 1組有5個聚類序列。每個聚類序列有一個百分比或“*”，百分比代表該序列與代表序列的相似度，“*”代表該序列即為代表序列。

4.去除冗餘的基本思路

首先對所有序列按照其長度進行排序，然後從最長的序列開始，形成第一個序列類，然後依次對序列進行處理，如果新的序列與已有的序列類的代表序列的相似性在cutoff以上則把該序列加到該序列類中，否則形成新的序列類。之所以快主要是兩個方面的原因：一個是使用了word過濾方法，即如果兩條序列之間的相似性在80%（假設序列長度為100），那麼它們至少有60個相同的長度為2的word，至少有40個相同的長度為3的word，至少有20個相同的長度為4的word。基於這個原則，在處理新的序列的時候，如果新的序列與已有序列的相同word的長度不能滿足這些要求則不需要進行比對了，這極大的降低了時間消耗；另外一個速度快的原因是使用了index table，可以很快的計算序列之間相同word的數目。

   #當序列相似性在80%時，有20個位點是有差異的，極端的情況就是這20個位點對應的長度為2的字串都不一樣，因此是40個不一樣，當有更多的不一樣時，兩條序列的相似性不可能在80%；同理，如果這20個位點對應的長度為4的字串都不一樣，則有80個不一樣。

5.使用方法和引數介紹

Cd-hit執行時用很多引數可以進行調整設定，其執行命令為（引數僅為示例）在剛才編譯的檔案路徑下執行：

cd-hit -i all.fasta -o new.fa -c 0.8 -aS 0.8 -d 0

下面簡單介紹一下重要的幾個引數：

-i：輸入檔案，fasta格式。

-o：輸出檔案字首，輸出檔案有兩個，分別為fasta格式序列檔案和以.clstr結尾的聚類資訊檔案。

-c：較短序列比對到長序列的bp與自身bp數的比值超過該數值則聚類為一組，預設為0.9。

-d：聚類資訊檔案中各個聚類組中序列名的長度，設為0則將取完整序列名。

-aL：控制代表序列比對嚴格程度的引數，預設為0，若設為0.8則表示比對區間要佔到代表（長）序列的80%。

-AL：控制代表序列比對嚴格程度的引數，預設為99999999，若設為40則表示代表序列的非比對區間要短於40bp。

-aS：控制短序列比對嚴格程度的引數，預設為0，若設為0.8則表示比對區間要佔到短序列的80%。

-AS：控制短序列比對嚴格程度的引數，預設為99999999，若設為40則表示短序列的非比對區間要短於40bp。

下圖詳解了-aL，-AL，-aS，-AS四個引數。

aL = Ra / R

AL = R - Ra

aS = Sa / S

AS = S - Sa

6.cdhit的缺點

1 它不能保證同一個序列類中的序列的相似性都在threshold之上，因為每次比對都是用新序列與序列類的代表序列進行，這就有可能使得序列類中除了代表序列外其他序列之間的相似性在threshold之下。比如A是代表序列，B與A的相似性大於0.95，C與A的相似性也大於0.95，但是這並不能保證B與C的相似性也大於0.95.

2 它不能保證一個序列類的病毒與另外一個序列類中的病毒的相似性也在threshold之上，原因還是在於用代表序列代表了整個序列類。

3 基於word filter的方法使得使用每個長度的word能夠處理的冗餘性水平有限，如使用長度為2的word只能夠得到相似性在50%以上的序列，長度為3的word只能夠得到相似性在66.7%以上的序列類，類似的，長度為5的word只能夠得到相似性在80%以上的序列。在實際應用的時候需要注意選擇的word長度與threshold的匹配。