vsearch(1)
β¨ NAME β¨
okay bestie, welcome to vsearch(1) β a versatile open-source tool for microbiome analysis, and when I say versatile I mean she does EVERYTHING π
Like, we're talking chimera detection, clustering, dereplication and rereplication, extraction, FASTA/FASTQ/SFF file processing, masking, orienting, pairwise alignment, restriction site cutting, searching, shuffling, sorting, subsampling, and taxonomic classification of amplicon sequences for metagenomics, genomics, and population genetics. She's the whole group chat in one binary fr fr π
π SYNOPSIS π
okay so the vibe here is: pick a mode, give her files, give her output flags, and she eats. every single one of these is a real invocation, no cap:
𧬠Chimera detection
vsearch (--uchime_denovo | --uchime2_denovo | --uchime3_denovo) fastafile (--chimeras | --nonchimeras | --uchimealns | --uchimeout) outputfile [options]
vsearch --uchime_ref fastafile (--chimeras | --nonchimeras | --uchimealns | --uchimeout) outputfile --db fastafile [options]
πΈ Clustering
vsearch (--cluster_fast | --cluster_size | --cluster_smallmem | --cluster_unoise) fastafile (--alnout | --biomout | --blast6out | --centroids | --clusters | --mothur_shared_out | --msaout | --otutabout | --profile | --samout | --uc | --userout) outputfile --id real [options]
π Dereplication and rereplication
vsearch --fastx_uniques (fastafile | fastqfile) (--fastaout | --fastqout | --tabbedout | --uc) outputfile [options]
vsearch (--derep_fulllength | --derep_id | --derep_prefix) fastafile (--output | --uc) outputfile [options]
vsearch --derep_smallmem (fastafile | fastqfile) --fastaout outputfile [options]
vsearch --rereplicate fastafile --output outputfile [options]
π¦ Extraction of sequences
vsearch --fastx_getseq fastafile (--fastaout | --fastqout | --notmatched | --notmatchedfq) outputfile --label label [options]
vsearch --fastx_getseqs fastafile (--fastaout | --fastqout | --notmatched | --notmatchedfq) outputfile (--label label --labels labelfile | --label_word label | --label_words labelfile) [options]
vsearch --fastx_getsubseq fastafile (--fastaout | --fastqout | --notmatched | --notmatchedfq) outputfile --label label [--subseq_start position] [--subseq_end position] [options]
π FASTA/FASTQ/SFF file processing
vsearch --fasta2fastq fastafile --fastqout outputfile [options]
vsearch --fastq_chars fastqfile [options]
vsearch --fastq_convert fastqfile --fastqout outputfile [options]
vsearch (--fastq_eestats | --fastq_eestats2) fastqfile --output outputfile [options]
vsearch --fastq_filter fastqfile [--reverse fastqfile] (--fastaout | --fastaout_discarded | --fastqout | --fastqout_discarded --fastaout_rev | --fastaout_discarded_rev | --fastqout_rev | --fastqout_discarded_rev) outputfile [options]
vsearch --fastq_join fastqfile --reverse fastqfile (--fastaout | --fastqout) outputfile [options]
vsearch --fastq_mergepairs fastqfile --reverse fastqfile (--fastaout | --fastqout | --fastaout_notmerged_fwd | --fastaout_notmerged_rev | --fastqout_notmerged_fwd | --fastqout_notmerged_rev | --eetabbedout) outputfile [options]
vsearch --fastq_stats fastqfile [--log logfile] [options]
vsearch --fastx_filter inputfile [--reverse inputfile] (--fastaout | --fastaout_discarded | --fastqout | --fastqout_discarded --fastaout_rev | --fastaout_discarded_rev | --fastqout_rev | --fastqout_discarded_rev) outputfile [options]
vsearch --fastx_revcomp inputfile (--fastaout | --fastqout) outputfile [options]
vsearch --sff_convert sff-file --fastqout outputfile [options]
πΈ Masking
vsearch --fastx_mask fastxfile (--fastaout | --fastqout) outputfile [options]
vsearch --maskfasta fastafile --output outputfile [options]
β Orienting
vsearch --orient fastxfile --db fastxfile (--fastaout | --fastqout | --notmatched | --tabbedout) outputfile [options]
π§ Pairwise alignment
vsearch --allpairs_global fastafile (--alnout | --blast6out | --matched | --notmatched | --samout | --uc | --userout) outputfile (--acceptall | --id real) [options]
βοΈ Restriction site cutting
vsearch --cut fastafile --cut_pattern pattern (--fastaout | --fastaout_rev | --fastaout_discarded | --fastaout_discarded_rev) outputfile [options]
π Searching
vsearch --search_exact fastafile --db fastafile (--alnout | --biomout | --blast6out | --mothur_shared_out | --otutabout | --samout | --uc | --userout | --lcaout) outputfile [options]
vsearch --usearch_global fastafile --db (fastafile | udbfile) (--alnout | --biomout | --blast6out | --mothur_shared_out | --otutabout | --samout | --uc | --userout | --lcaout) outputfile --id real [options]
π Shuffling and sorting
vsearch (--shuffle | --sortbylength | --sortbysize) fastafile --output outputfile [options]
π² Subsampling
vsearch --fastx_subsample fastafile (--fastaout | --fastqout) outputfile (--sample_pct real | --sample_size positive integer) [options]
π Taxonomic classification
vsearch --sintax fastafile --db (fastafile | udbfile) --tabbedout outputfile [--sintax_cutoff real] [options]
πΎ UDB database handling
vsearch --makeudb_usearch fastafile --output outputfile [options]
vsearch --udb2fasta udbfile --output outputfile [options]
vsearch (--udbinfo | --udbstats) udbfile [options]
π«Ά DESCRIPTION π«Ά
so here's the tea β β environmental or clinical molecular diversity studies generate large volumes of amplicons (e.g.; SSU-rRNA sequences) that need to be checked for chimeras, dereplicated, masked, sorted, searched, clustered or compared to reference sequences. The aim of vsearch is to offer a all-in-one open source tool to perform these tasks, using optimized algorithm implementations and harvesting the full potential of modern computers, thus providing fast and accurate data processing. iconique behaviour honestly π
Comparing nucleotide sequences is at the core of vsearch. To speed up comparisons, vsearch implements an extremely fast Needleman-Wunsch algorithm, making use of the Streaming SIMD Extensions (SSE2) of post-2003 x86-64 CPUs. If SSE2 instructions are not available, vsearch exits with an error message. On Power8 CPUs it will use AltiVec/VSX/VMX instructions, and on ARMv8 CPUs it will use Neon instructions. On other systems it can use the SIMD Everywhere (simde) library, if available.
Memory usage increases rapidly with sequence length: for example comparing two sequences of length 1 kb requires 8 MB of memory per thread, and comparing two 10 kb sequences requires 800 MB of memory per thread. π For comparisons involving sequences with a length product greater than 25 million (for example two sequences of length 5 kb), vsearch uses a slower alignment method described by Hirschberg (1975) and Myers and Miller (1988), with much smaller memory requirements.
π Input
vsearch accept as input fasta or fastq files containing one or several nucleotidic entries. In fasta files, each entry is made of a header and a sequence. The header is defined as the string comprised between the initial '>' symbol and the first space, tab or the end of the line, unless the --notrunclabels option is in effect, in which case the entire line is included. The header should contain printable ascii characters (33-126). The program will terminate with a fatal error if there are unprintable ascii characters. A warning will be issued if non-ascii characters (128-255) are encountered.
If the header matches the pattern '>[;]size=integer;label', the pattern '>label;size=integer;label', or the pattern '>label;size=integer[;]', vsearch will interpret integer as the number of occurrences (or abundance) of the sequence in the study. That abundance information is used or created during chimera detection, clustering, dereplication, sorting and searching. π―
The sequence is defined as a string of IUPAC symbols (ACGTURYSWKMDBHVN), starting after the end of the identifier line and ending before the next identifier line, or the file end. vsearch silently ignores ascii characters 9 to 13, and exits with an error message if ascii characters 0 to 8, 14 to 31, '.' or '-' are present. All other ascii or non-ascii characters are stripped and complained about in a warning message. π
In fastq files, each entry is made of sequence header starting with a symbol '@', a nucleotidic sequence (same rules as for fasta sequences), a quality header starting with a symbol '+' and a string of ASCII characters (offset 33 or 64), each one encoding the quality value of the corresponding position in the nucleotidic sequence.
vsearch operations are case insensitive, except when soft masking is activated. DUST masking is automatically applied during chimera detection, clustering, masking, pairwise alignment and searching. Soft masking is specified with the options '--dbmask soft' (for searching and chimera detection with a reference) or '--qmask soft' (for searching, de novo chimera detection, clustering and masking). When using soft masking, lower case letters indicate masked symbols, while upper case letters indicate regular symbols. Masked symbols are never included in the unique index words used for sequence comparisons, otherwise they are treated as normal symbols.
When comparing sequences during chimera detection, dereplication, searching and clustering, T and U are considered identical, regardless of their case. When aligning sequences, identical symbols will receive a positive match score (default +2). If two symbols are not identical, their alignment result in a negative mismatch score (default -4). Aligning a pair of symbols where at least one of them is an ambiguous symbol (BDHKMNRSVWY) will always result in a score of zero. Alignment of two identical ambiguous symbols (for example, R vs R) also receives a score of zero. When computing the amount of similarity by counting matches and mismatches after alignment, ambiguous nucleotide symbols will count as matching to other symbols if they have at least one of the nucleotides (ACGTU) they may represent in common. For example: W will match A and T, but also any of MRVHDN. When showing alignments (for example with the --alnout option) matches involving ambiguous symbols will be shown with a plus character (+) between them while exact matches between non-ambiguous symbols will be shown with a vertical bar character (|). π¦
okay bestie so vsearch isn't just a "read a file, write a file" kinda girl π she can read from standard files and write to standard files, but she can ALSO read from pipes and write to pipes! Like, multiple fasta files can be piped into vsearch for dereplication, no big deal. To do that, you just swap the file name out for:
- the symbol
-, which means/dev/stdinfor input files or/dev/stdoutfor output files (there's one exception for--db -, see the * below π), - a named pipe you made with the command
mkfifo, - a process substitution:
<(command)as input or>(command)as output. --db -is not accepted, and that's on purpose β it stops potential concurrent reads from stdin. For the advanced girlies, the workaround is to call--db /dev/stdindirectly. π
vsearch can automatically read compressed gzip or bzip2 files if the appropriate libraries are present during the compilation. She can also read pipes streaming compressed gzip or bzip2 data if the options --gzip_decompress or --bzip2_decompress are selected. Heads up tho: when reading from a pipe, the progress indicator is not updated. π
β¨ Options β¨ π
vsearch recognizes a large number of command-line commands and options. For easier navigation, options are grouped below by theme (chimera detection, clustering, dereplication and rereplication, FASTA/FASTQ file processing, masking, pairwise alignment, searching, shuffling, sorting, and subsampling). We start with the general options that apply to all themes. Options start with a double dash (--). A single dash (-) may also be used, except on NetBSD systems. Option names may be shortened as long as they are not ambiguous (e.g. --derep_f). iykyk β¨
Help and version commands π«Ά
--help,-hβ Display help text with brief information about all commands and options.--version,-vβ Output version information and a citation for the VSEARCH publication. Show the status of the support for gzip- and bzip2-compressed input files.
General options π
--bzip2_decompressβ When reading from a pipe streaming bzip2-compressed data, decompress the data. This option is not needed when reading from a standard bzip2-compressed file.--fasta_widthpositive integer β Fasta files produced by vsearch are wrapped (sequences are written on lines of integer nucleotides,80by default). Set the value to zero to eliminate the wrapping.--gzip_decompressβ When reading from a pipe streaming gzip-compressed data, decompress the data. This option is not needed when reading from a standard gzip-compressed file.--label_suffixstring β When writing FASTA or FASTQ files, add the suffix string to sequence headers.--logfilename β Write messages to the specified log file. Information written includes program version, amount of memory available, number of cores and command line options, and if need be, informational messages, warnings and fatal errors. The start and finish times are also recorded as well as the elapsed time and the maximum amount of memory consumed. The different vsearch commands can also write additional information to the log file.--maxseqlengthpositive integer β All vsearch operations discard sequences longer than integer (50,000nucleotides by default).--minseqlengthpositive integer β All vsearch operations discard sequences shorter than integer:1nucleotide by default for sorting or shuffling,32nucleotides for clustering and dereplication as well as the commands--makeudb_usearch,--sintax, and--usearch_global.--no_progressβ Do not show the gradually increasing progress indicator.--notrunclabelsβ Do not truncate sequence labels at first space or tab, but use the full header in output files. Turned off by default for all commands except the sintax command.--quietβ Suppress all messages to stdout and stderr except for warnings and fatal error messages.--samplestring β When writing FASTA or FASTQ files, add the the given sample identifier string to sequence headers. For instance, if the given string isABC, the text;sample=ABCwill be added to the header. Note that string will be truncated at the first;or blank character. Other characters (alphabetical, numerical and punctuations) are accepted.--threadspositive integer β Number of computation threads to use (1to1024). The number of threads should be less than or equal to the number of available CPU cores. The default is to use all available resources and to launch one thread per core. The following commands are multi-threaded:allpairs_global,cluster_fast,cluster_size,cluster_smallmem,cluster_unoise,fastq_mergepairs,fastx_mask,maskfasta,search_exact,sintax,uchime_ref, andusearch_global. Only one thread is used for the other commands.
π¦ Chimera detection options π¦
Chimera detection is based on a scoring function controlled by five options (--dn, --mindiffs, --mindiv, --minh, --xn). Sequences are first sorted by decreasing abundance, if available, and compared on their plus strand only (case insensitive).
Input sequences are masked as specified with the --qmask and --hardmask options. Masking of the database for reference based chimera detection is specified with the --dbmask option.
In de novo mode, input fasta file must present abundance annotations (i.e. a pattern [;]size=integer[;] in the fasta header). Input order matters for chimera detection, so we recommend to sort sequences by decreasing abundance (default of --derep_fulllength command). If your sequence set needs to be sorted, please see the --sortbysize command in the sorting section.
--abskewreal β When using--uchime_denovo, the abundance skew is used to distinguish in a three-way alignment which sequence is the chimera and which are the parents. The assumption is that chimeras appear later in the PCR amplification process and are therefore less abundant than their parents. For--uchime3_denovothe default value is16.0. For the other commands, the default value is2.0, which means that the parents should be at least 2 times more abundant than their chimera. Any positive value equal or greater than1.0can be used.--alignwidthpositive integer β When using--uchimealns, set the width of the three-way alignments (80nucleotides by default). Set to zero to eliminate wrapping.--borderlinefilename β Output borderline chimeric sequences to filename, in fasta format. Borderline chimeric sequences are sequences that have a high enough score but which are not sufficiently different from their closest parent.--chimerasfilename β Output chimeric sequences to filename, in fasta format. Output order may vary when using multiple threads.--dbfilename β When using--uchime_ref, detect chimeras using the reference sequences contained in filename. Reference sequences are assumed to be chimera-free. Chimeras cannot be detected if their parents, or sufficiently close relatives, are not present in the database. The file name must refer to a FASTA file or to a UDB file. If a UDB file is used, it should be created using the--makeudb_usearchcommand with the--dbmask dustoption.--dnstrictly positive real number β pseudo-count prior on the number of no votes, corresponding to the parameter n in the chimera scoring function (default value is1.4). Increasing--dnreduces the likelihood of tagging a sequence as a chimera (less false positives, but also more false negatives).--fasta_scoreβ Add the chimera score to the headers in the fasta output files for chimeras, non-chimeras and borderline sequences, using the format;uchime_denovo=float;.--lengthoutβ Write sequence length information to the output files in FASTA format by adding a;length=integerattribute in the header.--mindiffspositive integer β Minimum number of differences per segment (default value is3). The parameter is ignored with--uchime2_denovoand--uchime3_denovo.--mindivreal β Minimum divergence from closest parent (default value is0.8). The parameter is ignored with--uchime2_denovoand--uchime3_denovo.--minhreal β Minimum score (h). Increasing this value tends to reduce the number of false positives and to decrease sensitivity. Default value is0.28, and values ranging from0.0to1.0included are accepted. The parameter is ignored with--uchime2_denovoand--uchime3_denovo.--nonchimerasfilename β Output non-chimeric sequences to filename, in fasta format. Output order may vary when using multiple threads.--relabelstring β Relabel sequences using the prefix string and a ticker (1,2,3, etc.) to construct the new headers. Use--sizeoutto conserve the abundance annotations.--relabel_keepβ When relabelling, keep the old identifier in the header after a space.--relabel_md5β relabels every sequence by running the MD5 message digest algorithm on it. Bestie, heads up: the old headers get yeeted, gone, discarded, no notes π. Before the digest gets computed, the sequence gets uppercased and everyUgets swapped for aT. MD5 is a cryptographic hash function, and it's literally designed to make it super unlikely that two different inputs spit out the same output β even when those inputs are nearly identical but not quite. That said π there's still a teeny tiny non-zero chance two different inputs collide and give the same digest (yes, that's called a collision). MD5 drops a 128-bit (16-byte) digest, written out as 16 hexadecimal numbers (so 32 symbols picked from0123456789abcdef). Want to keep your abundance annotations? Use--sizeout. π«Ά--relabel_selfβ relabels sequences by using each sequence as its own label. The sequence is the label. Main character behaviour, honestly. β¨--relabel_sha1β same vibe as--relabel_md5but it runs the SHA1 message digest algorithm instead of MD5. SHA1 drops a 160-bit (20-byte) digest, written as 20 hexadecimal numbers (40 symbols). And the collision probability β two non-identical sequences ending up with the same digest β is smaller for SHA1 than for MD5. Upgrade moment fr fr π--selfβ when you're using--uchime_ref, this makes vsearch ignore a reference sequence whenever its label matches the label of the query sequence. Super useful for estimating the false-positive rate in your reference sequences. π--selfidβ when you're using--uchime_ref, ignore a reference sequence when its nucleotide sequence is strictly identical to the nucleotidic sequence of the query. Same energy, but it's checking the actual sequence instead of the label.--sizeinβ in de novo mode, abundance annotations (the pattern'[>;]size=integer[;]') that show up in sequence headers are taken into account by default (so--sizeinis always implied, bestie, it's just built in). This option gets ignored by--uchime_ref. π--sizeoutβ when relabelling, this adds abundance annotations to the fasta headers, using the format';size=integer;'. β¨--uchime_denovo filenameβ detects chimeras hanging out in the fasta-formatted filename, with no external references needed (that's the de novo part π ). It automatically sorts the sequences in filename by decreasing abundance first (peep the sorting section for the deets). Heads up though: multithreading is not supported. π«--uchime2_denovo filenameβ detects chimeras in the fasta-formatted filename using the UCHIME2 algorithm. This one's built for denoised amplicons (go look at--cluster_unoise). Just like the one above, it auto-sorts the sequences in filename by decreasing abundance beforehand (sorting section, you know the drill). And again: multithreading is not supported. π«--uchime3_denovo filenameβ also detects chimeras in the fasta-formatted filename using the UCHIME2 algorithm. The only difference from--uchime2_denovo? The default minimum abundance skew (--abskew) is set to16.0instead of2.0. That's it, that's the whole tea. π--uchime_ref filenameβ detects chimeras in the fasta-formatted filename by comparing them against reference sequences (that's the--dboption). And good news bestie: multithreading is supported here! π--uchimealns filenameβ writes the three-way global alignments (parentA,parentB,chimera) out to filename in a human-readable format. Wanna change the alignment length? Use--alignwidth. Note that output order can vary when you're running multiple threads, and every sequence gets converted to upper case before alignment. Lower case letters mean there's disagreement in the alignment. π--uchimeout filenameβ writes the chimera detection results to filename using an 18-field, tab-separated uchime-like format. Want the usearch v5-and-earlier compatible version instead? Use--uchimeout5. And yeah, row output order may vary when you're using multiple threads. Here's the full field breakdown, no skips:score: higher score means a more likely chimeric alignment.Q: query sequence label.A: parent A sequence label.B: parent B sequence label.T: top parent sequence label (i.e. the parent most similar to the query). That field gets removed when you use--uchimeout5.idQM: percentage of similarity of query (Q) and model (M) constructed as a part of parent A and a part of parent B.idQA: percentage of similarity of query (Q) and parent A.idQB: percentage of similarity of query (Q) and parent B.idAB: percentage of similarity of parent A and parent B.idQT: percentage of similarity of query (Q) and top parent (T).LY: yes votes in the left part of the model.LN: no votes in the left part of the model.LA: abstain votes in the left part of the model.RY: yes votes in the right part of the model.RN: no votes in the right part of the model.RA: abstain votes in the right part of the model.div: divergence, defined as(idQM - idQT).YN: query is chimeric (Y), or not (N), or is a borderline case (?).
--uchimeout5β when you're using--uchimeout, this writes the chimera detection results in a 17-field, tab-separated uchime-like format (it drops the 5th field of--uchimeout), which makes it compatible with usearch version 5 and earlier versions. π°οΈ--xlengthβ strips the header attribute";length=integer"from your input sequences. (Fun fact: that attribute is what--lengthoutadds to your output sequences.) βοΈ--xnstrictly positive real number β the weight of no votes, which corresponds to the parameterbetain the scoring function (default value is8.0). Bumping--xnup reduces the likelihood of tagging a sequence as a chimera β so fewer false positives, but also more false negatives. Pick your fighter. βοΈ--xsizeβ strips the abundance information out of the headers when writing the output file. π§Ή
π¦ Clustering options π¦
vsearch runs a single-pass, greedy centroid-based clustering algorithm, similar to the algorithms you'd find in usearch, DNAclust and sumaclust, for example. The important parameters here are the global clustering threshold (--id) and the pairwise identity definition (--iddef). Those two are the moment, lowkey. π
Input sequences get masked exactly as you specify with the --qmask and --hardmask options. π
--biomoutfilename β okay bestie, this one hands you an OTU table in the biom version 1.0 JSON file format, exactly as specified at https://biom-format.org/documentation/format_versions/biom-1.0.html π That format is literally about how to store a sparse matrix holding the abundances of the OTUs in the different samples.
And the tea? It is so much more efficient than the classic and mothur OTU table formats you get from the --otutabout and --mothur_shared_out options respectively, and it's the recommended one β at least for large tables. π
The OTUs are represented by the cluster centroids, and taxonomy info gets included for the OTUs if it's available. πΈ
Sample identifiers get extracted from the headers of all sequences in the input file. If the header contains ';sample=abc123;' or ';barcodelabel=abc123;' or a similar string somewhere, then that given sample identifier (here 'abc123') is what gets used. The semicolon is not mandatory at the beginning or end of the header. A sample identifier may contain any printable character except semicolons. And if no such sample label is found? Then the identifier in the initial part of the header gets used β but only letters, digits and underscores are allowed. π
OTU identifiers get extracted from the headers of the cluster centroid sequences. If the header contains ';otu=def789;' or a similar string somewhere, then that given OTU identifier (here 'def789') is used. Semicolon, again, not mandatory at the beginning or end of the header. An OTU identifier may contain any printable character except semicolons. If no such OTU label is found, the identifier in the initial part of the header is used, and all characters except semicolons are allowed there. Alternatively, OTU identifiers can be generated using the relabelling options (--relabel, --relabel_self, --relabel_sha1, or --relabel_md5). β¨
Taxonomy information, if present, also gets extracted from the headers of the centroid sequences. If the header contains ';tax=Homo_sapiens;' or a similar string somewhere, then that given taxonomy information (here 'Homo_sapiens') is used. Semicolon β you already know β not mandatory at the beginning or end of the header. The taxonomy information may contain any printable character except semicolons. π¦
And if what you actually need is an OTU table in the biom version 2.1 HDF5 file format? Then the biom utility may be used, as described at https://biom-format.org/documentation/biom_conversion.html. π
--centroidsfilename β write the cluster centroid sequences out to filename, in fasta format. The centroid is the sequence that seeded the cluster, i.e. the first sequence of the cluster. That's it, that's the whole vibe. π«Ά--clusterout_idβ adds cluster identifier information to your output files whenever you're using the--centroids,--consoutand--profileoptions. Little ID badge for every cluster. π·οΈ--clusterout_sortβ sorts some of the output files by decreasing abundance instead of input order. It applies to the--consout,--msaout,--profile,--centroids, and--ucoptions. And for--ucspecifically, the sorting applies only to the centroid information part (theClines). Noted. π--cluster_fastfilename β clusterize the fasta sequences in filename, and automatically sort by decreasing sequence length beforehand. It's giving "fast" for a reason. β--cluster_sizefilename β clusterize the fasta sequences in filename, and automatically sort by decreasing sequence abundance beforehand. π--cluster_smallmemfilename β clusterize the fasta sequences in filename without automatically modifying their order beforehand. Heads up: the sequences are expected to be sorted by decreasing sequence length, unless--usersortis used. No auto-sort, you're on your own, bestie. π--cluster_unoisefilename β perform denoising of the fasta sequences in filename according to the UNOISE version 3 algorithm by Robert Edgar, but without the de novo chimera removal step β that part may be performed afterwards with--uchime3_denovo. The options--minsize(default8) and--unoise_alpha(default2.0) may be specified. π₯
In this algorithm, the clustering of sequences depends on both the sequence distance and the abundance ratio. The abundance ratio (skew) is the abundance of a new sequence divided by the abundance of the centroid sequence. This skew must not be larger than beta if the sequences should be clustered together. Beta is calculated as 2 raised to the power of minus 1 minus alpha times the sequence distance. The sequence distance used is the number of mismatches in the alignment, ignoring gaps. Which means: the abundance must be exponentially lower as the distance increases from the centroid for a new sequence to be included in the cluster. Nearer sequences with higher abundances will form their own new clusters. Math girlies, this one's for you. π§
--clustersstring β output each cluster to a separate fasta file, using the prefix string and a ticker (0,1,2, etc.) to construct the path and filenames. One cluster, one file, no cap. π--consoutfilename β output cluster consensus sequences to filename. For each cluster, a center-star multiple sequence alignment is computed with the centroid as the center, using a fast algorithm β and ngl, it is not accurate when you're using low pairwise identity thresholds. A consensus sequence is constructed by taking the majority symbol (nucleotide or gap) from each column of the alignment. Columns containing a majority of gaps get skipped, except for terminal gaps. And if the--sizeinoption is specified, sequence abundances will be taken into account. π--cons_truncateβ this command is ignored. A warning is issued. That's the whole thing. Not me being dramatic, it literally does nothing. π--gapextstring β set penalties for a gap extension. See--gapopenfor a complete description of the penalty declaration system. The default is to initialize the six gap extending penalties using a penalty of2for extending internal gaps and a penalty of1for extending terminal gaps, in both query and target sequences (i.e.2I/1E). π--gapopenstring β set penalties for a gap opening. A gap opening can occur in six different contexts: in the query (Q) or in the target (T) sequence, at the left (L) or right (R) extremity of the sequence, or inside the sequence (I).
Sequence symbols (Q and T) can be combined with location symbols (L, I, and R), and numerical values, to declare penalties for all possible contexts: aQL/bQI/cQR/dTL/eTI/fTR, where abcdef are zero or positive integers, and '/' is used as a separator.
To simplify declarations, the location symbols (L, I, and R) can be combined, the symbol (E) can be used to treat both extremities (L and R) equally, and the symbols Q and T can be omitted to treat query and target sequences equally. For instance, the default is to declare a penalty of 20 for opening internal gaps and a penalty of 2 for opening terminal gaps (left or right), in both query and target sequences (i.e. 20I/2E). If only a numerical value is given, without any sequence or location symbol, then the penalty applies to all gap openings. To forbid gap-opening, an infinite penalty value can be declared with the symbol '*'. To use vsearch as a semi-global aligner, a null-penalty can be applied to the left (L) or right (R) gaps. β¨
vsearch always initializes the six gap opening penalties using the default parameters (20I/2E). The user is then free to declare only the values he/she wants to modify. The string is scanned from left to right, accepted symbols are (0123456789/LIREQT*), and later values override previous values. π
Please note that vsearch, in contrast to usearch, only allows integer gap penalties. Because the lowest gap penalties are 0.5 by default in usearch, all default scores and gap penalties in vsearch have been doubled to maintain equivalent penalties and to produce identical alignments. π«Ά
--idreal β do not add the target to the cluster if the pairwise identity with the centroid is lower than real (value ranging from0.0to1.0included). The pairwise identity is defined as the number of (matching columns) / (alignment length - terminal gaps). That definition can be modified by--iddef. π--iddef0|1|2|3|4β change the pairwise identity definition used in--id. Values accepted are:- CD-HIT definition: (matching columns) / (shortest sequence length).
- edit distance:
(matching columns) / (alignment length)β the basic math, no notes π - edit distance excluding terminal gaps β literally the same thing as
--id, bestie π«Ά - Marine Biological Lab definition β counts every gap opening (internal or terminal) as a single mismatch, whether or not that gap went and got extended. The formula is:
1.0 - [(mismatches + gap openings)/(longest sequence length)] - BLAST definition β equivalent to
--iddef 1when you're in a context of global pairwise alignment β¨
--lengthoutβ writes sequence length info into the output files in FASTA format by adding a;length=integerattribute in the header π--match integerβ the score you hand out for a match (aka identical nucleotides) in the pairwise alignment. Default value is2.--minsize positive integerβ sets the minimum abundance of sequences for denoising when you use--cluster_unoise. Default is8.--mismatch integerβ the score for a mismatch (aka different nucleotides) in the pairwise alignment. Default value is-4.--msaout filenameβ outputs a multiple sequence alignment and a consensus sequence for each cluster to filename, in fasta format. π§ Heads up though β vsearch computes center star multiple sequence alignments using a fast method, and its accuracy can drop significantly when you're working with low pairwise identity thresholds. The consensus sequence is built by grabbing the majority symbol (nucleotide or gap) from each column of the alignment. Columns that are mostly gaps get skipped, except for terminal gaps. And if you specified--sizein, sequence abundances get taken into account when computing the consensus π--mothur_shared_out filenameβ outputs an OTU table in the mothur 'shared' tab-separated plain text format, as described at (link) β¨https://www.mothur.org/wiki/Shared_fileβ© <https://www.mothur.org/wiki/Shared_file>. That format describes how a matrix holding the abundances of the OTUs across different samples gets stored. The first line starts with the strings'label','group'and'numOtus', then follows with a list of all OTU identifiers. The lines after that β one for each sample β start with the string'vsearch', then the sample identifier, then the total number of OTUs, then a list of abundances for each OTU in that sample, in the order given on the first line. OTU and sample identifiers get extracted from the FASTA headers of the sequences. The OTUs are represented by the cluster centroids. Peep the--biomoutoption for more deets β¨--otutabout filenameβ outputs an OTU table in the classic tab-separated plain text format, as a matrix holding the abundances of the OTUs across the different samples. First line starts with the string'#OTU ID', then follows a tab-separated list of all sample identifiers. The lines after that β one for each OTU β start with the OTU identifier, then follow with a tab-separated list of abundances for that OTU in each sample, in the order given on the first line. OTU and sample identifiers get extracted from the FASTA headers of the sequences (check the--sampleoption). The OTUs are represented by the cluster centroids. An extra column gets added to the right of the table if taxonomy info is available for at least one of the OTUs β it's labelled'taxonomy', and each row then holds the taxonomy info extracted for that OTU. See--biomoutfor more details πΈ--profile filenameβ outputs a sequence profile to a text file with the frequency of each nucleotide at each position in the multiple alignment, for each cluster. There's a FASTA-like header line for each cluster, followed by the profile info in a tab-separated format. The eight columns are: position (0-based), consensus nucleotide, number of As, number of Cs, number of Gs, number of Ts or Us, number of gap symbols, and finally the total number of ambiguous nucleotide symbols (B,D,H,K,M,N,R,S,Y,VorW). All numbers are integers. If--sizeinis specified, sequence abundances get taken into account π--qmask none|dust|softβ mask regions in sequences using thedustor thesoftmethods, or don't mask at all (none). β οΈ Warning bestie: when you use soft masking, clustering becomes case sensitive. The default is to mask usingdust.--qsegout filenameβ writes the aligned part of each query sequence to filename in FASTA format π--relabel stringβ relabels sequence identifiers in the output files produced by the--consout,--profileand--centroidsoptions. Go read the description of the same option under Chimera detection for the deets π--relabel_keepβ when relabelling, keeps the old identifier in the header, after a space.--relabel_md5β relabels sequence identifiers in the output files produced by the--consout,--profileand--centroidsoptions. Same drill β see the description of the same option under Chimera detection for details.--relabel_selfβ relabels sequence identifiers in the output files produced by the--consout,--profileand--centroidsoptions. Again, see the description of the same option under Chimera detection for details.--relabel_sha1β relabels sequence identifiers in the output files produced by the--consout,--profileand--centroidsoptions. Yep, same as the others β see the description of the same option under Chimera detection for details π--sizeinβ takes the abundance annotations present in the input fasta file into account (it searches for the pattern'[>;]size=integer[;]'in the sequence headers) π¦--sizeorderβ when an amplicon is close to 2 or more centroids, both of them within the distance specified with the--idoption, this resolves the ambiguity by clustering it with the centroid that has the highest abundance, not necessarily the closest one. The option only has an effect when the value specified with--maxacceptsis higher than one.--sizeorderturns on what's sometimes called abundance-based greedy clustering (AGC), as opposed to the default distance-based greedy clustering (DGC) π--sizeoutβ adds abundance annotations to the output fasta files (it adds the pattern';size=integer;'to sequence headers). If--sizeinis specified, abundance annotations get reported to the output files, and each cluster centroid receives a new abundance value equal to the total abundance of the amplicons included in the cluster (that's the--centroidsoption). If--sizeinis not specified, input abundances get set to1for amplicons, and to the number of amplicons per cluster for centroids π--strand plus|bothβ when comparing sequences with the cluster seed, check the plus strand only (the default) or check both strands π--tsegout filenameβ writes the aligned part of each target sequence to filename in FASTA format β¨--uc filenameβ outputs clustering results into filename using a tab-separated uclust-like format with10columns and3different types of entries (S,HorC). Each fasta sequence in the input file can be either a cluster centroid (S) or a hit (H) assigned to a cluster. Cluster records (C) summarize info (size, centroid label) for each cluster. In the context of clustering, the--uc_allhitsoption has no effect on the--ucoutput. Column content varies with the type of entry (S,HorC):- Record type:
S,H, orC. - Cluster number (zero-based).
- Centroid length (
S), query length (H), or cluster size (C). - Percentage of similarity with the centroid sequence (
H), or set to'*'(S,C). - Match orientation
+or-(H), or set to'*'(S,C). - Not used, always set to
'*'(S,C) or to zero (H). - Not used, always set to
'*'(S,C) or to zero (H). - Set to
'*'(S,C) or, forH, a compact representation of the pairwise alignment using the CIGAR format (Compact Idiosyncratic Gapped Alignment Report):M(match/mismatch),D(deletion) andI(insertion). The equal sign'='indicates that the query is identical to the centroid sequence (ignoring terminal gaps).
- Record type:
- Label of the query sequence (
H), or of the centroid sequence (S,C). - Label of the centroid sequence (
H), or set to'*'(S,C).
--unoise_alphareal β sets the alpha parameter for the--cluster_unoisecommand. Default is2.0. yes bestie that's a float, not an int, we do decimals here π--usersortβ when you're using--cluster_smallmem, this lets you feed sequences in any order you want, not just the decreasing-length ordering it normally expects. freedom!! π¦--xlengthβ strips the header attribute;length=integeroff your input sequences. (That's the same attribute the--lengthoutoption likes to add to output sequences, so this is basically the un-do button.)--xsizeβ yeets the abundance info out of the headers when writing the output file. bye bye size π«Ά...β okay so MOST of the searching options apply to clustering too, plus score filtering, gap penalties and masking (all defined back in the Searching section, go look π):--alnout,--blast6out,--fastapairs,--matched,--notmatched,--maxaccepts,--maxrejects,--samout,--userout,--userfields
πΈ Dereplication and rereplication options πΈ
VSEARCH can dereplicate your sequences with the commands --derep_fulllength, --derep_id, --derep_smallmem, --derep_prefix and --fastx_uniques. Heads up though: --derep_fulllength is depreciated and has been replaced by the new --fastx_uniques command, which can also handle FASTQ files on top of FASTA files. glow-up fr β¨
The --derep_fulllength, --derep_smallmem, and --fastx_uniques commands require strictly identical sequences of the same length, but they ignore upper/lower case and treat T and U as identical symbols. The --derep_id command requires both identical sequences and identical headers/labels. The --derep_prefix command groups sequences with a common prefix and does not require them to be equally long. The --derep_smallmem uses a much smaller amount of memory when dereplicating than the other files, and may be a bit slower and cannot read the input from a pipe. It takes both FASTA and FASTQ files as input but only writes FASTA output to the file specified with the --fastaout option. The --fastx_uniques command can write FASTQ output (specified with --fastqout) or FASTA output (specified with --fastaout) as well as a special tab-separated column text format (with --tabbedout). The other commands can write FASTA output to the file specified with the --output option. All dereplication commands, except --derep_smallmem, can write output to a special UCLUST-like file specified with the --uc option. The --rereplicate command can duplicate sequences in the input file according to the abundance of each input sequence. Other valid options are --fastq_ascii, --fastq_asciiout, --fastq_qmax, --fastq_qmaxout, --fastq_qmin, --fastq_qminout, --fastq_qout_max, --lengthout, --maxuniquesize, --minuniquesize, --relabel, --relabel_keep, --relabel_md5, --relabel_self, --relabel_sha1, --sizein, --sizeout, --strand, --topn, --xlength, and --xsize.
--derep_fulllengthfilename β merges strictly identical sequences found in filename. "Identical" here means same length AND same string of nucleotides (case insensitive, andTandUcount as the same thing). Peep--sizeinand--sizeoutif you wanna take abundance values into account and compute them. π«Ά This command does not support multithreading.--derep_idfilename β merges strictly identical sequences contained in filename, exactly like the--derep_fulllengthcommand does, BUT the sequence labels (identifiers) on the header line have to be identical too. same sequence AND same name, no exceptions π--derep_smallmemfilename β merges strictly identical sequences contained in filename, same as--derep_fulllength, but using much less memory. The output goes to a FASTA file specified with the--fastaoutoption. The output is written in the order the sequences first appear in the input, and not in descending abundance order like the other dereplication commands do. It can read FASTQ files but it cannot write them. It also cannot read from a pipe β it has to be a proper file, because it gets read twice. π Dereplication is performed with a 128 bit hash function and it is not verified that grouped sequences are actually identical, however the probability that two different sequences get grouped in a dataset of one billion unique sequences is approximately1e-21. Memory footprint is appr.24bytes times the number of unique sequence. Multithreading and the options--topn,--uc, or--tabbedoutare not supported.--derep_prefixfilename β merges sequences with identical prefixes contained in filename. A short sequence identical to an initial segment (prefix) of another sequence counts as a replicate of the longer sequence. If a sequence is identical to the prefix of two or more longer sequences, it gets clustered with the shortest of them. If they're equally long, it gets clustered with the most abundant one. Any remaining ties get solved using sequence headers and sequence input order. Sequence comparisons are case insensitive, andTandUare considered identical. This command does not support multithreading.--fastaoutfilename β writes the dereplicated sequences to filename, in fasta format and sorted by decreasing abundance. Identical sequences receive the header of the first sequence of their group. If--sizeoutis used, the number of occurrences (i.e. abundance) of each sequence is indicated at the end of their fasta header using the pattern';size=integer;'. This option is only valid for--fastx_uniquesand--derep_smallmem.--fastqoutfilename β writes the dereplicated sequences to filename, in fastq format and sorted by decreasing abundance. Identical sequences receive the header of the first sequence of their group. If--sizeoutis used, the number of occurrences (i.e. abundance) of each sequence is indicated at the end of their fastq header using the pattern';size=integer;'. This option is only valid for--fastx_uniques.--fastq_asciipositive integer β defines the ASCII character number used as the basis for the FASTQ quality score. Default is33, which is what the Sanger / Illumina 1.8+ FASTQ format uses (phred+33). The value64is used by the Solexa, Illumina 1.3+ and Illumina 1.5+ formats (phred+64). Only33and64are valid arguments. no other numbers, bestie, I don't make the rules π--fastq_asciioutpositive integer β when you're using--fastq_convert,--sff_convertor--fasta2fastq, this defines the ASCII character number used as the basis for the FASTQ quality score when writing FASTQ output files. Default is33. Only33and64are valid arguments.--fastq_qmaxpositive integer β specifies the maximum quality score accepted when reading FASTQ files. Default is41, which is the usual vibe for recent Sanger/Illumina 1.8+ files.--fastq_qmaxoutpositive integer β specifies the maximum quality score used when writing FASTQ files. Default is41, usual for recent Sanger/Illumina 1.8+ files. Older formats may use a maximum quality score of40.--fastq_qminpositive integer β specifies the minimum quality score accepted for FASTQ files. Default is0, usual for recent Sanger/Illumina 1.8+ files. Older formats may use scores between-5and2.--fastq_qminoutpositive integer β specifies the minimum quality score used when writing FASTQ files. Default is0, usual for recent Sanger/Illumina 1.8+ files. Older versions of the format may use scores between-5and2.--fastq_qout_maxβ for--fastx_uniques, this says that the new quality scores computed when dereplicating FASTQ files should be equal to the maximum (best) of the input quality scores for each position (corresponding to the lowest error probability). The default is to output a quality score corresponding to the average of the error probabilities for each position. so it's max vs. average, pick your fighter β--fastx_uniquesfilename β merges strictly identical sequences contained in FASTA or FASTQ file filename. Identical sequences are defined as having the same length and the same string of nucleotides (case insensitive,TandUare considered the same). Peep the options--sizeinand--sizeoutto take into account and compute abundance values. This command does not support multithreading. By default, the quality scores in FASTQ output files will correspond to the average error probability of the nucleotides in the each position. If the--fastq_qout_maxoption is given, the quality score will be the highest (best) quality score observed in each position. π--lengthoutβ writes sequence length information to the output files in FASTA and FASTQ format by adding a";length=integer"attribute in the header.--maxuniquesizepositive integer β discards sequences with a post-dereplication abundance value greater than integer.--minuniquesizepositive integer β bestie, this one yeets any sequence whose post-dereplication abundance is smaller than integer. Gone. Not invited to the party. π--outputfilename β writes the dereplicated sequences to filename, in fasta format, sorted by decreasing abundance (most popular first, obviously π). Identical sequences get the header of the first sequence of their group. And if you're using--sizeout, the number of occurrences (i.e. the abundance) of each sequence gets tacked onto the end of its fasta header using the pattern';size=integer;'. β οΈ Heads up though: this option is not allowed for--fastx_uniquesor--derep_smallmem. Not negotiable, sorry.--relabelstring β pls go read the description of the same option under Chimera detection for the deets. β¨--relabel_keepβ when relabelling, keep the old identifier in the header after a space. Preserving the lore π«Ά--relabel_md5β again, see the description of the same option under Chimera detection. Go look, iykyk.--relabel_selfβ yes bestie, same as above: see the description of the same option under Chimera detection.--relabel_sha1β you already know the vibe: check the description of the same option under Chimera detection.--rereplicatefilename β this one duplicates each sequence the number of times indicated by the abundance of each sequence in the specified file. (Option--sizeinis always implied here, no takebacksies.) The sequence labels are identical for the same sequence, unless you use--relabel,--relabel_self,--relabel_sha1or--relabel_md5to create unique labels. Output goes to the file you set with--output, in FASTA format. And no, the output file does not contain abundance information unless you specify--sizeoutβ and in that case, an abundance of1is used. One. Singular. π--sizeinβ takes into account the abundance annotations already present in the input fasta file (it searches for the pattern'[>;]size=integer[;]'in sequence headers). And FYI, this option is active by default when rereplicating. Slay, no setup needed.--sizeoutβ adds abundance annotations to the output fasta file (adds the pattern';size=integer;'to sequence headers). Here's the tea: if--sizeinis specified, each unique sequence receives a new abundance value corresponding to its total abundance (the sum of the abundances of its occurrences). But if--sizeinis not specified, input abundances are set to1, and each unique sequence gets a new abundance value corresponding to its number of occurrences in the input file. Math girlies rise up π--strand plus|bothβ when searching for strictly identical sequences, you check the plus strand only (default) or check both strands. Pick your fighter.--tabbedoutfilename β outputs clustering info to the specified tab-separated text file, with 6 columns and a row for each input sequence. Let's break it down, bestie, column by column:- Column 1 β the original label/header of the sequence.
- Column 2 β the label of the output sequence, which equals the label/header of the first sequence in each cluster, but potentially relabelled.
- Column 3 β the cluster number, starting from
0. - Column 4 β the sequence number within each cluster, starting at
0. - Column 5 β the number of sequences in the cluster.
- Column 6 β the original label/header of the first sequence in the cluster, before any potential relabelling.
π¨ This option is only valid for the --fastx_uniques command. Don't try it anywhere else, it's not gonna work.
--topnpositive integer β outputs only the top integer sequences (i.e. the most abundant ones). Only the main characters get screen time π--ucfilename β outputs full-length or prefix-dereplication results into filename using a tab-separated uclust-like format, with 10 columns and 3 different types of entries (S,HorC). Each fasta sequence in the input file can be either a cluster centroid (S) or a hit (H) assigned to a cluster. Cluster records (C) summarize information (size, centroid label) for each cluster. Also, in the context of dereplication, the option--uc_allhitshas no effect on the--ucoutput (rip to that option, it's doing nothing here π).
Column content varies with the type of entry (S, H or C):
- Record type:
S,H, orC. - Cluster number (zero-based).
- Sequence length (
S,H), or cluster size (C). - Percentage of similarity with the centroid sequence (
H), or set to'*'(S,C). - Match orientation
+or-(H), or set to'*'(S,C). - Not used, always set to
'*'(S,C) or0(H). - Not used, always set to
'*'(S,C) or0(H). - Not used, always set to
'*'. - Label of the query sequence (
H), or of the centroid sequence (S,C). - Label of the centroid sequence (
H), or set to'*'(S,C).
--xlengthβ strips the header attribute";length=integer"from input sequences. (Fun fact: this attribute is what the--lengthoutoption adds to output sequences. So it's giving β¨cleanup eraβ¨.)--xsizeβ strips abundance information from the headers when writing the output file. Bye bye size tags π¦
πΈ Extraction options πΈ
Sequences with headers matching certain criteria can be extracted from FASTA and FASTQ files using the --fastx_getseq, --fastx_getseqs and --fastx_getsubseq commands. π
The --fastx_getseq command requires the header to match a label specified with the --label option. If the --label_substr_match option is given, the label may be a substring located anywhere in the header β otherwise, the entire header must match the label. And no worries about caps, these matches are not case-sensitive. π€
The headers in the input file are truncated at the first space or tab character unless the --notrunclabels option is given. The matching sequences get written to the files specified with the --fastaout and --fastqout options, in FASTA and FASTQ format respectively. And the sequences that don't match? They get written to the files specified with the --notmatched and --notmatchedfq options, respectively. Nobody gets left out π«Ά
The --fastx_getsubseq command is similar to --fastx_getseq, but it will extract a subsequence of the matching sequences. The start position is specified with the --subseq_start option and the end position with the --subseq_end option. The positions are 1-based, meaning the first symbol of the sequence is at position 1. If the start or end position option isn't specified, the default is to start at the first position and end at the last position in the sequence. Easy β¨
The --fastx_getseqs command is similar to --fastx_getseq but allows more flexibility in specifying the label(s) to be matched. A single label may be specified using the --label option as described above. Alternatively, a file containing a list of labels to be matched may be specified with the --labels option β and that file must be a plain text file with one label on each line. The --label_word and --label_words options may be used to specify either a single word or a file containing a list of words, respectively, to be matched. Words are defined as character sequences delimited either by a character that is not alpha-numeric (A-Z, a-z, or 0-9) or by the beginning or end of the header. Word matching is case-sensitive (yes, unlike the label matching β pay attention bestie π). The --label_field option will limit the matching of words to a certain field in the header.
--fastaoutfilename β write the extracted sequences in FASTA format to the file with the given name. π--fastqoutfilename β write the extracted sequences in FASTQ format to the file with the given name. β οΈ This option is illegal if the input is in FASTA format. Can't fake the quality scores, sorry.--fastx_getseqfilename β extract sequences from the given FASTA or FASTQ file. Specify a label to match using the--labeloption. Output files are specified with the--fastaout,--fastqout,--notmatchedand--notmatchedfqoptions. π―--fastx_getseqsfilename β extract sequences from the given FASTA or FASTQ file. Specify the label or labels to match using one of the following options:--label,--labels,--label_word, or--label_words. Output files are specified with the--fastaout,--fastqout,--notmatchedand--notmatchedfqoptions. π«--fastx_getsubseqfilename β extract a certain part of some of the sequences in the given FASTA or FASTQ file. Specify labels to match using the--labeloption. Specify the subsequence range to be extracted with the--subseq_startand--subseq_endoptions. Output files are specified with the--fastaout,--fastqout,--notmatchedand--notmatchedfqoptions. βοΈ--labelstring β specify the label to match in the sequence header. Unless the--label_substr_matchoption is given, the label must match the entire header. And the comparison is not case-sensitive. π·οΈ--label_field stringβ bestie, this one is for when you're matching with--label_wordor--label_wordsand you need to tell it which field to look in π . The field name is just a string like"abc", and it has to come before the word you're matching with an equals sign (=) in between. The field also has to be fenced off by semicolons, or by the start/end of the header. So this header right here totally matches the label123in the fieldabc:"seq1;abc=123"β¨--label_substr_matchβ give it the--labelor--labelsoption and your labels are allowed to match anywhere in the header πΈ. Without this? Nope. A label has to match the ENTIRE header, no half measures, periodt.--label_word stringβ drop in a word you want to match in the sequence header. And what even is a word here? It's a string delimited by either the start or the end of the header, or by any symbol that is not a letter (A-Z,a-z) or a digit (0-9). Heads up though: the comparison is case-sensitive π--label_words filenameβ same vibe, but a whole file of words to match against the sequence headers. The plain text file must have one word on each line. Same word definition as above β delimited by the start/end of the header or by any symbol that isn't a letter (A-Z,a-z) or digit (0-9). And yes, still case-sensitive π--labels filenameβ a file of labels to match against the sequence headers, one label per line in a plain text file. Unless you also pass--label_substr_match, a label has to match the entire header. But good news: this comparison is not case-sensitive π«Ά--notmatched filenameβ write all the sequences that did not get extracted to the file with the given name, in FASTA format. The leftovers, the un-chosen ones, the ones that didn't make the cut π--notmatchedfq filenameβ same energy but in FASTQ format. > β οΈ This option is illegal if your input is in FASTA format. Don't do it, bestie, it will not go well.--subseq_end positive integerβ sets the end position when you're extracting subsequences with the--fastx_getsubseqcommand. Positions are 1-based, so sequences start at position 1, okay? If you don't specify this option, the default is to end at the end of the sequence π--subseq_start positive integerβ the starting position for subsequence extraction with--fastx_getsubseq. Again 1-based, sequences start at position 1. If you don't specify it, the default is to start at the beginning of the sequence (position 1) β¨
𧬠FASTA/FASTQ/SFF file processing options π§¬
Okay so this whole cluster is about analysing, trimming, filtering, converting, merging, joining or reverse-complementing sequences in FASTA, FASTQ or SFF files π
- Wanna analyse FASTQ files to figure out the quality encoding and the range of quality score values being used? That's
--fastq_chars, babes. - Converting between different FASTQ file variants?
--fastq_converthas you. - Statistical analysis of the quality and length of sequences in a FASTQ file? Say hi to
--fastq_stats,--fastq_eestats, and--fastq_eestats2πΈ - Trimming, filtering and converting sequences?
--fastq_filteror--fastx_filter, take it away. - SFF files to FASTQ?
--sff_convert. And--fasta2fastqwill convert a FASTA file to a FASTQ file with fake quality scores (iconic behaviour honestly). - Paired-end reads can be merged with
--fastq_mergepairsor joined with--fastq_joinπ - And
--fastx_revcompwill reverse-complement your sequences, slay. --eeoutβ when you're using--fastq_filter,--fastx_filteror--fastq_mergepairs, this includes the number of expected errors (ee) in the sequence header of FASTQ and FASTA output files π . It's a synonym of the--fastq_eeoutoption, so same thing, different name. Want that info gone from your headers? Use the--xeeoption instead.--eetabbedout filenameβ pair this with the--fastq_mergepairscommand and it writes statistics with expected errors of each merged read to the given file. It's a tab separated file with four columns, in this order:- The number of expected errors in the forward read
- The number of expected errors in the reverse read
- The number of observed errors in the forward read
- The number of observed errors in the reverse read
And just so we're clear on what "observed" means here: those are the number of differences in the overlap region of the merged sequence relative to each of the reads in the pair π
--fasta2fastq filenameβ takes the given FASTA file, slaps a fake nucleotide quality score onto the sequences, and writes them out to the FASTQ file you specified with--fastqoutβ¨. You can tweak that quality score with the--fastq_qmaxoutoption (default41), and the--fastq_asciioutoption adjusts the FASTQ output quality ASCII base character (default33).--fastaout filenameβ when you're using--fastq_filter,--fastq_mergepairsor--fastx_filter, this writes the sequences that pass the filter (or the merged sequences) to the given FASTA-formatted file πΈ--fastaout_rev filenameβ when using--fastq_filter, or--fastx_filter, write the reverse reads passing the filter to the given FASTA-formatted file π¦--fastaout_notmerged_fwd filenameβ when using--fastq_mergepairs, write forward reads that did not get merged to the specified FASTA file π--fastaout_notmerged_rev filenameβ when using--fastq_mergepairs, write reverse reads that did not get merged to the specified FASTA file π--fastaout_discarded filenameβ write sequences that do not pass the filter of the--fastq_filteror--fastx_filtercommand to the given FASTA-formatted file ποΈ--fastaout_discarded_rev filenameβ write reverse reads that do not pass the filter of the--fastq_filteror--fastx_filtercommand to the given FASTA-formatted file ποΈ--fastq_allowmergestaggerβ when using--fastq_mergepairs, this lets you merge staggered read pairs. And what's a staggered pair, you ask? It's a pair where the 3' end of the reverse read has an overhang to the left of the 5' end of the forward read. This can happen when a very short fragment gets sequenced. Important detail: the 3' overhang of the reverse read is not included in the merged sequence π. The opposite option is--fastq_nostagger, and the default is to discard staggered pairs entirely.--fastq_ascii positive integerβ defines the ASCII character number used as the basis for the FASTQ quality score. Default is33, which is what the Sanger / Illumina 1.8+ FASTQ format uses (phred+33). The value64is used by the Solexa, Illumina 1.3+ and Illumina 1.5+ formats (phred+64). And that's it β only33and64are valid arguments, no improvising π₯--fastq_asciiout positive integerβ when using--fastq_convert,--sff_convertor--fasta2fastq, this defines the ASCII character number used as the basis for the FASTQ quality score when writing FASTQ output files. Default is33. And same rule as above: only33and64are valid arguments π--fastq_chars filenameβ summarizes the composition of the sequence and quality strings in your input FASTQ file. For every sequence symbol,--fastq_charsgives you the number of occurrences of that symbol, its relative frequency, and the length of the longest run of that symbol π. For every character present in the quality strings, it gives you the ASCII value of the character, its relative frequency, and the number of times a k-mer of that character appears at the end of quality strings. The k-mer length is settable with--fastq_tail(default4). Oh and also β--fastq_charstries to auto-detect the quality encoding (Solexa, Illumina 1.3+, Illumina 1.5+ or Illumina 1.8+/Sanger) by analysing the range of observed quality score values π. If it succeeds, it suggests values for the--fastq_ascii(33or64),--fastq_qminand--fastq_qmaxoptions to use with the other commands that need a FASTQ input file. Literally so helpful.--fastq_convert filenameβ converts between the different variants of the FASTQ file format π . You must specify the quality encoding of the input file with the--fastq_asciioption (either33or64, default33), and the output quality encoding with the--fastq_asciioutoption (default33). You can limit the minimum and maximum output quality scores with the--fastq_qminoutand--fastq_qmaxoutoptions. The output file gets specified with the--fastqoutoption β¨--fastq_eeoutβ when using--fastq_filter,--fastx_filteror--fastq_mergepairs, this includes the number of expected errors (ee) in the sequence header of FASTQ and FASTA files π. It's a synonym of the--eeoutoption (yes, they're twins). Use the--xeeoption to remove this information from headers.--fastq_eestats filenameβ okay bestie, this one is a whole quality-score detective πβ¨. It analyzes a FASTQ file and reports statistics on the distributions of quality scores, error probabilities and expected accumulated errors. The report is a table of 21 tab-separated columns and it gets written to the file you named with the--outputoption. Column breakdown, no cap:- Column 1 = the position in the reads (
Pos). - Columns 2 and 3 = the number of reads (
Reads) and the percentage of reads (PctRecs) that include this position. - The remaining columns hold info about the distribution of quality scores at this position (
Q), the error probabilities at this position (Pe), and finally the expected number of accumulated errors from the beginning of the reads up until the current position (EE).
- Column 1 = the position in the reads (
For each of the Q, Pe and EE distributions you get all of these stats: minimum value (Min), lower quartile (Low), median (Med), mean (Mean), upper quartile (Hi), and maximum value (Max). And the quality encoding plus the range of quality values? You can specify those with --fastq_ascii, --fastq_qmin and --fastq_qmax. π
--fastq_eestats2 filenameβ the sequel π. It analyzes the FASTQ file you point at and reports stats on how many sequences would be retained at a combination of selected cutoffs for length truncation and maximum expected errors β the exact kinds of values you could potentially feed as arguments to the--fastq_trunclenand--fastq_maxeeoptions of the--fastq_filtercommand. The result, a table of two or more columns, goes to the file specified with--output. There's one line for each length truncation cutoff: column 1 on each line = the selected truncation length, and the following columns = the number of sequences and, in parenthesis, the percentage of sequences that would be retained at the selected EE levels.
The truncation length cutoffs are specified with --length_cutoffs, which requires a list of three comma-separated integers: the shortest cutoff, the longest cutoff, and the increment between cutoffs. The longest cutoff may be specified with a star (*), which means the limit is equal to the longest sequence in the input file. Default setting is "50,*,50", meaning truncation lengths of 50, 100, 150 and so on up to the longest sequence length should be used. π¦
The maximum expected error (EE) cutoffs are specified with --ee_cutoffs, which requires a comma-separated list of floating point numbers as its argument. Default setting is "0.5,1.0,2.0", which indicates that expected error levels of 0.5, 1.0 and 2.0 should be used.
--fastq_filter filenameβ trim and/or filter the sequences in the given FASTQ file. It's giving--fastx_filterenergy, but it works only on FASTQ files. See--fastx_filterfor the deets. β¨--fastq_join filenameβ join paired-end sequence reads into one sequence and add a gap between them using a padding sequence π«Ά. Important: the sequences are not merged like with thefastq_mergepairscommand β they're simply joined with a gap. The forward reads are the argument to this option and the reverse reads come from the--reverseoption. The resulting sequences = the forward read, the padding sequence, and the reverse complement of the reverse read. The padding sequence is set with--join_padgapand the padding quality is set with--join_padgapq. Default padding sequence string isNNNNNNNNand the default padding quality string isIIIIIIII, corresponding to a base quality score of 40 (a very high quality score with error probability 0.0001, iconic honestly π). The joined sequences are output to the file(s) specified with the--fastaoutor--fastqoutoptions.--fastq_maxdiffs positive integerβ when you're using--fastq_mergepairs, this sets the maximum number of non-matching nucleotides allowed in the overlap region. Heads up: this option has a strong influence on the merging success rate. Default value is10. π--fastq_maxdiffpct realβ when using--fastq_mergepairs, specify the maximum percentage of non-matching nucleotides allowed in the overlap region. Default value is100.0%. There are other, more sophisticated rules in the merging algorithm that will discard read pairs with a high fraction of mismatches, so don't think this is the whole story bestie π.--fastq_maxee realβ when using--fastq_filter,--fastq_mergepairsor--fastx_filter, this discards sequences with an expected error greater than the number you give it (value ranging from0.0to infinity). For a given sequence, the expected error is the sum of error probabilities for all the positions in the sequence. Since error probabilities can be small but not null, the expected error is always greater than zero, and at most equal to the length of the sequence when all positions in the sequence have an error probability of1.0.
Using the expected error as the lambda parameter in the Poisson distribution, it's possible to compute the probability of observing k errors. For instance, a read with an expected error of 1.0 has:
- 36.8% chance of having zero error,
- 36.8% chance of having one error,
- 18.4% chance of having two errors,
- 6.1% chance of having three errors,
- 1.5% chance of having four errors,
- 0.3% chance of having five errors,
- etc.
--fastq_maxee_rate realβ when using--fastq_filteror--fastx_filter, discard sequences with an average expected error greater than the specified number (value ranging from0.0to1.0included). For a given sequence, the average expected error is the sum of error probabilities for all the positions in the sequence, divided by the length of the sequence. Math girlies rise πΈ--fastq_maxlen positive integerβ when using--fastq_filter,--fastq_mergepairsor--fastx_filter, discard sequences with more than the specified number of bases. Too long? Not the vibe. π«--fastq_maxmergelen positive integerβ when using--fastq_mergepairs, specify the maximum length of the merged sequence (default is1,000,000).--fastq_maxns positive integerβ when using--fastq_filter,--fastq_mergepairsor--fastx_filter, discard sequences with more than the specified number of N's.--fastq_mergepairs filenameβ the main character of paired-end merging πβ¨. Merge paired-end sequence reads into one sequence. The forward reads are the argument to this option, and the reverse reads come from the--reverseoption. Reads with the same index/position in the forward and reverse files are considered to form a pair, even if their labels are different. So forward and reverse reads must appear in the same order and total number in both files. A warning is emitted if the forward and reverse files contain different numbers of reads.
The merged sequences are written to the file(s) specified with the --fastaout or --fastqout options. The non-merged reads can be output to the files specified with the --fastaout_notmerged_fwd, --fastaout_notmerged_rev, --fastqout_notmerged_fwd and --fastqout_notmerged_rev options. Statistics may be output to the file specified with the --eetabbedout option.
Sequences are truncated as specified with the --fastq_truncqual option to remove low-quality bases in the 3' end. Sequences shorter than specified with --fastq_minlen (after truncation) are discarded (1 by default). Sequences with too many ambiguous bases (N's), as specified with --fastq_maxns, are also discarded (no limit by default). Staggered reads are not merged unless the --fastq_allowmergestagger option is specified.
The minimum length of the overlap region between the reads may be specified with the --fastq_minovlen option (at least 5, default 10). The overlap region may not include more mismatches than specified with the --fastq_maxdiffs option (10 by default) or a higher percentage of mismatches than specified with the --fastq_maxdiffpct option (100.0% by default), otherwise the read pair is discarded. Additional rules will avoid merging of reads that cannot be aligned reliably and unambiguously. π§
The minimum and maximum length of the merged sequence may be specified with the --fastq_minmergelen and --fastq_maxmergelen options, respectively. The quality value limits for output files may be specified with the --fastq_qminout and --fastq_qmaxout options, but they apply only to the merged region. Other relevant options are: --fastq_ascii, --fastq_maxee, --fastq_nostagger, --fastq_qmax, --fastq_qmin, and --label_suffix. π
--fastq_minlenpositive integer β when you're running--fastq_filter,--fastq_mergepairsor--fastx_filter, any input sequence shorter than this many bases gets yeeted (default1). bye bestie π--fastq_minmergelenpositive integer β for--fastq_mergepairsonly: the minimum length the merged sequence is allowed to be. default is1. π--fastq_minovlenpositive integer β for--fastq_mergepairs: how much overlap the merged reads gotta have, minimum. default is10. and it must be at least 5, no negotiating β¨--fastq_minqualpositive integer β with--fastq_filteror--fastx_filter, any read that has even ONE base with a quality score under this value gets discarded. default is0, which means nothing gets tossed. zero drama πΈ--fastq_nostaggerβ with--fastq_mergepairs, this forbids merging staggered read pairs. fun fact: this is already the default behaviour of--fastq_mergepairs. if you wanna flip that, go look at the--fastq_allowmergestaggeroption π--fastq_qmaxpositive integer β the max quality score accepted when reading FASTQ files. default is41, which is the usual vibe for recent Sanger/Illumina 1.8+ files. π§--fastq_qmaxoutpositive integer β when using--fastq_mergepairs,--fastq_convert,--sff_convertor--fasta2fastq, this sets the max quality score used when writing FASTQ files. for the--fasta2fastqcommand specifically, the value you put here is literally the fake quality score used for the FASTQ output file. default is41β usual for recent Sanger/Illumina 1.8+ files. older formats might use a max quality score of40. and heads up: when you're using--fastq_mergepairs, this limit only applies to the merged region π--fastq_qminpositive integer β the min quality score accepted for FASTQ files. default is0, usual for recent Sanger/Illumina 1.8+ files. older formats may use scores between-5and2. β--fastq_qminoutpositive integer β with--fastq_mergepairs,--fastq_convertor--sff_convert, this is the min quality score used when writing FASTQ files. default is0, usual for recent Sanger/Illumina 1.8+ files. older versions of the format may use scores between-5and2. and same deal: the limit applies only to the merged region when you're using--fastq_mergepairsπ¦--fastq_statsfilename β analyzes a FASTQ file and reports how many reads are in it. you can specify the quality encoding and the range of quality values with--fastq_ascii,--fastq_qminand--fastq_qmax. this command requires the--logoption (not optional, bestie) and it outputs the following detailed stats on read length, quality score, length vs. quality distributions, and length / quality filtering:
Read length distribution:
L: read length.N: number of reads.Pct: fraction of reads with this length.AccPct: fraction of reads with this length or longer.
Quality score distribution:
ASCII: character encoding the quality score.Q: Phred quality score.Pe: probability of error associated with the quality score.N: number of bases with this quality score.Pct: fraction of bases with this quality score.AccPct: fraction of bases with this quality score or higher.
Length vs. quality distribution:
L: position in reads (starting from position2).PctRecs: fraction of reads with at least this length.AvgQ: average quality score at this position.P(AvgQ): error probability corresponding toAvgQ.AvgP: average error probability at this position.AvgEE: average expected error over all reads up to this position.Rate: growth rate ofAvgEEbetween this position and position- 1.RatePct:Rate(as explained above) expressed as a percentage.
Effect of expected error and length filtering: the first column is read lengths (L). the next four columns are the number of reads that would be retained by the --fastq_filter command if the reads were truncated at length L (option --fastq_trunclen L) and filtered to have a maximum expected error of 1.0, 0.5, 0.25 or 0.1 (with the option --fastq_maxee float). the last four columns are the fraction of reads that would be retained by --fastq_filter using those same length and maximum expected error parameters. π₯
Effect of minimum quality and length filtering: the first column is read lengths (Len). the next four columns are the fraction of reads that would be retained by the --fastq_filter command if the reads were truncated at length Len (option --fastq_trunclen Len) or at the first position with a quality Q equal to or lesser than 5, 10, 15 or 20 (option --fastq_truncqual Q).
--fastq_stripleftpositive integer β with--fastq_filteror--fastx_filter, strips this many bases off the left end of the reads. if the resulting read ends up null-length, the read gets discarded. no survivors π--fastq_striprightpositive integer β with--fastq_filteror--fastx_filter, strips this many bases off the right end of the reads. and again: if the resulting read is null-length, it's discarded. βοΈ--fastq_tailpositive integer β with--fastq_chars, counts how many times a series of characters of length k shows up at the end of quality strings. by default, k =4. π--fastq_trunceereal β with--fastq_filteror--fastx_filter, truncates sequences so their total expected error isn't higher than the value you give. π--fastq_truncee_ratereal β with--fastq_filteror--fastx_filter, truncates sequences so their average expected error per base isn't higher than the specified value. the truncation happens at the first occurrence. and fyi the average expected error per base is calculated as the total expected number of errors divided by the length of the sequence after truncation. math girlies rise π--fastq_trunclenpositive integer β with--fastq_filteror--fastx_filter, truncates sequences to the specified length. shorter sequences get discarded. sorry babe π--fastq_trunclen_keeppositive integer β with--fastq_filteror--fastx_filter, truncates sequences to the specified length β but this time shorter sequences are not discarded. keeping them, bestie π«Ά--fastq_truncqualpositive integer β with--fastq_filter,--fastq_mergepairsor--fastx_filter, truncates sequences starting from the first base whose base quality score is this value or lower. π--fastqoutfilename β with--fastq_filter,--fastq_mergepairs,--fastx_filteror--fasta2fastq, writes the sequences that passed the filter β or the merged or converted sequences β to the given FASTQ-formatted file. β¨--fastqout_revfilename β with--fastq_filteror--fastx_filter, writes the reverse reads that passed the filter to the given FASTQ-formatted file. π--fastqout_discardedfilename β with--fastq_filteror--fastx_filter, writes the sequences that did not pass the filter to the given FASTQ-formatted file. π--fastqout_discarded_revfilename β with--fastq_filteror--fastx_filter, writes the reverse reads that did not pass the filter to the given FASTQ-formatted file. π--fastqout_notmerged_fwdfilename β with--fastq_mergepairs, writes the forward reads that did not get merged to the specified FASTQ file. π¦--fastqout_notmerged_revfilename β with--fastq_mergepairs, writes the reverse reads that did not get merged to the specified FASTQ file. πΈ--fastx_filterfilename β bestie, this is the full glow-up station for your sequences π It trims and/or filters the sequences in the given FASTA or FASTQ file, then writes whatever survives to the FASTQ file you picked with--fastqoutand/or the FASTA file you picked with--fastaout. The sequences that got voted off the island go to the files named by--fastaout_discardedand--fastqout_discarded. The input format (FASTA or FASTQ) is automatically detected, no cap. If your input is paired, you can pass the reverse reads in with the--reverseoption, and the matching output lands in--fastqout_rev,--fastaout_rev,--fastqout_discarded_rev, and--fastaout_discarded_rev. Heads up though: output can not be written to FASTQ files if the input is in FASTA format β no quality scores, no FASTQ, periodt.
The sequences are trimmed first and then filtered based on the bases that are left. Trimming happens via --fastq_stripleft, --fastq_stripright, --fastq_truncee, --fastq_truncee_rate, --fastq_trunclen, --fastq_trunclen_keep and --fastq_truncqual. Filtering happens via --fastq_maxee, --fastq_maxee_rate, --fastq_maxlen, --fastq_maxns, --fastq_minlen (default 1), --fastq_minqual, --fastq_trunclen, --maxsize, and --minsize. Sequences that don't satisfy the requirements are discarded. For pairs of sequences, both sequences in a pair have to satisfy the requirements, otherwise both get discarded β ride or die, fr fr. If you give no shortening or filtering options at all, every sequence gets written to the output files, possibly after conversion from FASTQ to FASTA format. The --relabel option can relabel your output sequences, and --eeout can output the expected number of errors in each sequence.
Once all sequences have been processed, you'll see the number of kept and discarded sequences, plus how many of the kept ones were trimmed β¨
When the input is in FASTA format, these options are not accepted because quality scores literally don't exist here:--eeout,--fastq_ascii,--fastq_eeout,--fastq_maxee,--fastq_maxee_rate,--fastq_minqual,--fastq_out,--fastq_qmax,--fastq_qmin,--fastq_truncee,--fastq_truncee_rate,--fastq_truncqual,--fastqout_discarded,--fastqout_discarded_rev,--fastqout_rev.
--fastx_revcompfilename β reverse-complements the sequences in the given FASTA or FASTQ file into a file you name with--fastaoutand/or--fastqout. Mirror flip moment πͺ If the input file is in FASTA format, the output can not be written back to a FASTQ file, because the base quality scores just aren't there.--join_padgapstring β when you're running--fastq_join, this string is used as the sequence padding string. Default isNNNNNNNN(8 N's). The option accepts an empty string, or literally any other combination of ASCII characters.--join_padgapqstring β when you're running--fastq_join, this string is the quality padding string. Default is a string ofI's equal in length to the sequence padding string. The letterIcorresponds to a base quality score of 40, which means a very high quality base with an error probability of0.0001. The option accepts an empty string, or any other combination of ASCII characters.--lengthoutβ writes sequence length info to the output files in FASTA or FASTQ format by adding a;length=integerattribute in the header π--maxsizepositive integer β when using--fastq_filteror--fastx_filter, discard sequences with an abundance higher than the specified value.--minsizepositive integer β when using--fastq_filteror--fastx_filter, discard sequences with an abundance lower than the specified value.--outputfilename β when using--fastq_eestatsor--fastq_eestats2, write tabulated results to filename. Go read--fastq_eestats's and--fastq_eestats2's docs for the full tea on what that table contains.--relabel_keepβ when using--relabel, keep the old identifier in the header after a space.--relabelstring β please see the description of the same option under Chimera detection for details.--relabel_md5β please see the description of the same option under Chimera detection for details.--relabel_selfβ please see the description of the same option under Chimera detection for details.--relabel_sha1β please see the description of the same option under Chimera detection for details.--reversefilename β when using--fastq_filter,--fastx_filter,--fastq_mergepairsor--fastq_join, specify the FASTQ file containing the reverse reads.--sff_convertfilename β converts the given SFF file to FASTQ. The FASTQ output file is named with the--fastqoutoption. The sequence may be clipped as specified in the SFF file if you also pass--sff_clip; otherwise no clipping occurs. Bases that would have been clipped get converted to lower case, while everything else stays upper case. The output quality encoding can be set with--fastq_asciiout(default33). The minimum and maximum output quality scores can be limited using--fastq_qminoutand--fastq_qmaxout.--sff_clipβ says that the sequences converted by the--sff_convertcommand should be clipped at both ends as indicated in the SFF file. By default no clipping is performed.--xlengthβ strips the header attribute;length=integerfrom input sequences. (That attribute is what--lengthoutadds to output sequences.)--xsizeβ strips abundance information from the headers when writing the output file.--xeeβ strips information about expected errors (ee) from the output file headers. (That info is added by the--fastq_eeoutand--eeoutoptions.)
π Masking options π
An input sequence can be made of lower- or uppercase letters. When soft masking is specified, lower case letters are treated as symbols that should be masked. Otherwise, the case of the input sequences is simply ignored. iykyk π«Ά
DUST masking is performed automatically by the commands for chimera detection (uchime_denovo, uchime_ref), clustering (cluster_fast, cluster_smallmem, cluster_size), masking (maskfasta, fastx_mask), pairwise alignment (allpairs_global) and searching (search_exact, usearch_global).
Masking is usually specified with the --qmask option, while the --dbmask option is used for the database sequences specified with the --db option alongside the --usearch_global, --search_exact and --uchime_ref commands.
The argument to the --qmask and --dbmask options may be none, soft or dust:
noneβ no masking is performed.softβ the lower case symbols are masked.dustβ the sequence is masked using the DUST algorithm by Tatusov and Lipman to mask low-complexity regions.
If the --hardmask option is specified, all masked regions are converted to N's, otherwise masked regions are indicated by lower case letters.
If any sequence is masked, the masked version of the sequence (with lower case letters or N's) is used in all output files. Otherwise the sequence is unmodified. The exception is the sequences in the output file specified with the --uchimealns option, where the input sequences are converted to upper case first and lower case letters indicate disagreement between the aligned sequences.
The --qmask option (or --dbmask for database sequences) may be combined with the --hardmask option. The results of using the none, dust or soft argument to --qmask or --dbmask are presented below, assuming each input sequence contains both lower and uppercase symbols π
Results if the --hardmask option is off (default):
none: no masking, all symbols used, no change
dust: masked symbols lowercased, rest uppercased
soft: lowercase symbols masked, no case changes
Results if the --hardmask option is on:
none: no masking, all symbols used, no change
dust: masked symbols changed to Ns, rest unchanged
soft: lowercase symbols masked and changed to Ns
When a sequence region is masked, words in the region are not included in the indices used in the heuristic search algorithm. In all other aspects, the region is treated like any other region.
Regions in sequences that are hardmasked (with N's) have a zero alignment score and do not contribute to an alignment. π
--fastaoutfilename β write the masked sequences to filename, in fasta format. Applies only to the--fastx_maskcommand.--fastqoutfilename β write the masked sequences to filename, in fastq format. Applies only to the--fastx_maskcommand.--fastx_mask filenameβ mask regions in the sequences chilling inside the fasta or fastq file you point at. Default vibe is masking with DUST (peep--qmaskif you wanna change that behaviour). Output files get named with--fastaoutand--fastqout. And the min/max percentage of unmasked residues? Those are your--min_unmasked_pctand--max_unmasked_pctoptions respectively. π--hardmaskβ symbols in the masked regions get swapped out for N's. The default is literally just replacing the masked regions with lower case letters. That's it, that's the whole tea. π«--maskfasta filenameβ mask regions in the sequences in the fasta filefilename. Default masking is dust (hit up--qmaskto switch it up). The output file gets specified with--output. > β οΈ bestie this command is depreciated, please use--fastx_maskinstead, no cap.--max_unmasked_pct realβ discard sequences with MORE than the specified maximum percentage of unmasked residues. Only works with--fastx_mask, okurrr. π--min_unmasked_pct realβ discard sequences with LESS than the specified minimum percentage of unmasked residues. Also only works with--fastx_mask. πΈ--output filenameβ write the masked sequences tofilename, in fasta format. Only applies to the--mask_fastacommand, fyi.--qmask none|dust|softβ okay so, if the argument isdust, you mask regions in sequences using the DUST algorithm, which detects simple repeats and low-complexity regions. This is the default for chimera detection, clustering, masking, pairwise alignment, and searching. If the argument issoft, you mask the lower case letters in the input sequence. And if the argument isnone, you don't mask at all. Iconic trio honestly. β¨
π§ Orienting options π§
The --orient command can be used to orient the sequences in a given file in either the forward or the reverse complementary direction, based on a reference database you specify with the --db option. The two strands of each input sequence get compared to the reference database using nucleotide words. If one of the strands shares many more words with at least one sequence in the database than the other one does, then that strand gets chosen. π«Ά
The correctly oriented sequences may be written to a FASTA file specified with --fastaout, and to a FASTQ file specified with the --fastqout option (as long as the input was also in FASTQ format, obviously). If the result is uncertain β because the number of matching words is too similar β the original sequence gets written to the file specified with the --notmatched option. You can also write the results to a tab-delimited text file specified with --tabbedout. That file will contain the query label, the direction (+, - or ?), the number of matching words on the forward strand, and the number of matching words on the reverse complementary strand.
By default, a word length of 12 is used for this command, and you can adjust that with the --wordlength option. There has to be at least 4 times as many matches on one strand than the other for a strand to be selected. In addition to the common options, the following options may also be specified for this command: --dbmask, --qmask, --relabel, --relabel_keep, --relabel_md5, --relabel_self, --relabel_sha1, --sizein, and --sizeout. π
--db filenameβ read the reference database from the given file. It may be in FASTA, FASTQ or UDB format. If a UDB file is used it should have been created with a wordlength of 12. π--fastaout filenameβ write the correctly oriented sequences tofilename, in fasta format.--fastqout filenameβ write the correctly oriented sequences tofilename, in fastq format.--notmatched filenameβ write the sequences with undetermined direction tofilename, in the original format.--orient filenameβ orient the sequences in the given file. Simple as that. π--tabbedout filenameβ write the resuls to a tab-delimited text file with the specifiedfilename. This file will contain the query label, the direction (+,-or?), the number of matching words on the forward strand, and the number of matching words on the reverse complementary strand.
𧬠Pairwise alignment options π§¬
The results of the n * (n-1) / 2 pairwise alignments get written to the result files specified with --alnout, --blast6out, --fastapairs, --matched, --notmatched, --qsegout, --samout, --tsegout, --uc or --userout (peep the Searching section below for those). Specify either the --acceptall option to output ALL pairwise alignments, or specify an identity level with --id to discard the weak ones. Most other accept/reject options (again, see Searching options below) may also be used. Sequences are aligned on their plus strand only. Masking is performed as usual and specified with --qmask and --hardmask. π
--acceptallβ write the results of all alignments to output files. This option overrides all other accept/reject options (including--id). It's giving main character energy, it does NOT share the spotlight. π--allpairs_global filenameβ perform optimal global pairwise alignments of the fasta sequences contained infilename. Each sequence is compared to all sequencs that come after it in the file, resulting in a total of n * (n-1) / 2 pairwise alignments, where n is the total number of sequences. This command is multi-threaded. π₯--id realβ reject the sequence match if the pairwise identity is lower thanreal(value ranging from 0.0 to 1.0 included).--threads positive integerβ number of computation threads to use (1 to 1024). The number of threads should be lesser or equal to the number of available CPU cores. The default is to use all available resources and to launch one thread per logical core. Slay, it just goes feral by default. π--uc filenameβ output pairwise alignment results infilenameusing a tab-separated uclust-like format with 10 columns. Each sequence is compared to all other sequences, and all hits (--acceptall) or only some hits (--id float) are reported, with one pairwise comparison per line:- Record type, always set to
'H'. - Ordinal number of the target sequence (based on input order, starting from zero).
- Sequence length.
- Percentage of similarity with the target sequence.
- Match orientation, always set to
'+'. - Not used, always set to zero.
- Not used, always set to zero.
- Compact representation of the pairwise alignment using the CIGAR format (Compact Idiosyncratic Gapped Alignment Report):
M(match/mismatch),D(deletion) andI(insertion). The equal sign'='indicates that the query is identical to the centroid sequence (ignoring terminal gaps). - Label of the query sequence.
- Label of the target sequence.
- Record type, always set to
βοΈ Restriction site cutting options βοΈ
The input sequences in the file specified with the --cut command get cut into fragments at all restriction sites matching the pattern given with the --cut_pattern option. The fragments on the forward strand are written to the file specified with the --fastaout file, and the fragments on the reverse strand are written to the file specified with the --fastaout_rev option. Input sequences that do not match are written to the file specified with the option --fastaout_discarded, and their reverse complement is also written to the file specified with the --fastaout_discarded_rev option. The relabel options (--relabel, --relabel_self, --relabel_keep, --relabel_md5, and --relabel_sha1) may be used to relabel the output sequences. π§
--cut filenameβ specify the input file with sequences in FASTA format.--cut_pattern stringβ specify the restriction site cutting pattern and positions. The pattern is a string of lower- or uppercase letters specifying the nucleotides that must match, and may include ambiguous nucleotide symbols. The special characters"^"(circumflex) and"_"(underscore) are used to indicate the cutting position on the forward and reverse strand, respectively. For example, the pattern"G^AATT_C"is the pattern for the EcoRI restriction site. For such palindromic patterns (identical to its reverse complement) the command will output all possible fragments on both strands. For non-palindromic sites, it may be necessary to run the command also on the reverse complemented input sequences. Exactly one cutting site on each strand must be indicated. No cap, precision matters bestie. π--fastaout filenameβ specify the output file for the resulting fragments on the forward strand.--fastaout_rev filenameβ specify the output file for the resulting fragments on the reverse strand.--fastaout_discarded filenameβ specify the output file for the non-matching sequences.--fastaout_discarded_rev filenameβ specify the output file for the non-matching sequences, reverse complemented.
π Searching options π
Okay bestie, this is the section where vsearch decides how it spills the tea about what it found. Every flag here either writes your results somewhere or tweaks how the actual search behaves. There's a lot, so let's go one by one. π
--alnoutfilename β writes your pairwise global alignments to filename in a format that actual humans can read. Want to change how long those alignment lines are? That's what--rowlenis for. Heads up though: if you're running with multiple threads, the order things come out in might be all over the place. π--biomoutfilename β dumps your search results into an OTU table using the biom version 1.0 file format. The query file is what holds the samples, and the database file is what holds the OTUs. Both the sample IDs and the OTU IDs get pulled straight out of the sequence headers. Want the full lore? Go peek at the--biomoutoption over in the Clustering section. β¨--blast6outfilename β writes results to filename in a blast-ish, tab-separated, twelve-field format (the whole field list is right below), with one line for every query-target match β or for every non-match if you've got--output_no_hitsturned on.
β οΈ Warning bestie: vsearch uses global pairwise alignments, not blast's seed-and-extend algorithm. So some of those classic blast output values (alignment start and end, evalue, bit score) get reported differently here. Not a bug, just different.
Also, multi-threaded runs may shuffle the output order around. You can get something pretty similar with --userout filename and --userfields query+target+id+alnlen+mism+opens+qlo+qhi+tlo+thi+evalue+bits. The complete list and descriptions live in the 'Userfields' section of this manual. π
The twelve fields, in order:
query: the query label.target: the target (database sequence) label. This one gets set to*when there's no alignment at all.id: percentage of identity, a real value from 0.0 to 100.0. It's defined as 100 * (matching columns) / (alignment length - terminal gaps). Want other definitions? Check out fields id0 through id4.alnlen: length of the query-target alignment (the number of columns). Set to 0 when there's no alignment.mism: number of mismatches in the alignment (zero or a positive integer).opens: number of columns containing a gap opening (zero or a positive integer, and terminal gaps don't count).qlo: first nucleotide of the query that's aligned with the target. Always equal to 1 if there's an alignment, 0 otherwise (see qilo if you want to ignore initial gaps).qhi: last nucleotide of the query aligned with the target. Always equal to the length of the pairwise alignment, 0 otherwise (see qihi to ignore terminal gaps).tlo: first nucleotide of the target aligned with the query. Always equal to 1 if there's an alignment, 0 otherwise (see tilo to ignore initial gaps).thi: last nucleotide of the target aligned with the query. Always equal to the length of the pairwise alignment, 0 otherwise (see tihi to ignore terminal gaps).evalue: the expectancy-value (not computed for nucleotide alignments). Always set to -1.bits: the bit score (not computed for nucleotide alignments). Always set to 0.
--dbfilename β compares your query sequences (the ones you hand over with--usearch_global) against the target sequences sitting in filename, which can be FASTA or FASTQ, using global pairwise alignment. Alternatively you can point it at a preformatted UDB database name β the kind you'd make with themakeudb_usearchcommand (more on that below). π¦--dbmasknone|dust|softβ masks regions in your target database sequences using thedustmethod, thesoftmethod, or... doesn't mask at all (none).
β οΈ Warning: when you use soft masking, search commands become case sensitive. Big deal, lowkey.
The default is to mask using dust.
--dbmatchedfilename β writes every database target sequence that matched at least one query sequence out to filename, in fasta format. And if you've got--sizeouton, the number of queries that hit each target sequence gets shown using the pattern";size=integer;". π--dbnotmatchedfilename β writes the database target sequences that didn't match any query sequences out to filename, in fasta format.--fastapairsfilename β writes the pairwise alignments of your query and target sequences to filename, in fasta format.--fulldpβ okay so this one is literally just a dummy option for compatibility with usearch. The tea is: to max out search sensitivity, vsearch uses an 8-way 16-bit SIMD vectorized full dynamic programming algorithm (Needleman-Wunsch) whether or not you pass--fulldp. So it does nothing extra. It's giving "we have to have it because usearch does." π--gapextstring β sets the penalties for extending a gap. For the full breakdown of how the penalty declaration system works, go read--gapopen. The default here initializes the six gap-extending penalties with a penalty of 2 for extending internal gaps and a penalty of 1 for extending terminal gaps, in both the query and the target sequences (that's2I/1E). β¨--gapopenstring β sets the penalties for opening a gap. A gap opening can happen in six different contexts: in the query (Q) or in the target (T) sequence, at the left (L) or right (R) extremity of the sequence, or inside the sequence (I). You can combine the sequence symbols (QandT) with the location symbols (L,I, andR), plus numerical values, to declare penalties for every possible context:
aQL/bQI/cQR/dTL/eTI/fTR
where abcdef are zero or positive integers, and / is the separator. To keep your declarations simple, the location symbols (L, I, and R) can be combined, the symbol E can be used to treat both extremities (L and R) equally, and the symbols Q and T can be left out entirely to treat query and target sequences equally. For example, the default is to declare a penalty of 20 for opening internal gaps and a penalty of 2 for opening terminal gaps (left or right), in both query and target sequences β that's 20I/2E. If you give just a numerical value with no sequence or location symbol at all, then that penalty applies to all gap openings. To forbid gap-opening, you can declare an infinite penalty value using the symbol *. And if you want to use vsearch as a semi-global aligner, you can slap a null-penalty on the left (L) or right (R) gaps. π
Here's the important bit: vsearch always initializes the six gap opening penalties using the default parameters (20I/2E). So you're free to declare only the values you actually want to change. The string gets scanned from left to right, the accepted symbols are (0123456789/LIREQT*), and later values override earlier ones. No cap. π₯
Please note that vsearch, unlike usearch, only allows integer gap penalties. Because usearch's lowest gap penalties default to 0.5, all the default scores and gap penalties in vsearch have been doubled β this keeps the penalties equivalent and produces identical alignments. Math, but make it slay. π
--hardmaskβ masks sequence regions by swapping them out forNs, instead of just lowercasing them the way it does by default. Want the full details? Go see the Masking section. π§--idreal β rejects a sequence match whenever the pairwise identity dips below real (a value from 0.0 to 1.0, inclusive). Here's how the search actually works: it sorts target sequences by how many k-mers they share with the query, in decreasing order, using that as a proxy for sequence similarity. That efficient pre-filtering also blocks pairwise alignments with super short targets or weakly matching ones β because by default you need at least 12 shared k-mers before the pairwise alignment even starts, and at least one out of every 16 k-mers from the query has to match the target (see--wordlengthand--minwordmatchesif you wanna change that behaviour).
So yeah, using values lower than --id 0.5 probably won't scoop up more weakly matching targets, bestie. The pairwise identity is by default defined as (matching columns) / (alignment length - terminal gaps). You can change that definition with --iddef. π
--iddef0|1|2|3|4β changes the pairwise identity definition that--iduses. Accepted values are:- CD-HIT definition: (matching columns) / (shortest sequence length).
- edit distance: (matching columns) / (alignment length).
- edit distance excluding terminal gaps (this is the default definition for
--id).
- Marine Biological Lab definition, counting each gap opening (internal or terminal) as a single mismatch, whether or not the gap was extended:
1.0 - [(mismatches + gap openings)/(longest sequence length)] - BLAST definition, equivalent to
--iddef 1for global pairwise alignments.
Okay so real talk, the --userfields option is gonna let you ask for the fields id0 all the way to id4, on top of the plain old field id, and that's how you get the pairwise identity values for each of these different definitions. Obsessed honestly. π
--idprefixpositive integer β Reject the sequence match if the first integer nucleotides of the target don't match the query. Like, if the vibes don't start right, we're out. π«--idsuffixpositive integer β Reject the sequence match if the last integer nucleotides of the target don't match the query. Same energy but at the end. Byeeeee. π--lca_cutoffreal β This one adjusts the fraction of matching hits required for the last common ancestor (LCA) output with the--lcaoutoption during searches. The default value is1.0, which requires all hits to match at each taxonomic rank for that rank to be included. If you use a lower cutoff value, e.g.0.95, a small fraction of non-matching hits are allowed while that rank will still be reported. The argument to this option must be larger than 0.5, but not larger than 1.0. No cap fr. π§¬--lcaoutfilename β Outputs last common ancestor (LCA) information about the hits of each query to a text file in a tab-separated format. The first column contains the query id, while the second column contains the taxonomic information. The headers of the sequences in the database must contain taxonomic information in the same format as used with the--sintaxcommand, e.g."tax=k:Archaea,p:Euryarchaeota,c:Halobacteria". Only the initial parts of the taxonomy that are common to a large fraction of the hits of each query will be output. It is necessary to set the--maxacceptsoption to a value different from 1 for this information to be useful. The--top_hits_onlyoption may also be useful. The fraction of matching hits required may be adjusted by the--lca_cutoffoption (default1.0). π³--leftjustβ Reject the sequence match if the pairwise alignment begins with gaps. Not a single gap at the start allowed bestie. π--lengthoutβ Write sequence length information to the output files in FASTA format by adding a";length=integer"attribute in the header. It's giving receipts. π--matchinteger β Score assigned to a match (i.e. identical nucleotides) in the pairwise alignment. The default value is2. β¨--matchedfilename β Write query sequences matching database target sequences to filename, in fasta format. Besties only. π«Ά--maxacceptspositive integer β Maximum number of matching target sequences to accept before stopping the search for a given query. The default value is1. This option works in pair with--maxrejects. The search process sorts target sequences by decreasing number of k-mers they have in common with the query sequence, using that information as a proxy for sequence similarity. After pairwise alignments, if the first target sequence passes the acceptation criteria, it is accepted as best hit and the search process stops for that query. If--maxacceptsis set to a higher value, more matching targets are accepted. If--maxacceptsand--maxrejectsare both set to0, the complete database is searched. See--maxhitsoption for a control on the number of hits reported per query when search is done on both strands. π₯--maxdiffspositive integer β Reject the sequence match if the alignment contains more than integer substitutions, insertions or deletions. That's the whole budget, don't overspend. πΈ--maxgapspositive integer β Reject the sequence match if the alignment contains more than integer insertions or deletions. Gaps are lowkey the enemy. π--maxhitsnon-negative integer β Maximum number of hits to show once the search is terminated for a given query (hits are sorted by decreasing identity). When searching only on the plus strand (default situation, see--strand), the number of matching targets (--maxaccepts) and the number of hits (--maxhits) are the same. However, when searching on both strands, there could be two hits per target (one per strand):--maxhitsthen controls the overall number of reported hits per query. Unlimited by default or if the argument is zero. This option applies to--alnout,--blast6out,--fastapairs,--samout,--uc, or--useroutoutput files. π--maxidreal β Reject the sequence match if the percentage of identity between the two sequences is greater than real. Yes you read that right, greater than. Too similar, rejected, we gatekeep here. π--maxqsizepositive integer β Reject query sequences with an abundance greater than integer. πΈ--maxqtreal β Reject if the query/target sequence length ratio is greater than real. π¦--maxrejectspositive integer β Maximum number of non-matching target sequences to consider before stopping the search for a given query. The default value is32. This option works in pair with--maxaccepts. The search process sorts target sequences by decreasing number of k-mers they have in common with the query sequence, using that information as a proxy for sequence similarity. After pairwise alignments, if none of the first32examined target sequences pass the acceptation criteria, the search process stops for that query (no hit). If--maxrejectsis set to a higher value, more target sequences are considered. If--maxacceptsand--maxrejectsare both set to0, the complete database is searched. π―--maxsizeratioreal β Reject if the query/target abundance ratio is greater than real. π--maxslreal β Reject if the shorter/longer sequence length ratio is greater than real. β--maxsubspositive integer β Reject the sequence match if the pairwise alignment contains more than integer substitutions. π--midreal β Reject the sequence match if the percentage of identity is lower than real (ignoring all gaps, internal and terminal). Gaps? Don't know her. π§--mincolspositive integer β Reject the sequence match if the alignment length is shorter than integer. Too tiny, not the moment. π--minqtreal β Reject if the query/target sequence length ratio is lower than real. π--minsizeratioreal β Reject if the query/target abundance ratio is lower than real. π--minslreal β Reject if the shorter/longer sequence length ratio is lower than real. π--mintsizepositive integer β Reject target sequences with an abundance lower than integer. Tiny abundance = blocked. π«--minwordmatchesnon-negative integer β Minimum number of k-mers or word matches required for a sequence to be considered further. Default value is12for the default word length8. For word lengths3-15, the default minimum word matches are18,17,16,15,14,12,11,10,9,8,7,5and3, respectively. If the query sequence has fewer unique words than the number specified, all words in the query must match. If the argument is0, no word matches are required. π--mismatchinteger β Score assigned to a mismatch (i.e. different nucleotides) in the pairwise alignment. The default value is-4. Ouch, negative, as it should be. π--mothur_shared_outfilename β Write search results to an OTU table in the mothur 'shared' tab-separated plain text file format. The query file contains the samples, while the database file contains the OTUs. Sample and OTU identifiers are extracted from the header of these sequences. See the--otutaboutoption in the Clustering section for further details. π¦--notmatchedfilename β Write query sequences not matching database target sequences to filename, in fasta format. The flops, we love them anyway. π--otutaboutfilename β Write search results to an OTU table in the classic tab-separated plain text format. The query file contains the samples, while the database file contains the OTUs. Sample and OTU identifiers are extracted from the header of these sequences (--sampleoption). See the--mothur_shared_outoption in the Clustering section for further details. π--output_no_hitsβ Write both matching and non-matching queries to--alnout,--blast6out,--samoutor--useroutoutput files. Non-matching queries are labelled'No hits'in--alnoutfiles. Everyone gets included, no one left behind fr. π«Ά--patternstring β This option is ignored. It is provided for compatibility with usearch. It literally does nothing bestie, it's just here to be polite. π--qmask none|dust|softβ masks out the boring/low-complexity bits of your query sequences using either thedustalgo or thesoftalgo, or literally do nothing at all (none). β οΈ Heads up bestie: if you go withsoftmasking, your search commands become case sensitive, so like, watch out. Default vibe is masking withdust. π--qsegout filenameβ dumps just the aligned part of every query sequence intofilename, in FASTA format. Only the part that actually lined up, not the whole thing. β¨--query_cov realβ rejects a match if the fraction of the query that aligned to the target is lower thanreal(a value from0.0to1.0, both included). Query coverage =(matches + mismatches) / query sequence length. And no, internal or terminal gaps do not get counted here. π«Ά--rightjustβ if the pairwise alignment ends in gaps, the match gets rejected. Bye. π--rowlen positive integerβ sets how wide the alignment lines are in your--alnoutoutput. Default is64. Set it to0and wrapping is gone completely β one long line, main character energy. π--samheaderβ adds header lines to the SAM file when you've got--samoutgoing. The header includes lines starting with@HD,@SQand@PG, but no@RGlines (peep <https://github.com/samtools/hts-specs> for the lore). By default, zero header lines get written. π--samout filenameβ writes your alignment results tofilenamein SAM format (a tab-separated text file, iconic). If you're using--samheader, the SAM file kicks off with header lines. Every non-header line is a SAM record, and each one is either a query-target alignment or the fact that a query had no match at all (and yeah, the output order can wiggle around when you're running multiple threads). Each record has 11 mandatory fields plus optional ones β full description lives at <https://github.com/samtools/hts-specs>:- query sequence label.
- a combo of bitwise flags. Possible values:
0(top hit),4(no hit),16(reverse-complemented hit),256(secondary hit, i.e. every hit that isn't the top hit). - target sequence label.
- first position of a target aligned with the query (always
1for global pairwise alignments,0if there's no match). - mapping quality (ignored, always set to
'*'). - CIGAR string (set to
'*'if there's no match). - name of the target sequence matching the next read of the query (mate reads only β ignored, always set to
'*'). - position of the primary alignment of the next read of the query (mate reads only β ignored, always set to
0). - target sequence length (for multi-segment targets β ignored, always set to
0). - query sequence (complete, not just the segment aligned to the target like usearch does, we're not doing that here).
- quality string (ignored, always set to
'*').
Optional fields for query-target matches (number and order of fields may vary, no promises):
AS:i:?alignment score (aka percentage of identity).XN:i:?next best alignment score (always set to0).XM:i:?number of mismatches.XO:i:?number of gap openings (excluding terminal gaps).XG:i:?number of gap extensions (excluding terminal gaps).NM:i:?edit distance to the target (that's the sum ofXMandXG).MD:Z:?string for mismatching positions.YT:Z:UUstring representing the alignment type.
--search_exact filenameβ searches for exact full-length matches to the query sequences infilename, inside the fasta database of target sequences (--db). Only 100% exact matches get reported, and this command is way faster than--usearch_global(fr fr).--id,--maxacceptsand--maxrejectsare ignored, but all the other searching options can still be specified. π₯--selfβ rejects the match if the query and target labels are identical. Same name? No. π--selfidβ rejects the match if the query and target sequences are strictly identical. Same sequence? Also no. π«--sizeoutβ adds abundance annotations to the output of--dbmatched(using the pattern';size=integer;'), so it reports how many queries matched each target. Counted and cute. π--strand plus|bothβ when searching for similar sequences, check the plus strand only (default) or check both strands. Both is giving bi-directional β¨--target_cov realβ rejects the match if the fraction of the target sequence aligned to the query is lower thanreal. Target coverage =(matches + mismatches) / target sequence length. Internal or terminal gaps? Not taken into account, nope. πΈ--top_hits_onlyβ only the top hits with an equally high percentage of identity between the query and database sequence sets get written to whatever output you picked with--lcaout,--alnout,--samout,--userout,--blast6out,--uc,--fastapairs,--matchedor--notmatchedβ but not--dbmatchedand--dbnotmatched, they're exempt bestie. For each query, the top hit is the one with the highest percentage of identity (check--iddefif you wanna change how identity is measured). And if several top hits tie with exactly the same percentage of identity, how many matching targets get reported is controlled by--maxaccepts(default1), and how many hits total is controlled by--maxhits. π―--tsegout filenameβ writes the aligned part of each target sequence tofilenamein FASTA format. (Yes, the target version of--qsegout, we love symmetry.) π--uc filenameβ outputs searching results intofilenameusing a tab-separated uclust-like format with 10 columns. When you use the--search_exactcommand, the table layout is the same as with--allpairs_global. When you use--usearch_global, the table has two different types of entries: hit (H) or no hit (N). Every query sequence gets compared to all other sequences, and the best hit (--maxaccepts 1) or several hits (--maxaccepts > 1together with--uc_allhits) get reported (H). Output order may vary when using multiple threads. Column content changes depending on the entry type (HorN):- Record type:
HorN('hit' or 'no hit'). - Ordinal number of the target sequence (based on input order, starting from zero). Set to
'*'forN. - Sequence length. Set to
'*'forN. - Percentage of similarity with the target sequence. Set to
'*'forN. - Match orientation
+or-. Set to'.'forN. - Not used, always set to zero for
H, or'*'forN. - Not used, always set to zero for
H, or'*'forN. - Compact representation of the pairwise alignment using the CIGAR format (Compact Idiosyncratic Gapped Alignment Report):
M(match/mismatch),D(deletion) andI(insertion). The equal sign'='means the query is identical to the centroid sequence (terminal gaps ignored). Set to'*'forN. - Label of the query sequence.
- Label of the target centroid sequence. Set to
'*'forN.
- Record type:
--uc_allhitsβ with the commands--search_exactand--usearch_global, when you're using the--ucoption, this shows all hits, not just the top hit for each query. Show everything, gatekeep nothing. π--usearch_global filenameβ compares target sequences (--db) to the query sequences infilename(FASTA or FASTQ format), using global pairwise alignment. That's the whole thing. π--userfields stringβ when you're using--userout, this selects and orders the fields written to the output file. Fields get separated by+(likequery+target+id, for example). The 'Userfields' section has the complete list of fields, go look. π§--userout filenameβ writes user-defined tab-separated output tofilename. Pick your fields with the--userfieldsoption. Output order may vary when using multiple threads. And if--userfieldsis empty or just not there at all?filenameends up empty. Sad but true. π
--weak_idreal β shows you the hits with a percentage of identity of at least real, and it does it without terminating the search π«Ά. A normal search stops the second it's found enough hits (where "enough" is defined by--maxaccepts,--maxrejects, and--id). Since--weak_idreports the weak hits that are not deduced from--maxaccepts(but they do count towards--maxrejects), you can run high--idvalues and still preserve both speed and sensitivity. Iconic honestly. Logically though, real has to be smaller than the value indicated by--id, ok? no cap.--wordlengthpositive integer β the length of the words (i.e. k-mers π ) used for database indexing. The range of possible values goes from3to15, but values near8or9are generally recommended. Longer words may reduce the sensitivity/recall for weak similarities, but they can increase precision. Shorter words, on the other hand, may increase sensitivity or recall, but may reduce precision. Computation time generally increases with shorter words and decreases with longer words β but then it increases again for very long words π. Memory requirements for a part of the index increase with a factor of4each time the word length increases by one nucleotide, and this generally becomes significant for long words (12or more). The default value is8.--xlengthβ strips the header attribute;length=integerfrom input sequences. That attribute is the one that gets added to output sequences by the--lengthoutoption.
β¨ Shuffling options β¨
Fasta entries in the input file are outputted in a pseudo-random order π²
--lengthoutβ writes sequence length information to the output files in FASTA format by adding a;length=integerattribute in the header.--outputfilename β writes the shuffled sequences to filename, in fasta format.--randseedinteger β when shuffling sequence order, use integer as the seed. A given seed always produces the same output order (super useful for replicability π§). Set it to0to use a pseudo-random seed (that's the default behaviour).--relabelstring β relabels sequences using the prefix string plus a ticker (1,2,3, etc.) to construct the new headers. Use--sizeoutto conserve the abundance annotations.--relabel_keepβ when relabelling, keeps the old identifier in the header, after a space.--relabel_md5β relabels sequences using the MD5 message digest algorithm applied to each sequence. Former sequence headers get discarded π. The sequence is converted to upper case andUis replaced byTbefore the digest is computed. The MD5 digest is a cryptographic hash function designed to minimize the probability that two different inputs give the same output, even for very similar, but non-identical inputs. Still, there's always a very small, but non-zero probability that two different inputs give the same result β so like, don't bet your whole thesis on it. The MD5 digest generates a128-bit (16-byte) digest that is represented by16hexadecimal numbers (using32symbols among0123456789abcdef). Use--sizeoutto conserve the abundance annotations.--relabel_selfβ relabels sequences using the sequence itself as the label. main character energy, fr.--relabel_sha1β relabels sequences using the SHA1 message digest algorithm applied to each sequence. It's similar to the--relabel_md5option but uses the SHA1 algorithm instead of the MD5 algorithm. The SHA1 digest generates a160-bit (20-byte) result that is represented by20hexadecimal numbers (40symbols). The probability of a collision (two non-identical sequences having the same digest) is smaller for the SHA1 algorithm than it is for the MD5 algorithm. Use--sizeoutto conserve the abundance annotations.--sizeoutβ when you're using--relabel,--relabel_self,--relabel_md5or--relabel_sha1, this preserves and reports abundance annotations to the output fasta file (using the pattern';size=integer;').--shufflefilename β pseudo-randomly shuffles the order of the sequences contained in filename.--topnpositive integer β outputs only the first integer sequences after the pseudo-random reordering.--xlengthβ strips the header attribute;length=integerfrom input sequences. (Yep, the one--lengthoutadds to output sequences. You know the vibe.)--xsizeβ strips abundance information from the headers when writing the output file.
β¨ Sorting options β¨
Fasta entries get sorted by decreasing abundance (--sortbysize) or by sequence length (--sortbylength). To obtain a stable sorting order, ties are sorted by decreasing abundance (if present) and label in increasing alpha-numerical order (--sortbylength), or just by label in increasing alpha-numerical order (--sortbysize). Label sorting assumes that all sequences have unique labels π. The same applies to the automatic sorting performed during chimera checking (--uchime_denovo), dereplication (--derep_fulllength), and clustering (--cluster_fast and --cluster_size).
--lengthoutβ writes sequence length information to the output files in FASTA format by adding a;length=integerattribute in the header.--maxsizepositive integer β when using--sortbysize, discards sequences with an abundance value greater than integer.--minsizepositive integer β when using--sortbysize, discards sequences with an abundance value smaller than integer.--outputfilename β writes the sorted sequences to filename, in fasta format.--relabelstring β please see the description of the same option under Chimera detection for details.--relabel_keepβ when relabelling, keeps the old identifier in the header, after a space.--relabel_md5β please see the description of the same option under Chimera detection for details.--relabel_selfβ please see the description of the same option under Chimera detection for details.--relabel_sha1β please see the description of the same option under Chimera detection for details.--sizeoutβ when using--relabel, reports abundance annotations to the output fasta file (using the pattern';size=integer;').--sortbylengthfilename β sorts by decreasing length the sequences contained in filename. Peep the general options--minseqlengthand--maxseqlengthif you wanna eliminate short and long sequences.--sortbysizefilename β sorts by decreasing abundance the sequences contained in filename (missing abundance values are assumed to be';size=1'). Peep the options--minsizeand--maxsizeto eliminate rare and dominant sequences.--topnpositive integer β outputs only the top integer sequences (i.e. the longest or the most abundant).--xlengthβ strips the header attribute;length=integerfrom input sequences. (Added to output sequences by--lengthout, as always π)--xsizeβ strips abundance information from the headers when writing the output file.
β¨ Subsampling options β¨
Subsampling randomly extracts a certain number or a certain percentage of the sequences in the input file. If the --sizein option is in effect, the abundances of the input sequences are taken into account, and the sampling is performed as if the input sequences were rereplicated, subsampled and dereplicated before being written to the output file π. The extraction is performed as a random sampling with a uniform distribution among the input sequences, and it's performed without replacement. The input file is specified with the --fastx_subsample option, the output files are specified with the --fastaout and --fastqout options, and the amount of sequences to be sampled is specified with the --sample_pct or --sample_size options. The sequences that didn't get sampled may be written to files specified with the options --fasta_discarded and --fastq_discarded. The --fastq_ascii, --fastq_qmin and --fastq_qmax options are also available.
--fastaoutfilename β writes the sampled sequences to filename, in fasta format.--fastaout_discardedfilename β writes the sequences not sampled to filename, in fasta format.--fastq_asciipositive integer β defines the ASCII character number used as the basis for the FASTQ quality score. The default is33, which is used by the Sanger / Illumina 1.8+ FASTQ format (phred+33). The value64is used by the Solexa, Illumina 1.3+ and Illumina 1.5+ formats (phred+64). Only33and64are valid arguments, bestie β don't get creative.--fastq_qmaxpositive integer β specifies the maximum quality score accepted when reading FASTQ files. The default is41, which is usual for recent Sanger/Illumina 1.8+ files.--fastq_qminpositive integer β the minimum quality score you'll accept in FASTQ files, bestie. Default is0, which is the norm for recent Sanger/Illumina 1.8+ files. Older formats may use scores between-5and2, so likeβ¦ watch out for the ancients π--fastqoutfilename β writes the sampled sequences out to filename, in fastq format. Requires the input to be in fastq format. No fastq in, no fastq out, periodt.--fastqout_discardedfilename β writes the sequences that did not get sampled to filename, in fastq format. Also requires fastq input. The leftovers get their own moment πΈ--fastx_subsamplefilename β does subsampling from the sequences in the specified input file, which can be FASTA or FASTQ. Flexible queen β¨--lengthoutβ writes sequence length info into the output files, in FASTA format, by adding a;length=integerattribute in the header. Cute little metadata moment.--randseedpositive integer β when subsampling, use integer as the seed for the pseudo-random generator. A given seed always produces the same output, which is so useful for replicability (science besties, this is your best friend). Set it to0to use a pseudo-random seed, which is the default behaviour.--relabelstring β relabels sequences using the prefix string plus a ticker (1,2,3, etc.) to build the new headers. Use--sizeoutto conserve the abundance annotations, okurrr? π--relabel_keepβ when relabelling, keep the old identifier in the header, after a space. Nothing gets ghosted π»--relabel_md5β relabels sequences using the MD5 message digest algorithm applied to each sequence. Former sequence headers are discarded. The sequence is converted to upper case andUis replaced byTbefore the digest is computed. The MD5 digest is a cryptographic hash function designed to minimize the probability that two different inputs give the same output β even for very similar, but non-identical inputs. Still, there is always a very small, but non-zero probability that two different inputs give the same result. The MD5 digest generates a128-bit (16-byte) digest that is represented by16hexadecimal numbers (using32symbols among0123456789abcdef). Use--sizeoutto conserve the abundance annotations.
not me being dramatic but that "non-zero probability" line is the whole vibe check of hashing, fr fr
--relabel_selfβ relabels sequences using the sequence itself as the label. Main character energy: the sequence IS the name π«Ά--relabel_sha1β relabels sequences using the SHA1 message digest algorithm applied to each sequence. It's similar to the--relabel_md5option but uses the SHA1 algorithm instead of the MD5 algorithm. The SHA1 digest generates a160-bit (20-byte) result that is represented by20hexadecimal numbers (40symbols). The probability of a collision (two non-identical sequences having the same digest) is smaller for the SHA1 algorithm than it is for the MD5 algorithm. Use--sizeoutto conserve the abundance annotations.--sample_pctreal β subsample the given percentage of the input sequences. Accepted values range from0.0to100.0. Math girlies rise π§--sample_sizepositive integer β extract the given number of sequences. Exact count, no guessing π--sizeinβ take the abundance information of the input file into account. Otherwise, the abundance of each sequence is considered to be1. Lowkey important, ngl.--sizeoutβ write abundance information to the output file. Gotta keep the receipts π--xlengthβ strips the header attribute;length=integerfrom the input sequences. Heads up: that's the exact attribute added to output sequences by the--lengthoutoption. They're a duo π€--xsizeβ strips abundance information from the headers when writing the output file.
π Taxonomic classification options π
The vsearch command --sintax will classify the input sequences according to the Sintax algorithm, as described by Robert Edgar (2016) in SINTAX: a simple non-Bayesian taxonomy classifier for 16S and ITS sequences, BioRxiv, 074161. Preprint. doi: 10.1101/074161 (link https://doi.org/10.1101/074161).
The name of the fasta file containing the input sequences to be classified is given as an argument to the --sintax command. The reference sequence database is specified with the --db option. The results are written in a tab delimited text file whose name is specified with the --tabbedout option. The --sintax_cutoff option may be used to set a minimum level of bootstrap support for the taxonomic ranks to be reported. The --randseed option may be included to specify a seed for initialisation of the random number generator used by the algorithm.
Please note that when using multiple threads, the--randseedoption may not work as intended, because sequences may be processed in a random order by different threads. To ensure the same results each time, use a single thread (--threads 1) in combination with a fixed random seed specified with--randseed. π the multithreading betrayal
Multithreading is supported. Databases in UDB files are supported. The strand option may be specified.
The reference database must contain taxonomic information in the header of each sequence, in the form of a string starting with ;tax= and followed by a comma-separated list of up to nine taxonomic identifiers. Each taxonomic identifier must start with an indication of the rank, by one of the letters d (for domain), k (kingdom), p (phylum), c (class), o (order), f (family), g (genus), s (species), or t (strain). The letter is followed by a colon (:) and the name of that rank. Commas and semicolons are not allowed in the name of the rank. Non-ascii characters should be avoided in the names.
Example:
>X80725_S000004313;tax=d:Bacteria,p:Proteobacteria,
c:Gammaproteobacteria,o:Enterobacteriales,f:Enterobacteriaceae,
g:Escherichia/Shigella,s:Escherichia_coli,
t:str._K-12_substr._MG1655
The option --notrunclabels is turned on by default for this command, allowing spaces in the taxonomic identifiers. Spaces stay, bestie π¦
If two sequences in the reference database have equally many kmer matches with the query, the shortest sequence will be chosen by default. If they are equally long, the sequence appearing first in the database will be chosen. If the recommended option --sintax_random is specified, sequences with an equal number of kmer matches will instead be chosen by a random draw. Tie-breaker drama resolved π²
--dbfilename β read the reference sequences from filename, in FASTA, FASTQ or UDB format. These sequences need to be annotated with taxonomy.--randseedpositive integer β use integer as seed for the random number generator used in the Sintax algorithm. A given seed always produces the same output order (useful for replicability). Set to0to use a pseudo-random seed (default behaviour). Does not work correctly with multiple threads; please use--threads 1to ensure correct behaviour.--sintaxfilename β read the input sequences from filename, in FASTA or FASTQ format.--sintax_cutoffreal β specify a minimum level of bootstrap support for the taxonomic ranks that will be included in column4of the output file. For instance0.9, corresponding to90%. Slay at 90 and above π--sintax_randomβ break ties between sequences with equally many kmer matches by a random draw. This option is recommended and may be made the default in the future. Early adopter energy β--tabbedoutfilename β write the results to filename, in a tab-separated text format.- Column
1contains the query label. - Column
2contains the predicted taxonomy in the same format as for the reference data, with bootstrap support indicated in parentheses after each rank. - Column
3contains the strand. - If the
--sintax_cutoffoption is used, the predicted taxonomy will be repeated in column4while omitting the bootstrap values and including only the ranks with support at or above the threshold.
- Column
π§ UDB options π§
Databases to be used with the --usearch_global or --sintax commands may be prepared from FASTA files and stored to a binary UDB formatted file in order to speed up searching. This may be worthwhile when searching a large database repeatedly. The sequences are indexed and stored in a way that can be quickly loaded into memory. The commands and options below can be used to create and inspect UDB files. Pre-building your database? That's the glow-up β¨
--dbmasknone|dust|softβ specify the sequence masking method used with the--makeudb_usearchcommand, eithernone,dustorsoft. No masking is performed whennoneis specified. Whendustis specified, the DUST algorithm will be used for masking low complexity regions (short repeats and skewed composition). Lower case letters in the input file will be masked whensoftis specified (soft masking).--hardmaskβ mask sequences by replacing letters withNfor the--makeudb_usearchcommand. The default is to use lower case letters (soft masking). Choose your masking aesthetic π--makeudb_usearchfilename β create a UDB database file from the FASTA-formatted sequences in the file with the given filename. The UDB database is written to the file specified with the--outputoption.--outputfilename β specify the filename of a FASTA or UDB output file for the--makeudb_usearchor the--udb2fastacommand, respectively.--udb2fastafilename β bestie, this one reads the UDB database from the file with the given filename and writes the sequences out in FASTA format into the file you set with the--outputoption. π--udbinfofilename β shows you info about the UDB database sitting in the file with that filename. β¨--udbstatsfilename β reports stats about the indexed words in the UDB database inside the file with that filename. π--wordlengthpositive integer β sets the length of the words used when you create the UDB database index with the--makeudb_usearchcommand. Valid numbers go from3to15. The default is8. π
π
Userfields (fields accepted by the --userfields option) π
alnβ prints a string ofM(match/mismatch, meaning not a gap),D(delete, meaning a gap in the query) andI(insert, meaning a gap in the target) that represents the pairwise alignment. Empty field if there is no alignment. π¦alnlenβ prints the length of the query-target alignment (the number of columns). The field is0if there is no alignment. πΈbitsβ bit score (not computed for nucleotide alignments). Always set to0. πcalnβ compact representation of the pairwise alignment using the CIGAR format (Compact Idiosyncratic Gapped Alignment Report):M(match/mismatch),D(deletion) andI(insertion). The equal sign=indicates that the query is identical to the centroid sequence (ignoring terminal gaps). Empty field if there is no alignment. β¨evalueβ E-value (not computed for nucleotide alignments). Always set to-1. πextsβ number of columns containing a gap extension (zero or positive integer value). π§gapsβ number of columns containing a gap (zero or positive integer value, excluding terminal gaps). πidβ the percentage of identity, according to the identity definition specified by the--iddefoption. Equal toid0,id1,id2,id3orid4below. By default the same asid2. π₯id0β CD-HIT definition of the percentage of identity (real value ranging from0.0to100.0) using the length of the shortest sequence in the pairwise alignment as denominator:100 * (matching columns) / (shortest sequence length). πid1β the percentage of identity (real value ranging from0.0to100.0) is defined as the edit distance:100 * (matching columns) / (alignment length). πid2β the percentage of identity (real value ranging from0.0to100.0) is defined as the edit distance, excluding terminal gaps. βid3β Marine Biological Lab definition of the percentage of identity (real value ranging from0.0to100.0), counting each gap opening (internal or terminal) as a single mismatch, whether or not the gap was extended, and using the length of the longest sequence in the pairwise alignment as denominator:100 * (1.0 - [(mismatches + gaps) / (longest sequence length)]). πid4β BLAST definition of the percentage of identity (real value ranging from0.0to100.0), equivalent to--iddef 1in a context of global pairwise alignment. The fieldid4is always equal to the fieldid1. πidsβ number of matches in the alignment (zero or positive integer value). πmismβ number of mismatches in the alignment (zero or positive integer value). πopensβ number of columns containing a gap opening (zero or positive integer value, excluding terminal gaps). πpairsβ number of columns containing only nucleotides. That value corresponds to the length of the alignment minus the gap-containing columns (zero or positive integer value). π¦pctgapsβ number of columns containing gaps expressed as a percentage of the alignment length (real value ranging from0.0to100.0). πΈpctpvβ percentage of positive columns. When working with nucleotide sequences, this is equivalent to the percentage of matches (real value ranging from0.0to100.0). β¨pvβ number of positive columns. When working with nucleotide sequences, this is equivalent to the number of matches (zero or positive integer value). πqcovβ fraction of the query sequence that is aligned with the target sequence (real value ranging from0.0to100.0). The query coverage is computed as100.0 * (matches + mismatches) / query sequence length. Internal or terminal gaps are not taken into account. The field is set to0.0if there is no alignment. πqframeβ query frame (-3to+3). That field only concerns coding sequences and is not computed byvsearch. Always set to+0. π§qhiβ last nucleotide of the query aligned with the target. Always equal to the length of the pairwise alignment,0otherwise (seeqihito ignore terminal gaps). πqihiβ last nucleotide of the query aligned with the target (ignoring terminal gaps). Nucleotide numbering starts from1. The field is set to0if there is no alignment. π₯qiloβ first nucleotide of the query aligned with the target (ignoring initial gaps). Nucleotide numbering starts from1. The field is set to0if there is no alignment. πqlβ query sequence length (positive integer value). The field is set to0if there is no alignment. πqloβ first nucleotide of the query aligned with the target. Always equal to1if there is an alignment,0otherwise (seeqiloto ignore initial gaps). βqrowβ prints the sequence of the query segment as seen in the pairwise alignment (meaning with gap insertions if need be). Empty field if there is no alignment. πqsβ query segment length. Always equal to query sequence length. πqstrandβ query strand orientation (+or-for nucleotide sequences). Empty field if there is no alignment. πqueryβ query label. πrawβ raw alignment score (negative, null or positive integer value). The score is the sum of match rewards minus mismatch penalties, gap openings and gap extensions. The field is set to0if there is no alignment. πtargetβ target label. The field is set to'*'if there is no alignment. π¦tcovβ fraction of the target sequence that is aligned with the query sequence (real value ranging from0.0to100.0). The target coverage is computed as100.0 * (matches + mismatches) / target sequence length. Internal or terminal gaps are not taken into account. The field is set to0.0if there is no alignment. πΈtframeβ target frame (-3to+3). That field only concerns coding sequences and is not computed byvsearch. Always set to+0. β¨thiβ last nucleotide of the target aligned with the query. Always equal to the length of the pairwise alignment,0otherwise (seetihito ignore terminal gaps). πtihiβ last nucleotide of the target aligned with the query (ignoring terminal gaps). Nucleotide numbering starts from1. The field is set to0if there is no alignment. πtiloβ first nucleotide of the target aligned with the query (ignoring initial gaps). Nucleotide numbering starts from1. The field is set to0if there is no alignment. π§tlβ target sequence length (positive integer value). The field is set to0if there is no alignment. πtloβ first nucleotide of the target aligned with the query. Always equal to1if there is an alignment,0otherwise (seetiloto ignore initial gaps). π₯trowβ prints the sequence of the target segment as seen in the pairwise alignment (meaning with gap insertions if need be). Empty field if there is no alignment. πtsβ target segment length. Always equal to target sequence length. The field is set to0if there is no alignment. πtstrandβ target strand orientation (+or-for nucleotide sequences). Always set to'+', so reverse strand matches havetstrand '+'andqstrand '-'. Empty field if there is no alignment. β
π DELIBERATE CHANGES π
okay so if you're a usearch user, bestie, the whole goal here is to make you feel at home. π That's literally why vsearch was designed to behave like usearch β to some extent. But like any complex software, usearch is not free from quirks and inconsistencies. We decided not to reproduce some of them, and for complete transparency, we're documenting every deliberate change right here. no secrets, we're not that kind of girl π
- During a search with
usearch, when using the options--blast6outand--output_no_hits, for queries with no match the number of fields reported is 13, where it should be 12. This is corrected invsearch. (a little extra field just vibing in there, not anymore π«Ά) - The field
rawof the--userfieldsoption is not informative inusearch. This is corrected invsearch. - The fields
qlo,qhi,tlo,thinow have counterparts (qilo,qihi,tilo,tihi) reporting alignment coordinates ignoring terminal gaps. β¨ - In
usearch, when using the option--output_no_hits, queries that receive no match are reported in the--blast6outfile, but not in the alignment output file. This is corrected invsearch. vsearchintroduces a new--cluster_sizecommand that sorts sequences by decreasing abundance before clustering. πvsearchreintroduces--iddefalternative pairwise identity definitions that were removed fromusearch.vsearchextends the--topnoption to sorting commands.vsearchextends the--sizeinoption to dereplication (--derep_fulllength) and clustering (--cluster_fast).vsearchtreatsTandUas identical nucleotides during dereplication. π§¬vsearchsorting is stabilized by using sequence abundances or sequences labels as secondary or tertiary keys.vsearchby default uses the DUST algorithm for masking low-complexity regions. Masking behaviour is also slightly changed to be more consistent.
π NOVELTIES π
vsearch introduces new commands and new options not present in usearch 7. They're all described in the 'Options' section of this manual. Here's the short list, because you deserve a lil preview:
uchime2_denovo,uchime3_denovo,alignwidth,borderline,fasta_score(chimera checking) π¦cluster_size,cluster_unoise,clusterout_id,clusterout_sort,profile(clustering)fasta_width,gzip_decompress,bzip2_decompress(general option)iddef(clustering, pairwise alignment, searching)maxuniquesize(dereplication)relabel_md5,relabel_selfandrelabel_sha1(chimera detection, dereplication, FASTQ processing, shuffling, sorting)shuffle(shuffling)fastq_eestats,fastq_eestats2,fastq_maxlen,fastq_truncee(FASTQ processing)fastaout_discarded,fastqout_discarded(subsampling)rereplicate(dereplication/rereplication)
β¨ EXAMPLES β¨
bestie this is the part where it all clicks. let's go π«Ά
Align all sequences in a database with each other and output all pairwise alignments:
vsearch --allpairs_global database.fas --alnout results.aln --acceptall
Check for the presence of chimeras (de novo); parents should be at least 1.5 times more abundant than chimeras. Output non-chimeric sequences in fasta format (no wrapping):
vsearch --uchime_denovo queries.fas --abskew 1.5 --nonchimeras results.fas --fasta_width 0
Cluster with a 97% similarity threshold, collect cluster centroids, and write cluster descriptions using a uclust-like format:
vsearch --cluster_fast queries.fas --id 0.97 --centroids centroids.fas --uc clusters.uc
Dereplicate the sequences contained in queries.fas, take into account the abundance information already present, write unwrapped fasta sequences to queries_unique.fas with the new abundance information, discard all sequences with an abundance of 1:
vsearch --derep_fulllength queries.fas --sizein --fasta_width 0 --sizeout --output queries_unique.fas --minuniquesize 2
Mask simple repeats and low complexity regions in the input fasta file with the DUST algorithm (masked regions are lowercased), and write the results to the output file:
vsearch --maskfasta queries.fas --qmask dust --output queries_masked.fas
Search queries in a reference database, with a 80%-similarity threshold, take terminal gaps into account when calculating pairwise similarities, output pairwise alignments:
vsearch --usearch_global queries.fas --db references.fas --id 0.8 --iddef 1 --alnout results.aln
Search a sequence dataset against itself (ignore self hits), get all matches with at least 60% similarity, and collect results in a blast-like tab-separated format. Accept an unlimited number of hits (--maxaccepts 0), and compare each query to all other sequences, including unlikely candidates (--maxrejects 0):
vsearch --usearch_global queries.fas --db queries.fas --self --id 0.6 --blast6out results.blast6 --maxaccepts 0 --maxrejects 0
Shuffle the input fasta file (change the order of sequences) in a repeatable fashion (fixed seed), and write unwrapped fasta sequences to the output file:
vsearch --shuffle queries.fas --output queries_shuffled.fas --randseed 13 --fasta_width 0
Sort by decreasing abundance the sequences contained in queries.fas (using the size=integer information), relabel the sequences while preserving the abundance information (with --sizeout), keep only sequences with an abundance equal to or greater than 2:
vsearch --sortbysize queries.fas --output queries_sorted.fas --relabel sampleA_ --sizeout --minsize 2
π AUTHORS π
Implementation and documentation by TorbjΓΈrn Rognes, FrΓ©dΓ©ric MahΓ© and TomΓ‘s Flouri. absolute icons, no cap π
π CITATION π
Rognes T, Flouri T, Nichols B, Quince C, MahΓ© F. (2016) VSEARCH: a versatile open source tool for metagenomics. PeerJ 4:e2584 doi: 10.7717/peerj.2584 (link) β¨https://doi.org/10.7717/peerj.2584β©
π«Ά REPORTING BUGS π«Ά
Submit suggestions and bug-reports at (link) β¨https://github.com/torognes/vsearch/issuesβ© <https://github.com/torognes/vsearch/issues>, send a pull request on (link) β¨https://github.com/torognes/vsearchβ© <https://github.com/torognes/vsearch>, or compose a friendly or curmudgeont e-mail to TorbjΓΈrn Rognes (link) β¨torognes@ifi.uio.noβ© <torognes@ifi.uio.no>.
ngl "curmudgeont" is spelled exactly like that in the source, bestie β we're not touching it, that's β¨ heritage β¨
πΈ AVAILABILITY πΈ
Source code and binaries are available at <https://github.com/torognes/vsearch>.
π COPYRIGHT π
Copyright (C) 2014-2026, TorbjΓΈrn Rognes, FrΓ©dΓ©ric MahΓ© and TomΓ‘s Flouri
All rights reserved.
Contact: TorbjΓΈrn Rognes <torognes@ifi.uio.no>, Department of Informatics, University of Oslo, PO Box 1080 Blindern, NO-0316 Oslo, Norway
This software is dual-licensed and available under a choice of one of two licenses, either under the terms of the GNU General Public License version 3 or the BSD 2-Clause License.
π§ GNU General Public License version 3 π§
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see (link) β¨https://www.gnu.org/licenses/β© <https://www.gnu.org/licenses/>.
π The BSD 2-Clause License π
Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
- Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
- Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
We would like to thank the authors of the following projects for making their source code available:
vsearchincludes code from Google's CityHash project by Geoff Pike and Jyrki Alakuijala, providing some excellent hash functions available under a MIT license. π₯vsearchincludes code derived from Tatusov and Lipman's DUST program that is in the public domain.vsearchincludes public domain code written by Alexander Peslyak for the MD5 message digest algorithm.vsearchincludes public domain code written by Steve Reid and others for the SHA1 message digest algorithm.vsearchbinaries may include code from the zlib library, copyright Jean-Loup Gailly and Mark Adler.vsearchbinaries may include code from the bzip2 library, copyright Julian R. Seward.
π¦ SEE ALSO π¦
swipe, an extremely fast pairwise local (Smith-Waterman) database search tool by TorbjΓΈrn Rognes, available at (link) β¨https://github.com/torognes/swipeβ© <https://github.com/torognes/swipe>.
swarm, a fast and accurate amplicon clustering method by FrΓ©dΓ©ric MahΓ© and TorbjΓΈrn Rognes, available at (link) β¨https://github.com/torognes/swarmβ© <https://github.com/torognes/swarm>.
β VERSION HISTORY β
New features and important modifications of vsearch (short lived or minor bug releases may not be mentioned):
- v1.0.0 released November 28th, 2014
- First public release. π
- v1.0.1 released December 1st, 2014
- Bug fixes (
sortbysize, semicolon after size annotation in headers) and minor changes (labels as secondary sort key for most sorts, treatTandUas identical for dereplication, only output size in--dbmatchedfile if--sizeoutspecified).
- Bug fixes (
- v1.0.2 released December 6th, 2014
- Bug fixes (ssse3/sse4.1 requirement, memory leak).
- v1.0.3 released December 6th, 2014
- Bug fix (now writes help to
stdoutinstead ofstderr).
- Bug fix (now writes help to
v1.0.4 released December 8th, 2014 Added the --allpairs_global option. Also lowkey trimmed the memory requirements a lil and yeeted some memory leaks π
v1.0.5 released December 9th, 2014 Fixes a minor bug with the --allpairs_global and --acceptall options π
v1.0.6 released December 14th, 2014 Fixes a memory allocation bug in chimera detection (that's the --uchime_ref option bestie) π§
v1.0.7 released December 19th, 2014 Fixes a bug in the output from chimera detection when you're using the --uchimeout option β¨
v1.0.8 released January 22nd, 2015 Introduces several changes and bug fixes:
- a new linear memory aligner for alignment of sequences longer than 5,000 nucleotides π
- a new
--cluster_sizecommand that sorts sequences by decreasing abundance before clustering, - meaning of userfields
qlo,qhi,tlo,thichanged for compatibility with usearch, - new userfields
qilo,qihi,tilo,tihigive alignment coordinates ignoring terminal gaps, - in
--ucoutput files, a perfect alignment is indicated with a'='sign, - the option
--cluster_fastnow sorts sequences by decreasing length, then by decreasing abundance and finally by sequence identifier, - default
--maxseqlengthvalue set to 50,000 nucleotides, - fix for bug in alignment in rare cases,
- fix for lack of detection of under- or overflow in SIMD aligner.
v1.0.9 released January 22nd, 2015 Fixes a bug in the function sorting sequences by decreasing abundance (--sortbysize) π
v1.0.10 released January 23rd, 2015 Fixes a bug where the --sizein option was ignored and always treated as on, affecting clustering and dereplication commands π
v1.0.11 released February 5th, 2015 Introduces the possibility to output results in SAM format (for clustering, pairwise alignment and searching) π
v1.0.12 released February 6th, 2015 Temporarily fixes a problem with long headers in FASTA files π
v1.0.13 released February 17th, 2015 Fix a memory allocation problem when computing multiple sequence alignments with the --msaout and --consout options, as well as a memory leak. Also increased line buffer for reading FASTA files to 4MB.
v1.0.14 released February 17th, 2015 Fix a bug where the multiple alignment and consensus sequence computed after clustering ignored the strand of the sequences. Also decreased size of line buffer for reading FASTA files to 1MB again due to excessive stack memory usage.
v1.0.15 released February 18th, 2015 Fix bug in calculation of identity metric between sequences when using the MBL definition (--iddef 3) π₯
v1.0.16 released February 19th, 2015 Integrated patches from Debian for increased compatibility with various architectures π«Ά
v1.1.0 released February 20th, 2015 Added the --quiet option to suppress all output to stdout and stderr except for warnings and fatal errors. Added the --log option to write messages to a log file β¨
v1.1.1 released February 20th, 2015 Added info about --log and --quiet options to help text.
v1.1.2 released March 18th, 2015 Fix bug with large datasets. Fix format of help info.
v1.1.3 released March 18th, 2015 Fix more bugs with large datasets π
v1.2.0-1.2.19 released July 6th to September 8th, 2015 Several new commands and options added. Bugs fixed. Documentation updated.
v1.3.0 released September 9th, 2015 Changed to autotools build system π οΈ
v1.3.1 released September 14th, 2015 Several new commands and options. Bug fixes.
v1.3.2 released September 15th, 2015 Fixed memory leaks. Added '-h' shortcut for help. Removed extra 'v' in version number.
v1.3.3 released September 15th, 2015 Fixed bug in hexadecimal digits of MD5 and SHA1 digests. Added --samheader option.
v1.3.4 released September 16th, 2015 Fixed compilation problems with zlib and bzip2lib.
v1.3.5 released September 17th, 2015 Minor configuration/makefile changes to compile to native CPU and simplify makefile.
v1.4.0 released September 25th, 2015 Added --sizeorder option π
v1.4.1 released September 29th, 2015 Inserted public domain MD5 and SHA1 code to eliminate dependency on crypto and openssl libraries and their licensing issues π
v1.4.2 released October 2nd, 2015 Dynamic loading of libraries for reading gzip and bzip2 compressed files if available. Circumvention of missing gzoffset function in zlib 1.2.3 and earlier.
v1.4.3 released October 3rd, 2015 Fix a bug with determining amount of memory on some versions of Apple OS X π
v1.4.4 released October 3rd, 2015 Remove debug message.
v1.4.5 released October 6th, 2015 Fix memory allocation bug when reading long FASTA sequences.
v1.4.6 released October 6th, 2015 Fix subtle bug in SIMD alignment code that reduced accuracy.
v1.4.7 released October 7th, 2015 Fixes a problem with searching for or clustering sequences with repeats. In this new version, vsearch looks at all words occurring at least once in the sequences in the initial step. Previously only words occurring exactly once were considered. In addition, vsearch now requires at least 10 words to be shared by the sequences, previously only 6 were required. If the query contains less than 10 words, all words must be present for a match. This change seems to lead to slightly reduced recall, but somewhat increased precision, ending up with slightly improved overall accuracy. π§
v1.5.0 released October 7th, 2015 This version introduces the new option --minwordmatches that allows the user to specify the minimum number of matching unique words before a sequence is considered further. New default values for different word lengths are also set. The minimum word length is increased to 7.
v1.6.0 released October 9th, 2015 This version adds the relabeling options (--relabel, --relabel_md5 and --relabel_sha1) to the shuffle command. It also adds the --xsize option to the clustering, dereplication, shuffling and sorting commands β¨
v1.6.1 released October 14th, 2015 Fix bugs and update manual and help text regarding relabelling. Add all relabelling options to the subsampling command. Add the --xsize option to chimera detection, dereplication and fastq filtering commands. Refactoring of code π οΈ
v1.7.0 released October 14th, 2015 Add --relabel_keep option π
v1.8.0 released October 19th, 2015 Added --search_exact, --fastx_mask and --fastq_convert commands. Changed most commands to read FASTQ input files as well as FASTA files. Modified --fastx_revcomp and --fastx_subsample to write FASTQ files π
v1.8.1 released November 2nd, 2015 Fixes for compatibility with QIIME and older OS X versions π
v1.9.0 released November 12th, 2015 Added the --fastq_mergepairs command and associated options. This command has not been tested well yet. Included additional files to avoid dependency of autoconf for compilation. Fixed an error where identifiers in fasta headers where not truncated at tabs, just spaces. Fixed a bug in detection of the file format (FASTA/FASTQ) of a gzip compressed input file π
v1.9.1 released November 13th, 2015 Fixed memory leak and a bug in score computation in --fastq_mergepairs, and improved speed π₯
v1.9.2 released November 17th, 2015 Fixed a bug in the computation of some values with --fastq_stats.
v1.9.3 released November 19th, 2015 Workaround for missing x86intrin.h with old compilers.
v1.9.4 released December 3rd, 2015 Fixed incrementation of counter when relabeling dereplicated sequences.
v1.9.5 released December 3rd, 2015 Fixed bug resulting in inferior chimera detection performance π
v1.9.6 released January 8th, 2016 Fixed bug in aligned sequences produced with --fastapairs and --userout (qrow, trow) options.
v1.9.7 released January 12th, 2016 Masking behaviour is changed somewhat to keep the letter case of the input sequences unchanged when no masking is performed. Masking is now performed also during chimera detection. Documentation updated β¨
v1.9.8 released January 22nd, 2016 Fixed bug causing segfault when chimera detection is performed on extremely short sequences π±
v1.9.9 released January 22nd, 2016 Adjusted default minimum number of word matches during searches for improved performance.
v1.9.10 released January 25th, 2016 Fixed bug related to masking and lower case database sequences.
v1.10.0 released February 11th, 2016 Parallelized and improved merging of paired-end reads and adjusted some defaults. Removed progress indicator when stderr is not a terminal. Added --fasta_score option to report chimera scores in FASTA files. Added --rereplicate and --fastq_eestats commands. Fixed typos. Added relabelling to files produced with --consout and --profile options π
v1.10.1 released February 23rd, 2016 Fixed a bug affecting the --fastq_mergepairs command causing FASTQ headers to be truncated at first space (despite the bug fix release 1.9.0 of November 12th, 2015). Full headers are now included in the output (no matter if --notrunclabels is in effect or not) π«Ά
v1.10.2 released March 18th, 2016 Fixed a bug causing a segmentation fault when running --usearch_global with an empty query sequence. Also fixed a bug causing imperfect alignments to be reported with an alignment string of '=' in uc output files. Fixed typos in man file. Fixed fasta/fastq processing code regarding presence or absence of compression library header files π
v1.11.1 released April 13th, 2016 Added strand information in UC file for --derep_fulllength and --derep_prefix. Added expected errors (ee) to header of FASTA files specified with --fastaout and --fastaout_discarded when --eeout or --fastq_eeout option is in effect for fastq_filter and fastq_mergepairs. The options --eeout and --fastq_eeout are now equivalent πΈ
v1.11.2 released June 21st, 2016
okay so two bugs got absolutely obliterated here. bug numero uno was all about the --query_cov option, which was lowkey using a whole different coverage definition than the qcov userfield was. not the vibe. so now coverage is defined as the fraction of the whole query sequence length that's aligned with matching or mismatching residues in the target β and bestie, all gaps are ignored, periodt.
bug numero dos? it was living in the consensus sequences made during clustering, specifically when only N's were hanging out in some positions. before, those N's got converted to A's in the consensus (π not it). now the behaviour is changed so N's are produced in the consensus instead, and it should be way more compatible with usearch. slay.
v2.0.0 released June 24th, 2016
MAJOR version drop, this is a big one π
it now supports reading from pipes. two brand new options joined the chat: --gzip_decompress and --bzip2_decompress. one of these must be specified if you're reading compressed input from a pipe β but they're not required when reading from ordinary files, no cap.
also the vsearch header that used to get written to stdout is now written to stderr, which unlocks the ability to pipe results along for further processing β¨ and the file name - now represents standard input (/dev/stdin) or standard output (/dev/stdout) when reading or writing files, respectively. the code for reading FASTA and FASTQ files also got refactored. glow up fr.
v2.0.1 released June 30th, 2016
avoided a segmentation fault when masking very long sequences. we do not do segfaults here π
v2.0.2 released July 5th, 2016
avoided warnings when compiling with GCC 6. quiet, clean, iconic.
v2.0.3 released August 2nd, 2016
fixed bad compiler options that were causing Illegal instruction errors when running precompiled binaries. that was giving chaos energy, now it's giving stability π«Ά
v2.0.4 released September 1st, 2016
improved the error message for bad FASTQ quality values, and gave the manual a little makeover. self care πΈ
v2.0.5 released September 9th, 2016
added the options --fastaout_discarded and --fastqout_discarded, so you can output discarded sequences from subsampling to separate files. updated the manual too. organized queen behavior π
v2.1.0 released September 16th, 2016
new command alert π¨ --fastx_filter. new options: --fastq_maxlen, --fastq_truncee. and --minwordmatches can now go all the way down to 3.
v2.1.1 released September 23rd, 2016
fixed bugs in output to UC-files, plus improved help text and the manual. we love a documentation refresh β¨
v2.1.2 released September 28th, 2016
fixed incorrect abundance output from fastx_filter and fastq_filter when relabelling. numbers were lying, now they're not. π
v2.2.0 released October 7th, 2016
added OTU table generation options β --biomout, --mothur_shared_out and --otutabout β to the clustering and searching commands. tables for days π
v2.3.0 released October 10th, 2016
allowed zero-length sequences in FASTA and FASTQ files (yes, literally zero, we contain multitudes). added the --fastq_trunclen_keep option. and fixed a bug with output of OTU tables to pipes.
v2.3.1 released November 16th, 2016
fixed a bug where --minwordmatches 0 was being interpreted as the default minimum word matches for the given word length instead of actual zero. rude. now, when you use it in combination with --maxaccepts 0 and --maxrejects 0, it will allow complete bypass of kmer-based heuristics. full send πβοΈ
v2.3.2 released November 18th, 2016
fixed a bug where vsearch was reporting the ordinal number of the target sequence instead of the cluster number in column 2 on H-lines in the uc output file after clustering. bestie, for search and alignment commands both usearch and vsearch report the target sequence number here. consistency is the moment.
v2.3.3 released December 5th, 2016
a minor speed improvement. zoom zoom β‘
v2.3.4 released December 9th, 2016
fixed a bug in output of sequence profiles and updated the documentation. two birds, one slay π¦
v2.4.0 released February 8th, 2017
added support for Linux on Power8 systems (ppc64le) and Windows on x86_64. improved detection of pipes when reading FASTA and FASTQ files. and corrected the option for specifying output from the fastq_eestats command in the help text. global expansion era π
v2.4.1 released March 1st, 2017
fixed an overflow bug in fastq_stats and fastq_eestats that was affecting analysis of very large FASTQ files. also fixed maximum memory usage reporting on Windows. big files, no problem π«Ά
v2.4.2 released March 10th, 2017
default value for fastq_minovlen increased to 16, in accordance with the help text and for compatibility with usearch. plus minor changes for improved accuracy of paired-end read merging. accuracy is the aesthetic β¨
v2.4.3 released April 6th, 2017
fixed a bug with the progress bar for shuffling. and fixed missing N-lines in UC files with usearch_global, search_exact and allpairs_global when the output_no_hits option was not specified. nothing missing, everything accounted for π
v2.4.4 released August 28th, 2017
fixed a few minor bugs, improved error messages and updated the documentation. a lil maintenance glow-up π
v2.5.0 released October 5th, 2017
support for UDB database files ποΈ new commands: fastq_stripright, fastq_eestats2, makeudb_usearch, udb2fasta, udbinfo, and udbstats. new general option: no_progress. new options minsize and maxsize for fastx_filter. plus minor bug fixes, error message improvements and documentation updates. this release ATE.
v2.5.1 released October 25th, 2017
fixed a bug with a bad default value of 1 instead of 32 for minseqlength when using the makeudb_usearch command. 32, bestie. not 1. π
v2.5.2 released October 30th, 2017
fixed a bug where - as an argument to the fastq_eestats2 option was treated literally instead of being equivalent to stdin. little typo, big consequences π
v2.6.0 released November 10th, 2017
rewritten paired-end reads merger with improved accuracy π₯ decreased the default value for the fastq_minovlen option from 16 to 10, and the default value for the fastq_maxdiffs option is increased from 5 to 10. but heads up β there are now other, more important restrictions that will avoid merging reads that cannot be reliably aligned. safety first, always π‘οΈ
v2.6.1 released December 8th, 2017
improved parallelisation of paired end reads merging. multitasking icon β‘
v2.6.2 released December 18th, 2017
fixed the option xsize, which was partially inactive for the commands uchime_denovo, uchime_ref, and fastx_filter. it's awake now bestie π€
v2.7.0 released February 13th, 2018
added the commands cluster_unoise, uchime2_denovo and uchime3_denovo, contributed by Davide Albanese based on Robert Edgar's papers π also refactored the fasta and fastq print functions, as well as the code for extraction of abundance and other attributes from the headers.
v2.7.1 released February 16th, 2018
fix several bugs on Windows related to large files, use of - as a file name to mean stdin or stdout, alignment errors, missed kmers and corrupted UDB files. also added documentation of UDB-related commands. windows girlies, this one's for you π«Ά
v2.7.2 released April 20th, 2018
added the sintax command for taxonomic classification 𧬠and fixed a bug with incorrect FASTA headers of consensus sequences after clustering.
v2.8.0 released April 24th, 2018
added the fastq_maxdiffpct option to the fastq_mergepairs command. one tiny option, infinite possibilities β¨
v2.8.1 released June 22nd, 2018
fixes for compilation warnings with GCC 8. silence is golden π€«
v2.8.2 released August 21st, 2018
fix for wrong placement of semicolons in header lines in some cases when using the sizeout or xsize options. reduced memory requirements for full-length dereplication in cases with many duplicate sequences. improved wording of the fastq_mergepairs report. updated the manual regarding use of sizein and sizeout with dereplication. and changed a compiler option. a whole buffet of fixes π
v2.8.3 released August 31st, 2018
fix for segmentation fault for --derep_fulllength with --uc. we do not crash, we thrive π
v2.8.4 released September 3rd, 2018
further reduce memory requirements for dereplication when not using the uc option. also fixed output during subsampling when quiet or log options are in effect. leaner and cleaner πΈ
v2.8.5 released September 26th, 2018
fixed a bug in fastq_eestats2 that caused the values for large lengths to be much too high when the input sequences had varying lengths. the numbers were exaggerating. now they're honest π
v2.8.6 released October 9th, 2018
fixed a bug introduced in version 2.8.2 that caused derep_fulllength to include the full FASTA header in its output instead of stopping at the first space (unless the notrunclabels option is in effect). my bad, my bad, fixed now π
v2.9.0 released October 10th, 2018
added the fastq_join command. new tool, who dis π₯
v2.9.1 released October 29th, 2018
changed compiler options that select the target cpu and tuning, so the software can run on any 64-bit x86 system while still tuning for more modern variants. also avoided the illegal instruction error on some architectures, and updated documentation of the rereplicate command. everybody gets to play now π«Ά
v2.10.0 released December 6th, 2018
added the sff_convert command to convert SFF files to FASTQ π added some additional option argument checks. and fixed a segmentation fault bug that happened after some fatal errors when a log file was specified.
v2.10.1 released December 7th, 2018
improved the sff_convert command β it will now read several variants of the SFF format, and it's also able to read from a pipe β¨ warnings are given if there are minor problems, error messages have been improved, and there are minor speed and memory usage improvements. glow up complete.
v2.10.2 released December 10th, 2018
fixed a bug in sintax with reversed order of domain and kingdom. taxonomy stays in its lane now π
v2.10.3 released December 19th, 2018
ported to Linux on ARMv8 (aarch64) π¦Ύ and fixed a compilation warning with gcc version 8.1.0 and 8.2.0.
v2.10.4 released January 4th, 2019
fixed a serious bug in x86_64 SIMD alignment code introduced in version 2.10.3. added a link to BioConda in the README. fixed a bug in fastq_stats with sequence length 1. and fixed use of the equals symbol in UC files for identical sequences with cluster_fast. busy day, ate it all π
v2.11.0 released February 13th, 2019
added the ability to trim and filter paired-end reads using the reverse option with the fastx_filter and fastq_filter commands. added the --xee option to remove ee attributes from FASTA headers. plus a minor invisible improvement to the progress indicator. invisible but iconic β¨
v2.11.1 released February 28th, 2019
minor change to the handling of the weak_id and id options when using cluster_unoise. small tweak, big vibes.
v2.12.0 released March 19th, 2019
take sequence abundance into account when computing consensus sequences or profiles after clustering π warn when rereplicating sequences without abundance info. guess offset 33 in more cases with fastq_chars. and stricter checking of option arguments and option combinations. strict queen behavior, we love to see it π
π v2.13.0 β released April 11th, 2019
- Added the
--fastx_getseq,--fastx_getseqsand--fastx_getsubseqcommands, which pull sequences out of a FASTA or FASTQ file based on their labels π - Improved handling of ambiguous nucleotide symbols
- Corrected behaviour of the
--uchime_refcommand with and options--selfand--selfid - Strict detection of illegal options for each command
β¨ v2.13.1 β released April 26th, 2019
Minor changes to the allowed options for each command. All commands now allow the log, quiet and threads options. If more than 1 thread is specified for commands that are not multi-threaded, a warning will be issued. Minor changes to the manual.
πΈ v2.13.2 β released April 30th, 2019
Fixed bug related to improper handling of newlines on Windows. Allowed option strand plus to uchime_ref for compatibility.
π¦ v2.13.3 β released April 30th, 2019
Fixed bug in FASTQ parsing introduced in version 2.13.2. (the way we broke it and then un-broke it, iconic behaviour honestly π)
π v2.13.4 β released May 10th, 2019
Added information about support for gzip- and bzip2-compressed input files to the output of the version command. Adapted source code for compilation on FreeBSD and NetBSD systems.
π v2.13.5 β released July 2nd, 2019
- Added
cutcommand to fragment sequences at restriction sites - Silenced output from the
fastq_statscommand ifquietoption was given - Updated manual
π v2.13.6 β released July 2nd, 2019
Added info about cut command to output of help command.
β v2.13.7 β released September 2nd, 2019
Fixed bug in consensus sequence introduced in version 2.13.0.
π v2.14.0 β released September 11th, 2019
Added relabel_self option. Made fasta_width, sizein, sizeout and relabelling options valid for certain commands.
π«Ά v2.14.1 β released September 18th, 2019
Fixed bug with sequences written to file specified with fastaout_rev for commands fastx_filter and fastq_filter.
π§ v2.14.2 β released January 28th, 2020
Fixed some issues with the cut, fastx_revcomp, fastq_convert, fastq_mergepairs, and makeudb_usearch commands. Updated manual.
π v2.15.0 β released June 19th, 2020
Bestie this one is a whole glow-up, sit down:
- Update manual and documentation
- Turn on
notrunclabelsoption forsintaxcommand by default - Change
maxhits 0to mean unlimited hits, like the default - Allow non-ascii characters in headers, with a warning
- Sort centroids and uc too when
clusterout_sortspecified - Add cluster id to centroids output when
clusterout_idspecified - Improve error messages when parsing FASTQ files
- Add missing
fastq_qminoutoption and fixlabel_suffixoption forfastq_mergepairs - Add
derep_idcommand that dereplicates based on both label and sequence - Remove compilation warnings
π v2.15.1 β released October 28th, 2020
- Fix for dereplication when including reverse complement sequences and headers
- Make some extra checks when loading compression libraries and add more diagnostic output about them to the output of the
versioncommand - Report an error when
fastx_filteris used with FASTA input and options that require FASTQ input - Update manual
π v2.15.2 β released January 26th, 2021
No real functional changes, but some code and compilation changes. Compiles successfully on macOS running on Apple Silicon (ARMv8). Binaries available. Code updated for C++11. Minor adaptations for Windows compatibility, including the use of the C++ standard library for regular expressions. Minor changes for compatibility with Power8. Switch to C++ header files.
π₯ v2.16.0 β released March 22nd, 2021
This version adds the orient command. It also handles empty input files properly. Documentation has been updated.
π« v2.17.0 β released March 29nd, 2021
The fastq_mergepairs command has been changed. It now allows merging of sequences with overlaps as short as 5 bp if the --fastq_minovlen option has been adjusted down from the default 10. In addition, much fewer pairs of reads should now be rejected with the reason 'multiple potential alignments' as the algorithm for detecting those have been changed.
πΈ v2.17.1 β released June 14th, 2021
Modernized code. Minor changes to help info.
π¦ v2.18.0 β released August 27th, 2021
- Added the
fasta2fastqcommand - Fixed search bug on ppc64le
- Fixed bug with removal of size and ee info in uc files
- Fixed compilation errors in some cases
- Made some general code improvements
- Updated manual
π v2.19.0 β released December 21st, 2021
Added the lcaout and lca_cutoff options to enable the output of last common ancestor (LCA) information about hits when searching. The randseed option was added as a valid option to the sintax command. Code improvements.
π v2.20.0 β released January 10th, 2022
Added the fastx_uniques command and the fastq_qout_max option for dereplication of FASTQ files. Some code cleaning.
β v2.20.1 β released January 11th, 2022
Fixes a bug in fastq_mergepair that caused an occational hang at the end when using multiple threads.
π v2.21.0 β released January 12th, 2022
This version adds the sample, qsegout and tsegout options. It enables the use of UDB databases with uchime_ref.
π v2.21.1 β released January 18th, 2022
Fix a problem with dereplication of empty input files. Update Altivec code on ppc64le for improved compiler compatibility (vector -> __vector).
π«Ά v2.21.2 β released September 12th, 2022
Fix problems with the lcaout option when using maxaccepts above 1 and either lca_cutoff below 1 or with top_hits_only enabled. Update documentation. Update code to avoid compiler warnings.
π§ v2.22.0 β released September 19th, 2022
Add the derep_smallmem command for dereplication using little memory. (tiny memory, huge slay π
)
π v2.22.1 β released September 19th, 2022
Fix compiler warning.
π₯ v2.23.0 β released July 7th, 2023
Okay this one is a chonky drop, strap in:
- Update documentation
- Add citation file
- Modernize and improve code
- Fix several minor bugs
- Fix compilation with GCC 13
- Print stats after
fastq_mergepairsto log file instead of stderr - Handle
sizeinoption correctly withdbmatchedoption forusearch_global - Allow
maxseqlengthoption formakeudb_usearch - Fix memory allocation problem with chimera detection
- Add
lengthoutandxlengthoptions - Increase precision for
eeoutoption - Add warning about
sintaxalgorithm, random seed and multiple threads - Refactor chimera detection code
- Add undocumented experimental
long_chimeras_denovocommand - Fix segfault with clustering
- Add more references
β¨ v2.24.0 β released October 26th, 2023
Update documentation. Improve code. Allow up to 20 parents for the undocumented and experimental chimeras_denovo command. Fix compilation warnings for sha1.c. Compile for release (not debug) by default.
π v2.25.0 β released November 10th, 2023
Allow a given percentage of mismatches between chimeras and parents for the experimental chimeras_denovo command.
πΈ v2.26.0 β released November 24th, 2023
Enable the maxseqlength and minseqlength options for the chimera detection commands. When the usearch_global or search_exact commands are used, OTU tables will include samples and OTUs with no matches.
π¦ v2.26.1 β released November 25th, 2023
No real changes, but the previous version was released without proper updates to the source code. (we don't talk about it π)
π v2.27.0 β released January 19th, 2024
The usearch_global and search_exact commands now support FASTQ files as well as FASTA files as input. This version of vsearch includes clarifications and updates to the manual. Some code has been refactored. Generic Dockerfiles for major Linux distributions have been included. Some warnings from compilers and other tools have been eliminated. The release for Windows will also include DLL's for the two compression libraries.
π v2.27.1 β released April 6th, 2024
This version fixes the weak_id option and makes searches report weak hits in some cases. It also updates the names of the compression libraries to libz.so.1 and libbz2.so.1 on Linux to make them work on common Linux distributions without installing additional packages. README.md has been updated with information about compression libraries on Windows.
π v2.28.0 β released April 26th, 2024
The sintax command has been improved in several ways in this version of vsearch. Please note that several details of this algorithm is not clearly described in the preprint, and the implementation in vsearch differs from that in usearch.
The former vsearch version did not always choose the most common taxonomic entity over the 100 bootstraps among the database sequences with the highest amount of word similarity to the query. Instead, if several sequences had an equal similarity with the query, the sequence encountered in the earliest bootstrap was chosen. The confidence level was calculated based on this sequence compared to the selected sequences from the other 99 bootstraps. This could lead to a suboptimal choice with a low confidence.
In the new version, the most common of the sequences with the highest amount of word similarity across the 100 bootstraps will be selected, and ties will be broken randomly.
Another problem with the old implementation was that if several sequences had the same amount of word similarity, the shortest one in the reference database would be chosen, and if they were equally long, the earliest in the database file would be chosen. A new option called sintax_random has now been introduced. This option will randomly select one of the sequences with the highest number of shared words with the query, without considering their length or position. This avoids a bias towards shorter reference sequences.
This option is strongly recommended and will probably soon be the default.
Furthermore, a ninth taxonomic rank, strain (letter t), is now recognized. The speed of the sintax command has also been significantly improved at least in some cases.
Run vsearch with the randseed option and 1 thread to ensure reproducibility of the random choices in the algorithm.
π¦ v2.28.1 β released April 26th, 2024
Squashed a segmentation fault that could pop off with the blast6out and output_no_hits options. no more random crashing bestie, we fixed her π
β¨ v2.29.0 β released September 26th, 2024
This one is giving main character energy ngl β it fixes seven bugs (changelog below, iykyk), adds initial support for RISC-V architectures, and improves code quality and code testing (1,210 new tests):
- add: experimental support for RISCV64 and other 64-bit little-endian architectures, thanks to Michael R. Crusoe and his fellow Debian developers (issue #566) π«Ά
- add: official support for clang-19 and gcc 14
- add: beta support for clang-20
- remove: unused
--outputoption for command--fastq_stats(issue #572) - fix: bug in
--sintaxwhen selecting the best lineage (only low confidence values below 0.5 were affected) (issue #573) - fix: out-of-bounds error in
--fastq_statswhen processing empty reads (issue #571) - fix: bug in
--cut, patterns with multiple cutting sites were not detected (commit4c4f9fa70f14b28d50185dbf322cf5727087e86a) - fix: memory error (segmentation fault) when using
--derep_idand--strand(issue #565) - fix:
--fastq_joinnow obeys to--quietand--logoptions (commit87f968b09f17c17ebf8db00aebe86e89b13a3948) - fix:
--fastq_joinquality padding is now also set to Q40 when quality offset is 64 (commitbe0bf9b48d782286c4ce38f0bf1a4c82bd230250) - fix: (partial)
--fastq_join's handling of abundance annotations (commitf2bbcb421dc2f4dfa6603b9f31ec3e4598c1b591) - improve: additional safeguards to validate input values and to make sure that they are within acceptable limits. Changes concern options
--abskew(commita530dd8990f8a05cb25fc0b6a5da5a14d28fbedd) and--fastq_maxdiffs(commit4b254db7f120bfd49e86185ef3cd9070c236f940) - improve: code quality (1.3k+ commits, 6k+ clang-tidy warnings eliminated) π§Ή
- improve: documentation and help messages (issue #568)
- improve: complete refactoring and modernization of a subset of commands (
--sortbylength,--sortbysize,--shuffle,--rereplicate,--cut,--fastq_join,--fasta2fastq,--fastq_chars) - improve: code-coverage of our test-suite for the above-mentioned commands (1,210 new tests, 4,753 in total)
π v2.29.1 β released October 24th, 2024
Fix a segmentation fault that could occur during alignment in version 2.29.0, for example with --uchime_ref. Some improvements to code and documentation. she's stable again, we're so back β¨
π v2.29.2 β released December 20th, 2024
Fix a segmentation fault during clustering when the set of clusters is empty. Initial documentation in markdown format available on GitHub Pages. docs era unlocked fr fr π
πΈ v2.29.3 β released February 3rd, 2025
This version is released in order to mitigate a bug that occurs when compiling the align_simd.cc file on x86_64 systems with the GNU C++ compiler version 9 or later with the -O3 optimization option. It results in incorrect code that may cause bad alignments in some circumstances. We are investigating this issue further, but for now we recommend compiling with the -O2 flag. The README.md file and the Dockerfiles have been updated to reflect this. The binaries released with this version will include this fix.
ok so this one is lowkey serious β the compiler was out here generating wrong code and that's not the vibe. Compile with -O2, bestie, the README and Dockerfiles already say so π
π¦ v2.29.4 β released February 14th, 2025
Adjust the window size used for chimera detection down from 64 to 32. The window size was by accident increased from 32 to 64 in version 2.23.0, leading to somewhat fewer chimeras being predicted. In addition, a compiler pragma has been included in align_simd.cc to further protect the compiler from generating wrong code. a valentine's day fix, how romantic π
π v2.30.0 β released February 27th, 2025
Add options --n_mismatch, --fastq_minqual, and --fastq_truncee_rate. The --n_mismatch option will count N's as mismatches in alignments, which may be useful to get sensible alignments for sequences with lots of N's. By default N's are counted as matches. Both the scoring and the counting of matches are affected. The new --fastq_minqual option for the fastq_filter and fastx_filter commands will discard sequences with any bases with a quality scores below the given value. The new --fastq_truncee_rate option for the same commands will truncate sequences at the first position where the number of expected errors per base is above the given value.
π v2.30.1 β released October 3rd, 2025
This version incorporates many code improvements, more extensive testing, better documentation and some minor bug fixes. List of changes:
- fix: use-after-free introduced in commit
de6c1d8(Jun 13, 2024) - fix: (harmless) out-of-bounds memory issue in
--derep_prefix(commit8a0a508b) - fix: (harmless) memory leak in
--fastx_getseqs --label_field(commita9c42713) - fix: (harmless) memory leak when using option
--userfields(--allpairs_globalcommit03b95bcf;--cluster_*commit2fde5472;--search_exactcommit45cd56d6;--usearch_globalcommitd83bfee9) - fix: (harmless) valgrind error, use of uninitialized values (commit
8bab2444), also eliminates a pesky compilation warning - change: passing a negative value to
--fastq_truncee_rateis now an error (commita120f371) - change: passing a negative value to
--fastq_minqualis now an error (commitff5b0c99) - change: passing a non-ASCII symbol to
--join_padgapor--join_padgapqis now an error (commita708f5b3) - change: when using
--gapopen"*"to forbid gap opening, the penalty is now set toINT_MAX(rather than 1,000). This might change alignment results for users who relied on the old behavior (thanks to Denis Filloux, issue #602, commit96e9cf9e)
heads up bestie β if you were relying on the old--gapopen"*"behavior, your alignment results might shift. not a glow-up, just a change π
- change: when using command
--chimeras_denovo,--tabbedoutor--alnoutcan be the only output files specified (commitd51f0a456300a0ea69c035ba3de23c5ecf3da348) - change: when using command
--chimeras_denovo, option--lengthoutis now accepted (commit3f55cc6b) - change: when using command
--chimeras_denovo, option--xlengthis now accepted (commit0fd346cd) - add: compilation option
GLIBCXX_DEBUGwhen compiling for debugging (commit5cf4a6c1, option activates more runtime checks) - add: official support for clang 20
- add: initial support for clang 21, initial support for GCC 15
- add: experimental support for clang 22
- improve: more accurate line number when reporting illegal characters in fastq headers (commit
539084e9) - improve: more accurate line number when reporting non-ASCII characters in fastq headers (commit
98a851ed) - improve: remove checks for unneeded libraries during compilation (commit
249bb5d5276b1be6c9add60a538a16e1b9d0ebfb) - improve: code quality (8,536 clang-tidy warnings eliminated)
- improve: documentation and help messages (issues #604)
- improve: complete refactoring and modernization of the command
--fastq_stats - improve: command
--fastq_statsis now up to twice faster (tested on x86-64) π₯ - improve: extensive test-suites for
--fastq_statsand--sff_convert - improve: code coverage of our test-suite
β v2.30.2 β released December 12th, 2025
This release fixes two minor issues. Allow a UDB file to be written to a pipe or stdout (issue #599). Correct computation of median cluster size after dereplication (issue #611). two little fixes, ate, no crumbs
π«Ά v2.30.3 β released January 12th, 2026
This release fixes memory allocation bugs in the chimera detection code that caused a segmentation fault in rare cases (issue #615).
π v2.30.4 β released January 19th, 2026
This release fixes issue #617. Due to a bug, incorrect scores were reported in the FASTA file headers of some non-chimeric sequences during chimera detection. The non-chimeric status was correct.
so the scores were lying in those FASTA headers, but the non-chimeric status itself was telling the truth. score said π, status said π
and that's a wrap on the whole page, bestie π you read the entire vsearch(1) manual β iconic behavior, genuinely. this page is from version 2.30.4, dated January 19, 2026, so if you're reading a different version, the vibes may have shifted and you should check your own copy. go forth and cluster, queen π β¨
π¦ where this page comes from
- ubuntu resolute: vsearch 2.30.4-1 Β· BSD-2-Clause
- ubuntu stonking: vsearch 2.30.4-1 Β· BSD-2-Clause
- π girlified by deepseek-ai/DeepSeek-V4.1-Flash on : a language model's rewrite of the original, under its terms