基因/基因组分析最重要的一步就是序列比对算法，包括全局比对和局部比对（global and local alignment），这对找到特定相似度的序列特别重要，也是限制计算速度的步骤。32位的usearch由于内存限制的原因，速度很快的local alignment算法并没有被广泛使用，但是64位不限制内存的版本则非常受到青睐。本人所在的研究组购买了64位的USEARCH的多个版本，包括早期的7.0.1版本和最新的11.0.1版本。USEARCH支持多线程是大家都知道的，但是官方的网站中说不支持多进程，而且购买一次（one license）作者建议只能跑一个进程，不可以多个进程同时跑。不同用户使用usearch也不能同时使用，但是实际上，作者并没有对同时跑多个进程进行限制，即使是只有一个license，还是同样可以一个用户多次同时掉用usearch。笔者在Ubuntu Linux平台同时调用usearch的-usearch_local和-ublast两个算法，各自使用了40个核（线程），还同时掉用了另外一个-usearch_local进程，总共使用了120个核 （线程）。速度没有受到任何影响，将百万条序列和一个数据库进行比对，半个小时之内就完成了（ublast啊只需要2-3分钟）。usearch的这两个算法完全可以代替blastn算法和blastp算法，而且速度较快，对于多个样品，使用usearch将会节约很多时间。如果哪位读者对usearch用于宏基因组分析非常感兴趣，欢迎交流。
usearch -fastx_relabel otus.fa -prefix Otu -fastaout relabeled.fa -keep_annots
nohup time usearch11 -fastx_demux Hawaii_ITS_S0_L001_R1_001.fastq -reverse Hawaii_ITS_S0_L001_R2_001.fastq -index Hawaii_ITS_S0_L001_I1_001.fastq -barcodes Hawaii_barcode_ITS.fa -fastqout fwd_demux1.fq -output2 rev_demux1.fq -sample_delim . -threads 80 &
usearch11 -fastq_mergepairs *_R1.fq -fastqout all_samples_merged1.fq -relabel @ -threads 80
usearch11 -fastx_subsample all_samples_merged1.fq -sample_size 5000 -fastqout all_sub_for_primer_check1.fq -threads 80
usearch11 -search_oligodb all_sub_for_primer_check1.fq -db primer1.fa -strand both -userout primer_hits1.txt -userfields query+qlo+qhi+qstrand -threads 80
usearch11 -fastq_filter all_samples_merged1.fq -fastq_stripleft 22 -fastq_stripright 20 -fastq_maxee 1 -fastq_minlen 220 -fastaout QCd_merged.fa -threads 80
usearch11 -fastq_filter all_samples_merged.fq -fastq_stripleft 22 -fastq_stripright 19 -fastq_maxee 1 -fastq_minlen 220 -fastaout QCd_merged.fa -threads 80
usearch11 -fastx_uniques QCd_merged.fa -fastaout uniques.fasta -sizeout -relabel Uniq -threads 80
usearch11 -fastx_uniques pico.bac.qiime.fa -fastaout pico.bac.usearch.uniques.fasta -sizeout -relabel Uniq -threads 96
usearch11 -cluster_otus uniques.fasta -otus otus.fa -relabel Otu -threads 80
usearch11 -cluster_otus uniques.fasta -otus otus_98.5.fa -relabel Otu -id 98.5
usearch11 -unoise3 uniques.fasta -zotus ASVs1.fa -minsize 8 -threads 80
sed -i.tmp 's/Zotu/ASV_/' ASVs.fa
nohup time usearch11 -usearch_global QCd_merged.fa -db ASVs.fa -id 0.99 -otutabout ASV_counts.txt -threads 80 -strand both &
nohup time usearch11 -usearch_global QCd_merged.fa -db otus.fa -id 0.99 -otutabout otus_counts.txt -threads 80 -strand both &
sed -i.tmp 's/#OTU ID//' ASV_counts.txt
usearch11 -sintax otus.fa -db utax_reference_dataset_10.10.2017.fasta -tabbedout otus_tax_unite.txt -strand plus -sintax_cutoff 0.5 -threads 80
usearch11 -sintax ASVs.fa -db utax_reference_dataset_10.10.2017.fasta -tabbedout ASVs_tax_unite.txt -strand plus -sintax_cutoff 0.5 -threads 80
bash convert_usearch_tax.sh
usearch11 -nbc_tax ASVs.fa -db rdp_16s_v16.fa -strand plus -tabbedout ASV_tax_rdp1.txt
usearch11 -nbc_tax ASVs.fa -db rdp_its_v2.fa -strand plus -tabbedout ASV_tax_rdp.txt
usearch11 -nbc_tax ASVs.fa -db utax_reference_dataset_10.10.2017.fasta -strand plus -tabbedout ASV_tax_unite.txt -threads 80
usearch11 -nbc_tax ASVs.fa -db silva_16s_v123.fa -strand plus -tabbedout tax_rdp_silva.txt -threads 80
prodigal -i X10_16N_H_N_contig_long.fa -o X10_16N_H_N_contig_long.genes -a X10_16N_H_N_contig_long.faa -p meta
prodigal -i contig.fa -o contig_genes -a contig.faa -p meta
perl run_FragGeneScan.pl -genome contig.fa -out lab5_contigs_genepred -complete 0 -thread 96
usearch11 -cluster_fast X10_16N_H_N_contig_long.faa -id 0.9 -centroids X10_16N_H_N_contig_long.my.proteins.centroids.fasta -uc X10_16N_H_N_contig_long.my.clusters.uc -threads 80
usearch11 -cluster_fast uniques.fasta -id 0.97 -centroids uniques.centroids.fasta -uc uniques.clusters.uc -threads 80
usearch11 -makeudb_ublast uniprot_sprot.fasta -output uniprot_sprot.udb -threads 80
nohup time usearch11 -ublast X10_16N_H_N_contig_long.my.proteins.centroids.fasta -db uniprot_sprot.udb -evalue 1e-9 -accel 0.75 -id 0.9 -maxaccepts 1 -threads 8 -userout X10_16N_H_N_contig_long_protein_hits_ublast.m8 -userfields query+target+id+alnlen+mism+opens+qlo+qhi+tlo+thi+evalue+bits+ql+tl &
usearch11 -makeudb_usearch uniprot_sprot.fasta -output uniprot_sprot_usearch.udb -threads 80
nohup time usearch11 -usearch_local X10_16N_H_N_contig_long.my.proteins.centroids._0.7.fasta -db uniprot_sprot_usearch.udb -evalue 1e-9 -accel 0.75 -id 0.9 -maxaccepts 1 -threads 80 -userout X10_16N_H_N_contig_long_protein_hits_local.m8 -userfields query+target+id+alnlen+mism+opens+qlo+qhi+tlo+thi+evalue+bits+ql+tl &
time usearch11 -usearch_local X10_16N_H_N_contig_long.my.proteins.centroids.fasta -db uniprot_sprot_usearch.udb -evalue 1e-9 -accel 0.75 -id 0.9 -maxaccepts 1 -threads 8 -userout X10_16N_H_N_contig_long_protein_hits_local.m8 -userfields query+target+id+alnlen+mism+opens+qlo+qhi+tlo+thi+evalue+bits+ql+tl
makeblastdb -in uniprot_sprot.fasta -dbtype prot -out uniprot_sprot.fasta
nohup time blastp -db uniprot_sprot.fasta -query X10_16N_H_N_contig_long.my.proteins.centroids._0.7.fasta -num_threads 80 -outfmt 6 -out X10_16N_H_N_contig_long.my.proteins.m8_blast -max_target_seqs 1 -evalue 1e-9 &
usearch11 -makeudb_usearch 97_otus.fasta -wordlength 13 -output gg97.udb -threads 80
usearch11 -udb2bitvec gg97.udb -output gg97.bitvec -threads 80
usearch11 -search_16s X10_16N_H_N_contig_long.fa -bitvec gg97.bitvec -fastaout X10_16N_H_N_contig_16s.fa -tabbedout X10_16N_H_N_contig_results.txt -threads 80
usearch11 -search_16s S_Baltica_OS155.fna -bitvec gg97.bitvec -fastaout S_Baltica_OS155_16s.fa -tabbedout S_Baltica_OS155_results.txt -threads 80
usearch11 -fastx_split X10_16N_H_N_contig_long.fa -splits 50 -outname split@.fa
以上是常见的usearch分析流程用到的代码，包括扩增子分析和宏基因组分析。欢迎去本人github查看细节和下载相关数据：
https://github.com/jianshu93github.com局部比对和全局比对的区别：
局部比对用于从较长的基因/基因组中找到最优的局部片段，简单理解为随意位置开始和结束的全局比对，需要关注比对的长度和原始的query序列长度的差异，全局比对用于长度相似的序列，比如特定的功能基因（amoA，16S rRNA等），不需要关注长度，只需要关注identity，确定query序列和数据库序列的差异。
扫一扫
2445
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
