構造予測の対象はRNAへと拡大したが、現在、深層アラインメントの生成がボトルネックとなっている。本著者らは、配列を重複するダイマー(2塩基単位)として表現する、RNAおよびDNA向けの検索・アラインメントツールRiboseekを開発した。RiboseekはBLASTNやnhmmerよりも高い感度で相同性を検出でき、アラインメント生成速度はnhmmerやrMSAと比較してそれぞれ250倍、376倍高速である。構造を考慮した再アラインメントを行うことで、rMSAと同等の塩基対復元精度を実現しつつ、処理速度は29倍高速に保たれている。また、構造を持つと推定されるRNAについて、ウェブサーバーおよび170万件の事前計算済みアラインメントを提供している。Riboseekは https://github.com/steineggerlab/riboseekにて無償で利用できる。
インストール
Win11 WSLでテストした。
Github
# Linux AVX2 build (check using: cat /proc/cpuinfo | grep avx2)
wget https://mmseqs.com/riboseek/riboseek-linux-avx2.tar.gz; tar xvzf riboseek-linux-avx2.tar.gz; export PATH=$(pwd)/riboseek/bin/:$PATH
# Linux AVX2 & GPU build (req. glibc >= 2.17 and nvidia driver >= 525.60.13)
wget https://mmseqs.com/riboseek/riboseek-linux-gpu.tar.gz; tar xvzf riboseek-linux-gpu.tar.gz; export PATH=$(pwd)/riboseek/bin/:$PATH
#conda
mamba create -n riboseek -c conda-forge -c bioconda riboseek -y
conda activate riboseek
> riboseek
Riboseek enables fast RNA homology search using covariance models and structure-aware alignment.
riboseek Version: 1.0.0
©
usage: riboseek <command> [<args>]
Easy workflows for plain text input/output
easy-search Sensitive homology search
easy-cluster Slower, sensitive clustering
easy-linclust Fast linear time cluster, less sensitive clustering
easy-proteomecluster Cluster proteomes and identify reference proteomes
easy-proteomesearch Calculate pairwise proteome similarity
easy-taxonomy Taxonomic classification
easy-rbh Find reciprocal best hit
Main workflows for database input/output
search RNA homology search with dinucleotide profiles
map Map nearly identical sequences
rbh Reciprocal best hit search
linclust Fast, less sensitive clustering
cluster Slower, sensitive clustering
clusterupdate Update previous clustering with new sequences
taxonomy Taxonomic classification
Input database creation
databases List and download databases
createdb Convert FASTA/Q file(s) to a sequence DB
createindex Store precomputed index on disk to reduce search overhead
convertmsa Convert Stockholm/PFAM MSA file to a MSA DB
msa2profile Convert a MSA DB to a profile DB
Format conversion for downstream processing
convertalis Convert alignment DB to BLAST-tab, SAM or custom format
createtsv Convert result DB to tab-separated flat file
convert2fasta Convert sequence DB to FASTA format
taxonomyreport Create a taxonomy report in Kraken or Krona format
An extended list of all modules can be obtained by calling 'riboseek -h'.
> riboseek easy-search -h
usage: riboseek easy-search <i:queryFastaFile1[.gz|.bz2]> ... <i:queryFastaFileN[.gz|.bz2]>|<i:stdin> <i:targetFastaFile[.gz]>|<i:targetDB> <o:alignmentFile> <tmpDir> [options]
options:
-s FLOAT Sensitivity: 1.0 faster; 4.0 fast; 7.5 sensitive [5.700]
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [300]
--alignment-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment [3]
--alignment-output-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment
5: score only (output) cluster format [0]
--min-seq-id FLOAT List matches above this sequence identity (for clustering) (range 0.0-1.0) [0.000]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
--alt-ali INT Show up to this many alternative alignments [0]
-c FLOAT List matches above this fraction of aligned (covered) residues (see --cov-mode) [0.000]
--cov-mode INT 0: coverage of query and target
1: coverage of target
2: coverage of query
3: target seq. length has to be at least x% of query length
4: query seq. length has to be at least x% of target length
5: short seq. needs to be at least x% of the other seq. length [0]
--max-accept INT Maximum accepted alignments before alignment calculation for a query is stopped [2147483647]
--num-iterations INT Number of iterative profile search iterations [1]
--taxon-list STR Taxonomy ID, possibly multiple values separated by ','
--prefilter-mode INT prefilter mode: 0: kmer/ungapped 1: ungapped, 2: nofilter, 3: ungapped&gapped [0]
--search-type INT Search type 0: auto 1: amino acid, 2: translated, 3: nucleotide, 4: translated nucleotide alignment [0]
--start-sens FLOAT Start sensitivity [4.000]
--sens-steps INT Number of search steps performed from --start-sens to -s [1]
--translation-mode INT Translation AA seq from nucleotide by 0: ORFs, 1: full reading frames [0]
--format-mode INT Output format:
0: BLAST-TAB
1: SAM
2: BLAST-TAB + query/db length
3: Pretty HTML
4: BLAST-TAB + column headers
BLAST-TAB (0) and BLAST-TAB + column headers (4) support custom output formats (--format-output) [0]
--format-output STR Choose comma separated list of output columns from: query,target,evalue,gapopen,pident,fident,nident,qstart,qend,qlen
tstart,tend,tlen,alnlen,raw,bits,cigar,qseq,tseq,qheader,theader,qaln,taln,qframe,tframe,mismatch,qcov,tcov
qset,qsetid,tset,tsetid,taxid,taxname,taxlineage,qorfstart,qorfend,torfstart,torfend,ppos [query,target,fident,alnlen,mismatch,gapopen,qstart,qend,tstart,tend,evalue,bits]
--greedy-best-hits BOOL Choose the best hits greedily to cover the query [0]
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--gpu-server INT Use GPU server [0]
--gpu-server-wait-timeout INT Wait for GPU server for 0: don't wait -1: no wait limit: >0 this many seconds [600]
examples:
# Search multiple FASTA against FASTA (like BLASTP, TBLASTN, BLASTX, BLASTN --search-type 3, TBLASTX --search-type 2)
mmseqs easy-search examples/QUERY.fasta examples/QUERY.fasta examples/DB.fasta result.m8 tmp
# Iterative profile search from stdin (like PSI-BLAST)
cat examples/QUERY.fasta | mmseqs easy-search stdin examples/DB.fasta result.m8 tmp --num-iterations 2
# Profile search against small databases (e.g. PFAM, eggNOG)
mmseqs databases PFAM pfam_db tmp
mmseqs easy-search examples/QUERY.fasta pfam_db res.m8 tmp
# Exhaustive search against sequences or profiles (works for large DBs)
mmseqs easy-search examples/QUERY.fasta targetProfiles res.m8 tmp --exhaustive-search
# Increasing sensitivity search (from 2 to 7 in 3 steps)
mmseqs easy-search examples/QUERY.fasta examples/DB.fasta result.m8 tmp --start-sens 2 -s 7 --sens-steps 3
references:
- Steinegger M, Soding J: MMseqs2 enables sensitive protein sequence searching for the analysis of massive data sets. Nature Biotechnology, 35(11), 1026-1028 (2017)
- Mirdita M, Steinegger M, Soding J: MMseqs2 desktop and local web server app for fast, interactive sequence searches. Bioinformatics, 35(16), 2856-2858 (2019)
- Kallenborn F, Chacon A, Hundt C, Sirelkhatim H, Didi K, Cha S, Dallago C, Mirdita M, Schmidt B, Steinegger M: GPU-accelerated homology search with MMseqs2. Nature Methods (2025)
Show an extended list of options by calling 'riboseek easy-search -h'.
> riboseek search -h
usage: riboseek search <i:queryDB> <i:targetDB> <o:resultDB> <tmpDir> [options]
options:
-s FLOAT Sensitivity: 1.0 faster; 4.0 fast; 7.5 sensitive [7.500]
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [300]
-a BOOL Add backtrace string (convert to alignments with mmseqs convertalis module) [0]
--alignment-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment [2]
--alignment-output-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment
5: score only (output) cluster format [0]
-e DOUBLE List matches below this E-value (range 0.0-inf) [1.000E-03]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
--alt-ali INT Show up to this many alternative alignments [0]
--max-rejected INT Maximum rejected alignments before alignment calculation for a query is stopped [2147483647]
--max-accept INT Maximum accepted alignments before alignment calculation for a query is stopped [2147483647]
--e-profile DOUBLE Include sequences matches with < E-value thr. into the profile (>=0.0) [1.000E-01]
--num-iterations INT Number of iterative profile search iterations [1]
--taxon-list STR Taxonomy ID, possibly multiple values separated by ','
--prefilter-mode INT prefilter mode: 0: kmer/ungapped 1: ungapped, 2: nofilter, 3: ungapped&gapped [0]
--rescore-mode INT Rescore diagonals with:
0: Hamming distance
1: local alignment (score only)
2: local alignment
3: global alignment
4: longest alignment fulfilling window quality criterion [0]
--allow-deletion BOOL Allow deletions in a MSA [0]
--min-length INT Minimum codon number in open reading frames [30]
--max-length INT Maximum codon number in open reading frames [32734]
--max-gaps INT Maximum number of codons with gaps or unknown residues before an open reading frame is rejected [2147483647]
--contig-start-mode INT Contig start can be 0: incomplete, 1: complete, 2: both [2]
--contig-end-mode INT Contig end can be 0: incomplete, 1: complete, 2: both [2]
--orf-start-mode INT Orf fragment can be 0: from start to stop, 1: from any to stop, 2: from last encountered start to stop (no start in the middle) [1]
--forward-frames STR Comma-separated list of frames on the forward strand to be extracted [1,2,3]
--reverse-frames STR Comma-separated list of frames on the reverse strand to be extracted [1,2,3]
--translate INT Translate ORF to amino acid [0]
--use-all-table-starts BOOL Use all alternatives for a start codon in the genetic table, if false - only ATG (AUG) [0]
--id-offset INT Numeric ids in index file are offset by this value [0]
--search-type INT Search type 0: auto 1: amino acid, 2: translated, 3: nucleotide, 4: translated nucleotide alignment [3]
--start-sens FLOAT Start sensitivity [4.000]
--sens-steps INT Number of search steps performed from --start-sens to -s [1]
--translation-mode INT Translation AA seq from nucleotide by 0: ORFs, 1: full reading frames [0]
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--gpu INT Use GPU (CUDA) if possible [0]
--gpu-server INT Use GPU server [0]
--gpu-server-wait-timeout INT Wait for GPU server for 0: don't wait -1: no wait limit: >0 this many seconds [600]
examples:
# Search RNA query against RNA target database
riboseek search queryDB targetDB resultDB tmp
# Iterative profile search
riboseek search queryDB targetDB resultDB tmp --num-iterations 3
references:
- Steinegger M, Soding J: MMseqs2 enables sensitive protein sequence searching for the analysis of massive data sets. Nature Biotechnology, 35(11), 1026-1028 (2017)
Show an extended list of options by calling 'riboseek search -h'.
> riboseek map -h
usage: riboseek map <i:queryDB> <i:targetDB> <o:alignmentDB> <tmpDir> [options]
By Milot Mirdita <milot@mirdita.de> & Martin Steinegger <martin.steinegger@snu.ac.kr>
options: prefilter:
--seed-sub-mat TWIN Substitution matrix file for k-mer generation [aa:dinuc.out,nucl:dinuc.out]
-s FLOAT Sensitivity: 1.0 faster; 4.0 fast; 7.5 sensitive [2.000]
-k INT k-mer length (0: automatically set to optimum) [0]
--target-search-mode INT target search mode (0: regular k-mer, 1: similar k-mer) [0]
--k-score TWIN k-mer threshold for generating similar k-mer lists [seq:2147483647,prof:2147483647]
--alph-size TWIN Alphabet size (range 2-21) [aa:25,nucl:5]
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [300]
--split INT Split input into N equally distributed chunks. 0: set the best split automatically [0]
--split-mode INT 0: split target db; 1: split query db; 2: auto, depending on main memory [2]
--split-memory-limit BYTE Set max memory per split. E.g. 800B, 5K, 10M, 1G. Default (0) to all available system memory [0]
--comp-bias-corr INT Correct for locally biased amino acid composition (range 0-1) [0]
--comp-bias-corr-scale FLOAT Correct for locally biased amino acid composition (range 0-1) [1.000]
--diag-score BOOL Use ungapped diagonal scoring during prefilter [1]
--exact-kmer-matching INT Extract only exact k-mers for matching (range 0-1) [0]
--mask INT Mask sequences in prefilter stage with tantan: 0: w/o low complexity masking, 1: with low complexity masking [0]
--mask-prob FLOAT Mask sequences is probablity is above threshold [0.900]
--mask-lower-case INT Lowercase letters will be excluded from k-mer search 0: include region, 1: exclude region [0]
--mask-n-repeat INT Repeat letters that occur > threshold in a rwo [0]
--min-ungapped-score INT Accept only matches with ungapped alignment score above threshold [15]
--add-self-matches BOOL Artificially add entries of queries with themselves (for clustering) [0]
--spaced-kmer-mode INT 0: use consecutive positions in k-mers; 1: use spaced k-mers [1]
--spaced-kmer-pattern STR User-specified spaced k-mer pattern
--local-tmp STR Path where some of the temporary files will be created
align:
-c FLOAT List matches above this fraction of aligned (covered) residues (see --cov-mode) [0.950]
--cov-mode INT 0: coverage of query and target
1: coverage of target
2: coverage of query
3: target seq. length has to be at least x% of query length
4: query seq. length has to be at least x% of target length
5: short seq. needs to be at least x% of the other seq. length [2]
--wrapped-scoring BOOL Double the (nucleotide) query sequence during the scoring process to allow wrapped diagonal scoring around end and start [0]
-e DOUBLE List matches below this E-value (range 0.0-inf) [1.000E-03]
-a BOOL Add backtrace string (convert to alignments with mmseqs convertalis module) [0]
--min-seq-id FLOAT List matches above this sequence identity (for clustering) (range 0.0-1.0) [0.900]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
profile:
--pca Pseudo count admixture strength
--pcb Pseudo counts: Neff at half of maximum admixture (range 0.0-inf)
misc:
--taxon-list STR Taxonomy ID, possibly multiple values separated by ','
--rescore-mode INT Rescore diagonals with:
0: Hamming distance
1: local alignment (score only)
2: local alignment
3: global alignment
4: longest alignment fulfilling window quality criterion [2]
--min-length INT Minimum codon number in open reading frames [10]
--max-length INT Maximum codon number in open reading frames [32734]
--max-gaps INT Maximum number of codons with gaps or unknown residues before an open reading frame is rejected [2147483647]
--contig-start-mode INT Contig start can be 0: incomplete, 1: complete, 2: both [2]
--contig-end-mode INT Contig end can be 0: incomplete, 1: complete, 2: both [2]
--orf-start-mode INT Orf fragment can be 0: from start to stop, 1: from any to stop, 2: from last encountered start to stop (no start in the middle) [1]
--forward-frames STR Comma-separated list of frames on the forward strand to be extracted [1,2,3]
--reverse-frames STR Comma-separated list of frames on the reverse strand to be extracted [1,2,3]
--translation-table INT 1) CANONICAL, 2) VERT_MITOCHONDRIAL, 3) YEAST_MITOCHONDRIAL, 4) MOLD_MITOCHONDRIAL, 5) INVERT_MITOCHONDRIAL, 6) CILIATE
9) FLATWORM_MITOCHONDRIAL, 10) EUPLOTID, 11) PROKARYOTE, 12) ALT_YEAST, 13) ASCIDIAN_MITOCHONDRIAL, 14) ALT_FLATWORM_MITOCHONDRIAL
15) BLEPHARISMA, 16) CHLOROPHYCEAN_MITOCHONDRIAL, 21) TREMATODE_MITOCHONDRIAL, 22) SCENEDESMUS_MITOCHONDRIAL
23) THRAUSTOCHYTRIUM_MITOCHONDRIAL, 24) PTEROBRANCHIA_MITOCHONDRIAL, 25) GRACILIBACTERIA, 26) PACHYSOLEN, 27) KARYORELICT, 28) CONDYLOSTOMA
29) MESODINIUM, 30) PERTRICH, 31) BLASTOCRITHIDIA [1]
--translate INT Translate ORF to amino acid [0]
--use-all-table-starts BOOL Use all alternatives for a start codon in the genetic table, if false - only ATG (AUG) [0]
--id-offset INT Numeric ids in index file are offset by this value [0]
--start-sens FLOAT Start sensitivity [4.000]
--sens-steps INT Number of search steps performed from --start-sens to -s [1]
common:
--sub-mat TWIN Substitution matrix file [aa:dinuc.out,nucl:dinuc.out]
--max-seq-len INT Maximum sequence length [65535]
--db-load-mode INT Database preload mode 0: auto, 1: fread, 2: mmap, 3: mmap+touch [0]
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--mpi-runner STR Use MPI on compute cluster with this MPI command (e.g. "mpirun -np 42")
--force-reuse BOOL Reuse tmp filse in tmp/latest folder ignoring parameters and version changes [0]
--remove-tmp-files BOOL Delete temporary files [0]
expert:
--filter-hits BOOL Filter hits by seq.id. and coverage [0]
--sort-results INT Sort results: 0: no sorting, 1: sort by E-value (Alignment) or seq.id. (Hamming) [1]
--create-lookup INT Create database lookup file (can be very large) [0]
references:
- Steinegger M, Soding J: MMseqs2 enables sensitive protein sequence searching for the analysis of massive data sets. Nature Biotechnology, 35(11), 1026-1028 (2017)
- Steinegger M, Mirdita M, Soding J: Protein-level assembly increases protein sequence recovery from metagenomic samples manyfold. Nature Methods, 16(7), 603-606 (2019)
> riboseek rbh -h
sage: riboseek rbh <i:queryDB> <i:targetDB> <o:alignmentDB> <tmpDir> [options]
By Eli Levy Karin
options: prefilter:
--comp-bias-corr INT Correct for locally biased amino acid composition (range 0-1) [0]
--comp-bias-corr-scale FLOAT Correct for locally biased amino acid composition (range 0-1) [1.000]
--add-self-matches BOOL Artificially add entries of queries with themselves (for clustering) [0]
--seed-sub-mat TWIN Substitution matrix file for k-mer generation [aa:dinuc.out,nucl:dinuc.out]
-s FLOAT Sensitivity: 1.0 faster; 4.0 fast; 7.5 sensitive [4.000]
-k INT k-mer length (0: automatically set to optimum) [0]
--target-search-mode INT target search mode (0: regular k-mer, 1: similar k-mer) [0]
--k-score TWIN k-mer threshold for generating similar k-mer lists [seq:2147483647,prof:2147483647]
--alph-size TWIN Alphabet size (range 2-21) [aa:25,nucl:5]
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [300]
--split INT Split input into N equally distributed chunks. 0: set the best split automatically [0]
--split-mode INT 0: split target db; 1: split query db; 2: auto, depending on main memory [2]
--split-memory-limit BYTE Set max memory per split. E.g. 800B, 5K, 10M, 1G. Default (0) to all available system memory [0]
--diag-score BOOL Use ungapped diagonal scoring during prefilter [1]
--exact-kmer-matching INT Extract only exact k-mers for matching (range 0-1) [0]
--mask INT Mask sequences in prefilter stage with tantan: 0: w/o low complexity masking, 1: with low complexity masking [0]
--mask-prob FLOAT Mask sequences is probablity is above threshold [0.900]
--mask-lower-case INT Lowercase letters will be excluded from k-mer search 0: include region, 1: exclude region [0]
--mask-n-repeat INT Repeat letters that occur > threshold in a rwo [0]
--min-ungapped-score INT Accept only matches with ungapped alignment score above threshold [15]
--spaced-kmer-mode INT 0: use consecutive positions in k-mers; 1: use spaced k-mers [1]
--spaced-kmer-pattern STR User-specified spaced k-mer pattern
--local-tmp STR Path where some of the temporary files will be created
--disk-space-limit BYTE Set max disk space to use for reverse profile searches. E.g. 800B, 5K, 10M, 1G. Default (0) to all available disk space in the temp folder [0]
align:
-a BOOL Add backtrace string (convert to alignments with mmseqs convertalis module) [0]
--alignment-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment [3]
--alignment-output-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment
5: score only (output) cluster format [0]
--wrapped-scoring BOOL Double the (nucleotide) query sequence during the scoring process to allow wrapped diagonal scoring around end and start [0]
-e DOUBLE List matches below this E-value (range 0.0-inf) [1.000E-03]
--min-seq-id FLOAT List matches above this sequence identity (for clustering) (range 0.0-1.0) [0.000]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
--alt-ali INT Show up to this many alternative alignments [0]
-c FLOAT List matches above this fraction of aligned (covered) residues (see --cov-mode) [0.000]
--cov-mode INT 0: coverage of query and target
1: coverage of target
2: coverage of query
3: target seq. length has to be at least x% of query length
4: query seq. length has to be at least x% of target length
5: short seq. needs to be at least x% of the other seq. length [0]
--max-rejected INT Maximum rejected alignments before alignment calculation for a query is stopped [2147483647]
--max-accept INT Maximum accepted alignments before alignment calculation for a query is stopped [2147483647]
--score-bias FLOAT Score bias when computing SW alignment (in bits) [0.000]
--realign BOOL Compute more conservative, shorter alignments (scores and E-values not changed) [0]
--realign-score-bias FLOAT Additional bias when computing realignment [-0.200]
--realign-max-seqs INT Maximum number of results to return in realignment [2147483647]
--corr-score-weight FLOAT Weight of backtrace correlation score that is added to the alignment score [0.000]
--gap-open TWIN Gap open cost [aa:11,nucl:5]
--gap-extend TWIN Gap extension cost [aa:1,nucl:2]
--zdrop INT Maximal allowed difference between score values before alignment is truncated (nucleotide alignment only) [40]
--exhaustive-search-filter INT Filter result during search: 0: do not filter, 1: filter [0]
profile:
--pca Pseudo count admixture strength
--pcb Pseudo counts: Neff at half of maximum admixture (range 0.0-inf)
--mask-profile INT Mask query sequence of profile using tantan [0,1] [1]
--e-profile DOUBLE Include sequences matches with < E-value thr. into the profile (>=0.0) [1.000E-03]
--wg BOOL Use global sequence weighting for profile calculation [0]
--filter-msa INT Filter msa: 0: do not filter, 1: filter [1]
--filter-min-enable INT Only filter MSAs with more than N sequences, 0 always filters [0]
--max-seq-id FLOAT Reduce redundancy of output MSA using max. pairwise sequence identity [0.0,1.0] [0.900]
--qid STR Reduce diversity of output MSAs using min.seq. identity with query sequences [0.0,1.0]
Alternatively, can be a list of multiple thresholds:
E.g.: 0.15,0.30,0.50 to defines filter buckets of ]0.15-0.30] and ]0.30-0.50] [0.0]
--qsc FLOAT Reduce diversity of output MSAs using min. score per aligned residue with query sequences [-50.0,100.0] [-20.000]
--cov FLOAT Filter output MSAs using min. fraction of query residues covered by matched sequences [0.0,1.0] [0.000]
--diff INT Filter MSAs by selecting most diverse set of sequences, keeping at least this many seqs in each MSA block of length 50 [1000]
--pseudo-cnt-mode INT use 0: substitution-matrix or 1: context-specific pseudocounts [0]
--profile-output-mode INT Profile output mode: 0: binary log-odds 1: human-readable frequencies [0]
--num-iterations INT Number of iterative profile search iterations [1]
--exhaustive-search BOOL For bigger profile DB, run iteratively the search by greedily swapping the search results [0]
--lca-search BOOL Efficient search for LCA candidates [0]
misc:
--taxon-list STR Taxonomy ID, possibly multiple values separated by ','
--prefilter-mode INT prefilter mode: 0: kmer/ungapped 1: ungapped, 2: nofilter, 3: ungapped&gapped [0]
--rescore-mode INT Rescore diagonals with:
0: Hamming distance
1: local alignment (score only)
2: local alignment
3: global alignment
4: longest alignment fulfilling window quality criterion [0]
--allow-deletion BOOL Allow deletions in a MSA [0]
--min-length INT Minimum codon number in open reading frames [10]
--max-length INT Maximum codon number in open reading frames [32734]
--max-gaps INT Maximum number of codons with gaps or unknown residues before an open reading frame is rejected [2147483647]
--contig-start-mode INT Contig start can be 0: incomplete, 1: complete, 2: both [2]
--contig-end-mode INT Contig end can be 0: incomplete, 1: complete, 2: both [2]
--orf-start-mode INT Orf fragment can be 0: from start to stop, 1: from any to stop, 2: from last encountered start to stop (no start in the middle) [1]
--forward-frames STR Comma-separated list of frames on the forward strand to be extracted [1,2,3]
--reverse-frames STR Comma-separated list of frames on the reverse strand to be extracted [1,2,3]
--translation-table INT 1) CANONICAL, 2) VERT_MITOCHONDRIAL, 3) YEAST_MITOCHONDRIAL, 4) MOLD_MITOCHONDRIAL, 5) INVERT_MITOCHONDRIAL, 6) CILIATE
9) FLATWORM_MITOCHONDRIAL, 10) EUPLOTID, 11) PROKARYOTE, 12) ALT_YEAST, 13) ASCIDIAN_MITOCHONDRIAL, 14) ALT_FLATWORM_MITOCHONDRIAL
15) BLEPHARISMA, 16) CHLOROPHYCEAN_MITOCHONDRIAL, 21) TREMATODE_MITOCHONDRIAL, 22) SCENEDESMUS_MITOCHONDRIAL
23) THRAUSTOCHYTRIUM_MITOCHONDRIAL, 24) PTEROBRANCHIA_MITOCHONDRIAL, 25) GRACILIBACTERIA, 26) PACHYSOLEN, 27) KARYORELICT, 28) CONDYLOSTOMA
29) MESODINIUM, 30) PERTRICH, 31) BLASTOCRITHIDIA [1]
--translate INT Translate ORF to amino acid [0]
--use-all-table-starts BOOL Use all alternatives for a start codon in the genetic table, if false - only ATG (AUG) [0]
--id-offset INT Numeric ids in index file are offset by this value [0]
--sequence-overlap INT Overlap between sequences [0]
--sequence-split-mode INT Sequence split mode 0: copy data, 1: soft link data and write new index, [1]
--headers-split-mode INT Header split mode: 0: split position, 1: original header [0]
--search-type INT Search type 0: auto 1: amino acid, 2: translated, 3: nucleotide, 4: translated nucleotide alignment [0]
--start-sens FLOAT Start sensitivity [4.000]
--sens-steps INT Number of search steps performed from --start-sens to -s [1]
--translation-mode INT Translation AA seq from nucleotide by 0: ORFs, 1: full reading frames [0]
common:
--sub-mat TWIN Substitution matrix file [aa:dinuc.out,nucl:dinuc.out]
--max-seq-len INT Maximum sequence length [65535]
--db-load-mode INT Database preload mode 0: auto, 1: fread, 2: mmap, 3: mmap+touch [0]
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--gpu INT Use GPU (CUDA) if possible [0]
--gpu-server INT Use GPU server [0]
--gpu-server-wait-timeout INT Wait for GPU server for 0: don't wait -1: no wait limit: >0 this many seconds [600]
--mpi-runner STR Use MPI on compute cluster with this MPI command (e.g. "mpirun -np 42")
--force-reuse BOOL Reuse tmp filse in tmp/latest folder ignoring parameters and version changes [0]
--remove-tmp-files BOOL Delete temporary files [0]
expert:
--filter-hits BOOL Filter hits by seq.id. and coverage [0]
--sort-results INT Sort results: 0: no sorting, 1: sort by E-value (Alignment) or seq.id. (Hamming) [0]
--create-lookup INT Create database lookup file (can be very large) [0]
--chain-alignments INT Chain overlapping alignments [0]
--merge-query INT Combine ORFs/split sequences to a single entry [1]
--strand INT Strand selection only works for DNA/DNA search 0: reverse, 1: forward, 2: both [1]
references:
- Steinegger M, Soding J: MMseqs2 enables sensitive protein sequence searching for the analysis of massive data sets. Nature Biotechnology, 35(11), 1026-1028 (2017)
> riboseek cluster -h
usage: riboseek cluster <i:sequenceDB> <o:clusterDB> <tmpDir> [options]
By Martin Steinegger <martin.steinegger@snu.ac.kr> & Lars von den Driesch
options: prefilter:
--seed-sub-mat TWIN Substitution matrix file for k-mer generation [aa:dinuc.out,nucl:dinuc.out]
-s FLOAT Sensitivity: 1.0 faster; 4.0 fast; 7.5 sensitive [4.000]
-k INT k-mer length (0: automatically set to optimum) [0]
--target-search-mode INT target search mode (0: regular k-mer, 1: similar k-mer) [0]
--k-score TWIN k-mer threshold for generating similar k-mer lists [seq:2147483647,prof:2147483647]
--alph-size TWIN Alphabet size (range 2-21) [aa:25,nucl:5]
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [20]
--split INT Split input into N equally distributed chunks. 0: set the best split automatically [0]
--split-mode INT 0: split target db; 1: split query db; 2: auto, depending on main memory [2]
--split-memory-limit BYTE Set max memory per split. E.g. 800B, 5K, 10M, 1G. Default (0) to all available system memory [0]
--comp-bias-corr INT Correct for locally biased amino acid composition (range 0-1) [1]
--comp-bias-corr-scale FLOAT Correct for locally biased amino acid composition (range 0-1) [1.000]
--diag-score BOOL Use ungapped diagonal scoring during prefilter [1]
--exact-kmer-matching INT Extract only exact k-mers for matching (range 0-1) [0]
--mask INT Mask sequences in prefilter stage with tantan: 0: w/o low complexity masking, 1: with low complexity masking [0]
--mask-prob FLOAT Mask sequences is probablity is above threshold [0.900]
--mask-lower-case INT Lowercase letters will be excluded from k-mer search 0: include region, 1: exclude region [0]
--mask-n-repeat INT Repeat letters that occur > threshold in a rwo [0]
--min-ungapped-score INT Accept only matches with ungapped alignment score above threshold [15]
--add-self-matches BOOL Artificially add entries of queries with themselves (for clustering) [0]
--spaced-kmer-mode INT 0: use consecutive positions in k-mers; 1: use spaced k-mers [1]
--spaced-kmer-pattern STR User-specified spaced k-mer pattern
--local-tmp STR Path where some of the temporary files will be created
align:
-c FLOAT List matches above this fraction of aligned (covered) residues (see --cov-mode) [0.800]
--cov-mode INT 0: coverage of query and target
1: coverage of target
2: coverage of query
3: target seq. length has to be at least x% of query length
4: query seq. length has to be at least x% of target length
5: short seq. needs to be at least x% of the other seq. length [0]
-a BOOL Add backtrace string (convert to alignments with mmseqs convertalis module) [0]
--alignment-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment [3]
--alignment-output-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment
5: score only (output) cluster format [0]
--wrapped-scoring BOOL Double the (nucleotide) query sequence during the scoring process to allow wrapped diagonal scoring around end and start [0]
-e DOUBLE List matches below this E-value (range 0.0-inf) [1.000E-03]
--min-seq-id FLOAT List matches above this sequence identity (for clustering) (range 0.0-1.0) [0.000]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
--alt-ali INT Show up to this many alternative alignments [0]
--max-rejected INT Maximum rejected alignments before alignment calculation for a query is stopped [2147483647]
--max-accept INT Maximum accepted alignments before alignment calculation for a query is stopped [2147483647]
--score-bias FLOAT Score bias when computing SW alignment (in bits) [0.000]
--realign BOOL Compute more conservative, shorter alignments (scores and E-values not changed) [0]
--realign-score-bias FLOAT Additional bias when computing realignment [-0.200]
--realign-max-seqs INT Maximum number of results to return in realignment [2147483647]
--corr-score-weight FLOAT Weight of backtrace correlation score that is added to the alignment score [0.000]
--gap-open TWIN Gap open cost [aa:11,nucl:5]
--gap-extend TWIN Gap extension cost [aa:1,nucl:2]
--zdrop INT Maximal allowed difference between score values before alignment is truncated (nucleotide alignment only) [40]
clust:
--cluster-mode INT 0: Set-Cover (greedy)
1: Connected component (BLASTclust)
2,3: Greedy clustering by sequence length (CDHIT) [0]
--max-iterations INT Maximum depth of breadth first search in connected component clustering [1000]
--similarity-type INT Type of score used for clustering. 1: alignment score 2: sequence identity [2]
--set-mode BOOL 0: Cluster by each entry
1: Cluster by set [0]
--single-step-clustering BOOL Switch from cascaded to simple clustering workflow [0]
--cluster-steps INT Cascaded clustering steps from 1 to -s [3]
--cluster-reassign BOOL Cascaded clustering can cluster sequence that do not fulfill the clustering criteria.
Cluster reassignment corrects these errors [0]
kmermatcher:
--weights STR Weights used for cluster priorization
--cluster-weight-threshold FLOAT Weight threshold used for cluster priorization [0.900]
--cluster-version INT Cluster version: 1: Cluster1, 2: Cluster2 [1]
--kmer-per-seq INT k-mers per sequence [21]
--kmer-per-seq-scale TWIN Scale k-mer per sequence based on sequence length as kmer-per-seq val + scale x seqlen [aa:0.000,nucl:0.200]
--adjust-kmer-len BOOL Adjust k-mer length based on specificity (only for nucleotides) [0]
--hash-shift INT Shift k-mer hash initialization [67]
--include-only-extendable BOOL Include only extendable [0]
--ignore-multi-kmer BOOL Skip k-mers occurring multiple times (>=2) [0]
--include-count-table BOOL Include count table [1]
--num-count-table INT Number of count table based center swapping [2]
--include-adjacency BOOL Include adjacency based center swapping [1]
--num-adjacency INT Number of adjacency based center swapping [3]
--use-parallelism BOOL Enable or disable parallel execution for group assignment and related k-mer processing steps [0]
--linclust-version INT Linclust version: 1: Linclust1, 2: Linclust2 [2]
--clust-hash BOOL Use clusthash before kmermatcher in linclust [0]
profile:
--pca Pseudo count admixture strength
--pcb Pseudo counts: Neff at half of maximum admixture (range 0.0-inf)
misc:
--taxon-list STR Taxonomy ID, possibly multiple values separated by ','
--rescore-mode INT Rescore diagonals with:
0: Hamming distance
1: local alignment (score only)
2: local alignment
3: global alignment
4: longest alignment fulfilling window quality criterion [0]
common:
--sub-mat TWIN Substitution matrix file [aa:dinuc.out,nucl:dinuc.out]
--max-seq-len INT Maximum sequence length [65535]
--db-load-mode INT Database preload mode 0: auto, 1: fread, 2: mmap, 3: mmap+touch [0]
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--remove-tmp-files BOOL Delete temporary files [0]
--force-reuse BOOL Reuse tmp filse in tmp/latest folder ignoring parameters and version changes [0]
--mpi-runner STR Use MPI on compute cluster with this MPI command (e.g. "mpirun -np 42")
expert:
--filter-hits BOOL Filter hits by seq.id. and coverage [0]
--sort-results INT Sort results: 0: no sorting, 1: sort by E-value (Alignment) or seq.id. (Hamming) [0]
--filter-cludb-file STR Specify a cluDB file to filter
--filter-seqdb-file STR Specify a seqDB of cluDB file to filter
examples:
# Cascaded clustering of FASTA file
mmseqs cluster sequenceDB clusterDB tmp
# --cov-mode
# Sequence 0 1 2
# Q: MAVGTACRPA 60% IGN 60%
# T: -AVGTAC--- 60% 100% IGN
# Cutoff -c 0.7 - + -
# -c 0.6 + + +
# Cascaded clustering with reassignment
# - Corrects criteria-violations of cascaded merging
# - Produces more clusters and is a bit slower
mmseqs cluster sequenceDB clusterDB tmp --cluster-reassign
references:
- Steinegger M, Soding J: MMseqs2 enables sensitive protein sequence searching for the analysis of massive data sets. Nature Biotechnology, 35(11), 1026-1028 (2017)
- Hauser M, Steinegger M, Soding J: MMseqs software suite for fast and deep clustering and searching of large protein sequence sets. Bioinformatics, 32(9), 1323-1330 (2016)
- Steinegger M, Soding J: Clustering huge protein sequence sets in linear time. Nature Communications, 9(1), 2542 (2018)
> riboseek linclust -h
usage: riboseek linclust <i:sequenceDB> <o:clusterDB> <tmpDir> [options]
By Martin Steinegger <martin.steinegger@snu.ac.kr>
options: prefilter:
--comp-bias-corr INT Correct for locally biased amino acid composition (range 0-1) [1]
--comp-bias-corr-scale FLOAT Correct for locally biased amino acid composition (range 0-1) [1.000]
--add-self-matches BOOL Artificially add entries of queries with themselves (for clustering) [0]
--alph-size TWIN Alphabet size (range 2-21) [aa:25,nucl:5]
--spaced-kmer-mode INT 0: use consecutive positions in k-mers; 1: use spaced k-mers [0]
--spaced-kmer-pattern STR User-specified spaced k-mer pattern
--mask INT Mask sequences in prefilter stage with tantan: 0: w/o low complexity masking, 1: with low complexity masking [0]
--mask-prob FLOAT Mask sequences is probablity is above threshold [0.900]
--mask-lower-case INT Lowercase letters will be excluded from k-mer search 0: include region, 1: exclude region [0]
--mask-n-repeat INT Repeat letters that occur > threshold in a rwo [0]
-k INT k-mer length (0: automatically set to optimum) [0]
--split-memory-limit BYTE Set max memory per split. E.g. 800B, 5K, 10M, 1G. Default (0) to all available system memory [0]
align:
-a BOOL Add backtrace string (convert to alignments with mmseqs convertalis module) [0]
--alignment-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment [3]
--alignment-output-mode INT How to compute the alignment:
0: automatic
1: only score and end_pos
2: also start_pos and cov
3: also seq.id
4: only ungapped alignment
5: score only (output) cluster format [0]
--wrapped-scoring BOOL Double the (nucleotide) query sequence during the scoring process to allow wrapped diagonal scoring around end and start [0]
-e DOUBLE List matches below this E-value (range 0.0-inf) [1.000E-03]
--min-seq-id FLOAT List matches above this sequence identity (for clustering) (range 0.0-1.0) [0.900]
--min-aln-len INT Minimum alignment length (range 0-INT_MAX) [0]
--seq-id-mode INT 0: alignment length 1: shorter, 2: longer sequence [0]
--alt-ali INT Show up to this many alternative alignments [0]
-c FLOAT List matches above this fraction of aligned (covered) residues (see --cov-mode) [0.800]
--cov-mode INT 0: coverage of query and target
1: coverage of target
2: coverage of query
3: target seq. length has to be at least x% of query length
4: query seq. length has to be at least x% of target length
5: short seq. needs to be at least x% of the other seq. length [0]
--max-rejected INT Maximum rejected alignments before alignment calculation for a query is stopped [2147483647]
--max-accept INT Maximum accepted alignments before alignment calculation for a query is stopped [2147483647]
--score-bias FLOAT Score bias when computing SW alignment (in bits) [0.000]
--realign BOOL Compute more conservative, shorter alignments (scores and E-values not changed) [0]
--realign-score-bias FLOAT Additional bias when computing realignment [-0.200]
--realign-max-seqs INT Maximum number of results to return in realignment [2147483647]
--corr-score-weight FLOAT Weight of backtrace correlation score that is added to the alignment score [0.000]
--gap-open TWIN Gap open cost [aa:11,nucl:5]
--gap-extend TWIN Gap extension cost [aa:1,nucl:2]
--zdrop INT Maximal allowed difference between score values before alignment is truncated (nucleotide alignment only) [40]
clust:
--cluster-mode INT 0: Set-Cover (greedy)
1: Connected component (BLASTclust)
2,3: Greedy clustering by sequence length (CDHIT) [0]
--max-iterations INT Maximum depth of breadth first search in connected component clustering [1000]
--similarity-type INT Type of score used for clustering. 1: alignment score 2: sequence identity [2]
--set-mode BOOL 0: Cluster by each entry
1: Cluster by set [0]
kmermatcher:
--weights STR Weights used for cluster priorization
--cluster-weight-threshold FLOAT Weight threshold used for cluster priorization [0.900]
--kmer-per-seq INT k-mers per sequence [21]
--kmer-per-seq-scale TWIN Scale k-mer per sequence based on sequence length as kmer-per-seq val + scale x seqlen [aa:0.000,nucl:0.200]
--adjust-kmer-len BOOL Adjust k-mer length based on specificity (only for nucleotides) [0]
--hash-shift INT Shift k-mer hash initialization [67]
--include-only-extendable BOOL Include only extendable [0]
--ignore-multi-kmer BOOL Skip k-mers occurring multiple times (>=2) [0]
--include-count-table BOOL Include count table [1]
--num-count-table INT Number of count table based center swapping [2]
--include-adjacency BOOL Include adjacency based center swapping [1]
--num-adjacency INT Number of adjacency based center swapping [3]
--use-parallelism BOOL Enable or disable parallel execution for group assignment and related k-mer processing steps [0]
--linclust-version INT Linclust version: 1: Linclust1, 2: Linclust2 [2]
--clust-hash BOOL Use clusthash before kmermatcher in linclust [0]
profile:
--pca Pseudo count admixture strength
--pcb Pseudo counts: Neff at half of maximum admixture (range 0.0-inf)
misc:
--rescore-mode INT Rescore diagonals with:
0: Hamming distance
1: local alignment (score only)
2: local alignment
3: global alignment
4: longest alignment fulfilling window quality criterion [0]
common:
--threads INT Number of CPU-cores used (all by default) [10]
--compressed INT Write compressed output [0]
-v INT Verbosity level: 0: quiet, 1: +errors, 2: +warnings, 3: +info [3]
--sub-mat TWIN Substitution matrix file [aa:dinuc.out,nucl:dinuc.out]
--max-seq-len INT Maximum sequence length [65535]
--db-load-mode INT Database preload mode 0: auto, 1: fread, 2: mmap, 3: mmap+touch [0]
--remove-tmp-files BOOL Delete temporary files [0]
--force-reuse BOOL Reuse tmp filse in tmp/latest folder ignoring parameters and version changes [0]
--mpi-runner STR Use MPI on compute cluster with this MPI command (e.g. "mpirun -np 42")
expert:
--filter-hits BOOL Filter hits by seq.id. and coverage [0]
--sort-results INT Sort results: 0: no sorting, 1: sort by E-value (Alignment) or seq.id. (Hamming) [0]
--filter-cludb-file STR Specify a cluDB file to filter
--filter-seqdb-file STR Specify a seqDB of cluDB file to filter
examples:
# Linear-time clustering of FASTA file
mmseqs linclust sequenceDB clusterDB tmp
--cov-mode
# Sequence 0 1 2
# Q: MAVGTACRPA 60% IGN 60%
# T: -AVGTAC--- 60% 100% IGN
# Cutoff -c 0.7 - + -
# -c 0.6 + + +
# Cluster nucleotide sequences
他のコマンドは割愛。1記事あたりの文字数制限(5万文字以内)があるため。
実行方法
QUERY.fastaに含まれる配列を、DB.fastaの配列データベースに対して検索する。
riboseek easy-search query.fasta DB.fasta output_aln.m8 tmp
easy-searchは、入力されたFASTAから内部データベースを自動作成して、そのまま検索、読みやすい表形式出力まで行う便利なコマンドとなっている。
NCBI ntのように、同じDB.fastaを何度も検索する場合は、riboseek createdbコマンドを使って事前にデータベースを作成しておく方が効率的となる。
#1.DB作成
riboseek createdb DB.fasta target_db
riboseek createindex target_db tmp_index
#2.ここではクエリもDB作成
riboseek createdb query.fasta query_db
#3.Search
riboseek search query_db target_db search_result tmp_search -a
#4.検索結果を表形式に変換、convertalisコマンドを使う
riboseek convertalis query_db target_db search_result output.m8
#5.MSAに変換。 result2msaコマンドを使う
riboseek result2msa query_db target_db search_result msa_db --msa-format-mode 6
#6.MSAファイルを取り出す、unpackdbコマンドを使う
riboseek unpackdb msa_db msa_a3m --unpack-suffix a3m --unpack-name-mode 0
このあとも、通常検索の結果から、RNAの二次構造情報を考慮した再検索を行うこともできる。この場合、MSAを展開したmsa_a3mではなく、3で作成した検索結果データベースsearch_resultを使う。
#7.riboseek cmbuildで検索結果から共分散モデル(CM)を構築
riboseek cmbuild query_db target_db search_result cm_model --cmlite-msa-eval 1e-3 --threads 4
#8.riboseek cmsearchを使い、作成したCMから候補配列を再検索
riboseek cmsearch cm_model target_db search_result cm_search_result --cm-region 3.0 --threads 4
#9.CM検索の結果からMSAを作成。result2msaコマンドを使う
riboseek result2msa query_db target_db cm_search_result cm_msa_db --msa-format-mode 6
#10.MSAファイルを取り出す、unpackdbコマンドを使う
riboseek unpackdb cm_msa_db cm_msa_a3m --unpack-suffix a3m --unpack-name-mode 0
#11.CM検索の結果を表形式に変換
riboseek convertalis query_db target_db cm_search_result output_cm.m8
テストラン
クエリの配列QUERY.fastaを、配列データベースDB.fastaに対してサーチする。ここではmmseq系でお馴染みのeasy-searchコマンドを使う。easy-searchは、FASTA配列を入力して、相同性検索と検索結果の表形式への変換までを一括で行う。
git clone https://github.com/steineggerlab/riboseek.git
cd riboseek/
riboseek easy-search query.fasta example/DB.fasta output_aln.m8 tmp
出力例
出力ファイルには、querとtargetの配列名、配列一致率、アラインメント長、ミスマッチ数、ギャップ数が記録される。続いて、queryとtarget上の一致位置、E-value、bit scoreが出力される。


