一篇不知道拖了多久的筆記……欠下的債早晚還是要還的……
0.軟件下載及安裝
官網(wǎng)：
參考教程：
其實之前已經(jīng)有兩位前輩寫過教程了的，我這只是再給需要的朋友舉個例子打個樣。詳情請戳下面兩個鏈接
下載軟件：
解壓：
unzip MCScanX.zip
cd到MCScanX后，按照
介紹的方法，給
msa.h
dissect_multiple_alignment.h
detect_collinear_tandem_arrays.h
這三個文件前面添加上#include <unistd.h>
即，用vim將文件打開后在開頭加上上述的參數(shù)段。這個錯誤的原因是，MCScanX 不支持64位系統(tǒng)。
然后make一下就安裝好了
1.文件準備
從TAIR上下載TAIR10的pep和gff3數(shù)據(jù)，要先用notepad++的正則替換功能（vim的批量替換如果會用的話也是可以的）將蛋白質(zhì)的基因ID后面的注釋信息給刪掉，否則會影響后續(xù)的工作。接下來處理gff3文件。
先看一眼原始的數(shù)據(jù)長啥樣：
我們需要的是第1、4、5列和第9列的ID部分
用awk、grep和vim的全局替換功能可以得到所需要的簡化版gff3文件：
cat TAIR10.gff3 |awk '/gene/ {print $1 "\t" $3 "\t" $4 "\t" $5 "\t" $9 }'> tmp.gff3
$3是用來標記用的，我們只取含有g(shù)ene的列
但是會有一些問題，比如：
在搜索gene的時候只要包含gene的詞就會被檢索到，于是有以下三個混進來的：
pseudogene
transposable_element_gene
mRNA_TE_gene
于是我就：
grep -v 'mRNA_TE_gene' tmp2.gff3 > tmp3.gff3
grep -v 'pseudogene' tmp3.gff3 > tmp4.gff3
grep -v 'transposable_element_gene' tmp4.gff3 > tmp5.gff3
由于當時學(xué)藝不精比較菜就只能先暫時這樣用一下了……
現(xiàn)在想想好像在上一步awk的時候可以：
cat TAIR10.gff3 |awk '/^gene/ {print $2 "\t" $0}' > tmp.gff3
加上一個“^”表示以gene開頭，應(yīng)該就可以免去上面那繁瑣的幾步了,最后生成的時候再倒一下回來就好了。
可以用
grep -v "^gene" tmp3.gff3
來檢測是否還包含有非gene的行。
經(jīng)過一系列的awk操作之后獲得最終的結(jié)果：
另一個物種如法炮制即可。
updates:也可以這么寫，三步搞定：
awk '$3 == "gene" {print $1 "\t" $3 "\t" $4 "\t" $5 "\t" $9}' TAIR10.gff3 > new.gff3
用vim全局替換處理ID列的冗余信息
用awk將信息位置調(diào)換正確
2.建庫
makeblastdb -in AT.pep -dbtype prot -parse_seqids -out ATdb
因為我需要做的是maca對擬南芥的共線性比較，所以只需要對擬南芥建庫即可。（如果需要物種內(nèi)也進行共線性比較的話這里需要把兩個物種的蛋白質(zhì)文件cat到一起來建庫）
3.比對
nohup blastp -query maca.pep -db ATdb -out AT_maca.blast -evalue 1e-10 -num_threads 30 -outfmt 6 -num_alignments 5 &
參數(shù)稍微解釋一下。
-query 被比對的物種的蛋白質(zhì)文件
-db 上一步建的庫
-out 輸出結(jié)果
-evalue 設(shè)置輸出結(jié)果的e-value值，值越小相似度越高。
-num_threads 比對所用的線程數(shù)
-outfmt 輸出的格式
-num_alignments 輸出比對上的序列的最大值條目數(shù)
這一步還是要點時間的。記得掛nohup。
4.運行MCScanX
將之前的兩個物種的gff3文件cat成一個文件后命名為AT_maca.gff
將上一步生成的AT_maca.blast和AT_maca.gff放到新建的AT_maca文件夾內(nèi)
MCScanX ./AT_maca
之后就會生成一個叫AT_maca.collinearity的文件和一個AT_maca.html文件夾。
AT_maca.html里面是這樣滴：
5.下游分析
java circle_plotter -g AT_maca.gff -s AT_maca.collinearity -c circle.ctl -o AT_maca_circle.png
這一步需要修改修改circle.ctl文件
只要把你希望畫出來的序列的ID給寫進去就好了
可以畫dot圖，bar圖和circle圖，個人覺得circle圖最好看了~
大概就是這樣了。如果有什么錯誤的地方歡迎討論指出喲。
