1、準(zhǔn)備數(shù)據(jù)
mkdir mscanx
cd mscanx
mkdir data
cp 10\ Garlic.pep ~/mscanx/data/garlic_pep.fasta
cp Af.pep ~/mscanx/data/Afi_pep.fasta
cd ~/mscanx/data/
ll
將來(lái)兩個(gè)基因的蛋白序列放一起
cat Afi_pep.fasta garlic_pep.fasta >Af_As.fasta
perl -p -i -e 's/*//; s/*/X/g' Af_As.fasta 替換掉* //當(dāng)有特殊字符($, ^)要替換時(shí)，要用單引號(hào)把s///g引起來(lái)
-p 標(biāo)志要求Perl循環(huán)處理看到的所有行， -i 標(biāo)志要求Perl原地更新匹配的行內(nèi)容，而 -e 標(biāo)志則要求Perl把傳入的字符串當(dāng)作一個(gè)Perl程序處理。-p 標(biāo)志要求Perl循環(huán)處理看到的所有行， -i 標(biāo)志要求Perl原地更新匹配的行內(nèi)容，而 -e 標(biāo)志則要求Perl把傳入的字符串當(dāng)作一個(gè)Perl程序處理。g 修飾符則確保所有匹配都被更新，而不只是任何給定行里的第一個(gè)匹配。
# 2、diamond 超快蛋白序列比對(duì)
cd /opt/biosoft/
wget https://github.com/bbuchfink/diamond/releases/download/v2.1.5/diamond-linux64.tar.gz
tar xzf diamond-linux64.tar.gz
或者用conda 安裝
首先建立數(shù)據(jù)庫(kù)
diamond makedb --in Af_As.fasta --db afas
產(chǎn)生一個(gè)afas.dmnd 的文件
diamond blastp --db afas.dmnd --query Af_As.fasta --out AF_AS.balst --outfmt 6 --se
nsitive --max-target-seqs 10 --evalue 1e-5
主要參數(shù)說(shuō)明
--db/-d 輸入比對(duì)數(shù)據(jù)庫(kù)
--query/-q 比對(duì)序列
--threads/-p 線程數(shù)
--out/-o 輸出文件
--outfmt/-f 輸出文件格式，默認(rèn)6（表格）
--evalue/-e 比對(duì)的最大evalue值（默認(rèn)0.001)
--max-target-seqs/-k 比對(duì)到的最大序列數(shù)，默認(rèn)值是25
其他參數(shù)：
--top 百分?jǐn)?shù)的形式表示--max-target-seqs
--min-score 最小評(píng)分
--id 給出指定百分比的數(shù)據(jù)
--subject-cover 最小覆蓋度
--unal (0,1) 是否輸出未比對(duì)上的reads（0=no, 1=yes）
--sensitive 建議對(duì)齊較長(zhǎng)的序列
--more-sensitive 比對(duì)準(zhǔn)確度更高
--block-size/b，一次處理的十億堿基的大小，主要控制內(nèi)存使用，默認(rèn)為2（預(yù)計(jì)使用此內(nèi)存數(shù)量的大約六倍，即默認(rèn)內(nèi)存使用將到達(dá)12G），轉(zhuǎn)錄流程使用0.2
--salltitles 將全長(zhǎng)標(biāo)題包含在DAA格式中，默認(rèn)DAA文件僅包含縮短序列ID（直到第一個(gè)空白字符）
準(zhǔn)備AF_AS.gff 文件