導(dǎo)讀

用query基因、DNA序列、16S擴(kuò)增子序列、蛋白序列比對基因注釋、物種注釋 、基因功能注釋數(shù)據(jù)庫。Blast準(zhǔn)確性高，速度慢。Diamond速度更快，常用于二代測序短序列比對。

一、Diamond blastp

Diamond blastx功能可用于DNA序列比對【蛋白庫】

【核酸】比對【蛋白庫】

diamond blastx \
--db /path/to/protein.dmnd \
--query /path/to/genome.fna \
-e 1e-5 --outfmt 6 --more-sensitive --max-target-seqs 1 --threads 52 --quiet --id 80 --subject-cover 50 --query-cover 50 \ 
--out /path/to/anno.txt

【蛋白】比對【蛋白庫】

• diamond格式化蛋白庫

diamond makedb -p 52 \
--in /path/to/protein.fasta \
--db /path/to/protein.dmnd

• query比對格式化蛋白庫
• paraments

blastp 比蛋白
blastx 比核酸
blastp  Align amino acid query sequences against a protein reference database
blastx  Align DNA query sequences against a protein reference database
--query-cover            minimum query cover% to report an alignment
--subject-cover          minimum subject cover% to report an alignment
--out (-o) 輸出
--query (-q) 輸入
--db (-d) 格式化數(shù)據(jù)庫
--in Fasta數(shù)據(jù)庫
--quiet 安靜模式
--verbose (-v) 啰嗦模式
--outfmt (-f) 輸出格式 
        5   = BLAST XML
        6   = BLAST tabular
        100 = DIAMOND alignment archive (DAA)
        101 = SAM
--top / --max-target-seqs 保留個數(shù)
--evalue (-e) 最大evalue
--id 最小一致性
--sensitive 敏感模式 (default: fast)
--more-sensitive 更敏感模式 (default: fast)
--block-size (-b) sequence block size in billions of letters (default=2.0)
--index-chunks (-c) number of chunks for index processing

• command

diamond blastp \
--db /path/to/protein.dmnd \
--query /path/to/genome.faa \
-e 1e-5 --outfmt 6 --more-sensitive --max-target-seqs 1 --threads 52 --quiet --id 60 \
--out /path/to/anno.txt

Diamond比對結(jié)果

1. qseqid |  query  sequence id 
2. sseqid |  subject sequence id 
3. pident |  percentage of identical matches 
4. length |  alignment length 
5. mismatch |  number of mismatches 
6. gapopen |  number of gap openings 
7. q start |  start of alignment in query 
8. q end |  end of alignment in query 
9. s start |  start of alignment in subject 
10. s end |  end of alignment in subject 
11. evalue
12. bitscore

e value解讀

二、Blastn & Blastp

Blastn【基因】比對【基因庫】

注釋信息可以是物種，也可以是功能，比對就O了

• blast格式化數(shù)據(jù)庫

makeblastdb -in database.fna \
-input_type fasta -dbtype nucl -title prefix -parse_seqids \
 -out /path/to/prefix 

-dbtype <String, nucl, prot>
-input_type <String, asn1_bin, asn1_txt, blastdb, fasta>
-title <String> Title for BLAST database
-parse_seqids Option to parse seqid for FASTA input
-out <String> Name of BLAST database to be created

• query比對格式化數(shù)據(jù)庫

blastn -db /path/to/prefix \
-query /path/to/query.fasta \
-outfmt 6 \
-out /path/to/qcov80_id80_blast.out \
-num_threads 64 -perc_identity 80 -qcov_hsp_perc 80

-max_target_seqs: Maximum number of aligned sequences to keep. Default = 500
-evalue: E值就是Score值可靠性的評價。
它表明在隨機(jī)的情況下，其它序列與目標(biāo)序列相似度要大于S值的可能性。所以它的分值越低越好。
-perc_identity: 序列一致性
-qcov_hsp_perc <Real, 0..100> Percent query coverage per hsp。HSP是high scoring pair

Blastp【蛋白】比對【蛋白庫】

• blast格式化數(shù)據(jù)庫

makeblastdb -in /path/to/pro.fas \
-input_type fasta -dbtype prot 
-title blastn_pro \
-out /path/to/blastn_pro -parse_seqids

• query比對格式化數(shù)據(jù)庫

blastp -db /path/to/blastn_pro \
-query query.faa \
-out blast_pro.out \
-num_threads 52 -qcov_hsp_perc 80 -outfmt 6 \
-evalue 1e-5 -num_alignments 1

Blast比對結(jié)果

1、Query id：查詢序列ID標(biāo)識
2、Subject id：比對上的目標(biāo)序列ID標(biāo)識
3、% identity：序列比對的一致性百分比
4、alignment length：符合比對的比對區(qū)域的長度
5、mismatches：比對區(qū)域的錯配數(shù)
6、gap openings：比對區(qū)域的gap數(shù)目
7、q start：Query id起始位點(diǎn)
8、q end：Query id終止位點(diǎn)
9、s start：Subject id起始位點(diǎn)
10、s end：Subject id終止位點(diǎn)
11、e value：比對結(jié)果的期望值
12、bit score：比對結(jié)果的bit score值

The sensitivity can be adjusted using the options: 
--mid-sensitive, 
--sensitive, 
--more-sensitive,
--very-sensitive,
--ultra-sensitive.

Diamond 比對軟件使用
The default e-value cutoff of DIAMOND is 0.001 while that of BLAST is 10, so by default the program will search a lot more stringently than BLAST and not report weak hits.

三、比對軟件默認(rèn)取值

prokka

https://github.com/tseemann/prokka/blob/master/README.md

--evalue [n.n]    Similarity e-value cut-off (default '1e-06')

default '1e-06'

diamond

default e-value cutoff of DIAMOND is 0.001 while that of BLAST is 10

1e-3

blast

default e-value cutoff of DIAMOND is 0.001 while that of BLAST is 10

10

kofamscan

-E, --e-value <e_value>    Largest E-value required of the hits

親自測試是 1e-3

????