寫在前面
快一年了，寫了重測序三兄弟插件，幾乎所有用戶都可以用 TBtools 輕松完成 BSAseq 數(shù)據分析。當然，我自己基本也沒怎么折騰過，畢竟當時只是為了上個課。前幾天將家人送回老家后，開始進入渾渾噩噩的狀態(tài)。大體是每天也不知道干啥或者不干啥，盡管確實也在干活。索性，就還是簡單整理一下一些插件功能，也再測試測試看看。于是決定，干脆先重現(xiàn)一下一些論文分析結果。比如這篇 Nature Plant 論文，定位了番茄分支數(shù)目的決定位點論文（分支越多，果實越多嘛，很好理解）。
于是下載了數(shù)據
SRR8307487 suppressed_A
SRR8307489 branched_A
于是我下載了基因組和測序數(shù)據
看起來不錯，數(shù)據挺大的?？梢蚤_干，整體步驟如下：
- 讀段回帖
- 比對結果排序
- 標記PCR重復
1.讀段回帖
首先，需要測序得到的reads比對到基因組上，使用 TBtools 的 「BWA-mem2 插件」即可。具體界面如下：
過了大概 24 小時，第一個樣品比對結束了，并產生了一個 10Gb 的 BAM 文件。這里似乎速度比較慢。有點懷疑數(shù)據是否是有接頭。回頭我可以跑一個 FastQC 看看。當然，邏輯上沒啥問題，畢竟一共是 32Gb+ 的reads。
整體時間差不多。也要考量本地電腦 6 個線程確實慢。手上我的CPU也是 2016年的。
2. 比對結果排序
這一步很快，大體20分鐘不到。
3. 標記PCR重復
重測序數(shù)據分析，最大的問題就是PCR重復會影響SNP的檢測，所以需要標記出來
這一步也很快
4. 檢測基因組變異 Call SNPs
開始檢測變異
點擊開始，然后報錯
調整后，重新 Start，于是開始跑
我們用的是 bcftools，速度還是很快的
5. 過濾低質量 SNP
幾乎所有變異檢測的軟件或流程出來的結果都是需要過濾的。因為每個人認為可靠的 SNP 的標準是不同的。比如有些人覺得測序深度要10X，有些人覺得5X就足夠了。在TBtools的這個Pipeline中，我們直接用默認的。邏輯上對BSAseq影響不大。
這一步邏輯上非?？臁Ｊ聦嵣?，整體流程瓶頸還是在比對。等待幾分鐘即可。
6. 進行QTL位點檢測
使用相對簡單，不過還是有兩三個參數(shù)要手動給一下，具體在設置輸入的 VCF 文件時，會有預覽窗口，通過預覽文本窗口的信息，就可以直接復制黏貼設置，整體如下
不會很久，不過一般也要幾分鐘
復現(xiàn)結果如下
與原文結果比較
注意：圖片美化上自然還有很多辦法....；所以結果挺好
具體輸出的表格信息也可以看一下定位的區(qū)間
尷尬了，區(qū)間沒對應上。當然，第一反應是他沒有問題，同時咱們的流程也沒問題。
隨后查看了下 2019年05月，過去三年多了，基因組版本更新了，現(xiàn)在我們用的是新版本?？纯丛瓉淼淖⑨屧谖覀冇玫陌姹镜淖⑨寧缀危?br>
結果沒問題，偏離一點點。當然了。精準的結果我們完全可以通過指定某個染色體，重新跑一次，結合 deltaSNP 和 G’value 來看。此外，我們其實只用了論文中的一部分數(shù)據，每個混池都只選了10+個材料，完全可以再下載兩組數(shù)據，合并好然后跑一下。另外還注意到，原文作者對群體做了分層，
OK，復現(xiàn)結果。我相信，只要合并數(shù)據，完全就可以解決這個問題，感興趣的朋友可以試試
SRR8307487 suppressed_A
SRR8307488 suppressed_B
# 下面兩組
SRR8307489 branched_A
SRR8307490 branched_B
如果你復現(xiàn)了，有結果了，歡迎投稿，有稿酬
寫在最后
至于為什么要來復現(xiàn)這個圖稿，自然我仍然希望 TBtools的「插件模式」可以給更多朋友打開數(shù)據分析大門，自己動手，豐衣足食。邏輯上現(xiàn)在只要你有群體，其實可能2000塊錢不到，就可以搞定兩個F2或者F1群體的混池測序，為何不試試呢？也沒多少錢。
最后，大伙用系列插件時，請引用Wrapper對應的原軟件，上述包括但不僅限于：
bwa-mem2, samtools, bcftools, qtlseqr
至于 TBtools 最近我們觀測的最高引用頻率是 每天在Web of Science 被引 10 次，所以是否引用，大伙可以看心情。