(codexより) カラムについてはMMseqのマニュアルで説明されている(リンク)。
helpをみると、easyコマンドは以下のような種類が用意されている。
- easy-search:配列データの読み込みから、データベース作成、相同性検索、検索結果の表形式出力まで一括して行う。
- easy-cluster:配列データの読み込みから、感度重視の配列クラスタリング、クラスタ結果の出力まで行う。
- easy-linclust:配列データの読み込みから、高速な配列クラスタリング、クラスタ結果の出力まで行う。
- easy-proteomecluster:複数のプロテオーム配列データの読み込みから、プロテオーム間のクラスタリングと代表プロテオームの選定まで行う。
- easy-proteomesearch:プロテオーム配列データの読み込みから、プロテオーム間の類似性検索と比較結果の出力まで行う。
- easy-taxonomy:配列データの読み込みから、相同性検索、分類群の割り当て、分類結果のレポート作成まで行う。
- easy-rbh:2つの配列データセットの検索から、相互ベストヒット(RBH)の抽出まで行う。
RiboseekはMMseqs2を基盤として開発されており、MMseqs2由来の汎用コマンドも同じ実行ファイルに含まれている。easy-proteomesearchやeasy-proteomeclusterが表示さているが、Riboseekの中心機能がタンパク質向けという意味ではない。RNAを中心としたヌクレオチド配列向けのツールとなる。
その他
- ヌクレオチド相同性検索に使われるBLASTNのような高速なツールは、ヒューリスティックなワードヒットのシード化と伸長に依存しているため、その感度はシード長、伸長パラメータ、および十分に保存された局所的な一致領域が存在するかどうかに左右される(4、5)。その結果、置換、挿入・欠失、あるいは配列マスキングによってシードとして利用できる領域が少なくなると、遠縁のヌクレオチドホモログが見逃される可能性がある。
- ヌクレオチド相同性検索の別の手法として、単一のクエリ配列を配列ファミリーの統計モデルに置き換えることによって感度を高めるプロファイルベースの手法がある。現在の標準的なツールであるnhmmerは、この枠組みをDNA/RNAの相同性検索へ拡張し、遠縁のヌクレオチドホモログの検出能力を向上させた(6)。しかし、プロファイル隠れマルコフモデル(profile HMM)を用いるヌクレオチド検索は、最速のシードベース手法よりも2桁以上遅く、大規模なデータベースや多数のクエリへの拡張性が制限されている。
- 高感度側の別のアプローチは、多くのRNAファミリーでは相補的置換を通じて二次構造が保存されるため、塩基対形成相互作用を維持したまま一次配列が多様化し得るという事実に基づいている(7、8)。共分散モデルは、配列保存性と構造的共変動を単一の確率的枠組みのなかで同時にモデル化することにより、このシグナルを捉える(7)。Infernal(9)に実装されたことで、共分散モデルはRfam(10)などのリソースにおける高感度なRNAファミリー検出の標準手法となった。Infernalのフィルターパイプラインと、HMMによってバンド制限された共分散モデルアラインメントは、共分散モデル検索の速度を大幅に向上させたが、それでもなお計算負荷は非常に高い。
- RNA構造予測では、多数の相同配列を使った深いMSAが精度に関わる。現在では、RNA構造予測における深いアラインメントの作成が大きなボトルネックになっている (アラインメントの深さがモデルの精度に影響する)。
- AlphaFold 3(11)とProtenix-v2(12)は、いずれもRNA MSAを入力として受け付ける。AlphaFold 3のパイプラインでは、nhmmerを用いてRNAcentral(13)およびNT(14)を検索することによりRNA MSAを構築する。さらに低速な方法としてrMSA(15)があり、nhmmerおよびBLASTNによる検索に、Rfamに基づくフィルタリングと構造を考慮したフィルタリングを組み合わせている。大規模データベースを用いたMSA生成はクエリごとに実行され、数時間を要し主要なボトルネックとなっている。
- MMseqs2(16)のようなツールは、タンパク質相同性検索における速度と感度のトレードオフを大きく変え、AlphaFold2(17)やColabFold(18)によるタンパク質構造予測を含むさまざまな応用を可能にした。しかし、ヌクレオチド配列については、これに相当する進歩は達成されていない。本研究では、プロファイルベースのダイマー採点、高速なプレフィルタリング、およびギャップを含むアラインメントを組み合わせ、大規模に相同配列を同定する配列ベースの検索エンジンRiboseekを提示する(論文図1)。
- 本研究のベンチマーク全体において、RiboseekはBLASTN、nhmmer、およびペアワイズSmith–Waterman検索よりも高い感度を示しながら、nhmmerおよびInfernalよりも2桁高速だった。
- Riboseekで作成したMSAをProtenix-v2に入力すると、nhmmerのMSAを使った場合より、構造予測の指標が改善した。また、約173万本のRNAについて事前計算済みMSAを作成し、Web APIやColabFoldから利用できるようにした。
- (レポジトリより) Riboseekは、大規模なRNA配列セットに対する高速かつ高感度な検索を可能にする。本ツールは、事前フィルタリングにジヌクレオチド(2塩基)単位のアルファベットを用い、NTデータベース規模の塩基配列データベース検索においてGPUによる高速化をサポートしている。検索機能に加え、Riboseekは後続のRNA解析に用いるための多重配列アラインメントや共分散モデルの生成も行う(つまり、簡単に言うと、Riboseekは、RNAやDNAの大規模な配列データベースから検索対象の配列と相同な配列を高速に見つけ、多重配列アラインメント(MSA)の作成や、RNAの構造情報を考慮した再検索まで行えるツール)
- テストデータに含まれる8,791配列を対象にCM検索を実行したところ、メモリ不足で強制終了した。 (テスト環境:Win11、64 GBメモリ)。CM検索は通常の配列検索よりも多くのメモリと計算時間を必要とするため、大規模なデータベースを検索する場合は、十分な計算資源が必要。
- 普通の検索(easy-searchやsearch)では、DNA配列もFASTA形式で使用できる。ただしRiboseekはRNA向けに設計されたヌクレオチド配列検索ツールのため、CM検索では、内部的にTをUとして扱い、RNAの二次構造を考慮した検索となる。
- easy-searchは、FASTA入力から通常の相同性検索を行い、結果をm8形式で出力するところまで行う。現時点では、CM構築や構造を考慮した再検索までを一括実行するeasy-cmsearchのようなラッパーコマンドは用意されていない。
引用
Fast remote nucleotide sequence alignment with Riboseek
Sukhwan Park, Kieran Didi, Andrew Favor, Anton Bushuiev, Soohyun Kim, View ORCID ProfileMilot Mirdita, Martin Steinegger
bioRxiv, Posted July 31, 2026.
関連
参考
Fast Remote Nucleotide Sequence Alignment With Riboseek
— Biology+AI Daily (@BiologyAIDaily) August 3, 2026
1 Riboseek is a fast homology search + alignment tool for RNA and DNA aimed at the current bottleneck in RNA structure prediction pipelines: generating deep, high-quality nucleotide MSAs at scale (RNAcentral/NT-sized… pic.twitter.com/YCwPc1ufU9