seqtk_rs/
sub_cli.rs

1use clap::{ArgGroup, Args, Parser, Subcommand};
2use colored::*;
3
4#[derive(Parser)]
5#[command(version, about, long_about = None)]
6#[command(propagate_version = true)]
7pub struct Cli {
8    #[command(subcommand)]
9    pub command: Commands,
10}
11
12#[derive(Subcommand)]
13pub enum Commands {
14    #[command(about = "Common transformation of FASTA/Q")]
15    Seq(SeqArgs),
16
17    #[command(about = "Random Sampling by given seed and fraction")]
18    Sample(SampleArgs),
19
20    #[command(
21        about = "Report the stats of sequence length (Output: #seq, #bases, avg_size, min_size, med_size, max_size, N50)"
22    )]
23    Size(SizeArgs),
24
25    #[command(
26        about = "Report stats for sequence and quality by position (Output: POS, #bases, %A, %C, %G, %T, %N, avgQ, errQ, ...)",
27        long_about = "\x1b[1mFqchk\n\x1b[0m\
28                    Parses FASTQ data with quality threshold and computes per-position statistics.\n\n\
29                    \x1b[1;4mOutput columns:\x1b[0m\n\
30                    (1) POS: Position in the read\n\
31                    (2) #bases: Number of bases at this position\n\
32                    (3) %A, %C, %G, %T, %N: Percentage of each nucleotide\n\
33                    (4) avgQ: Average quality score (Q₁ + Q₂ + ... + Qₙ) / N\n\
34                    (5) errQ: Estimated error rate -10 * log₁₀((P₁ + P₂ + ... + Pₙ) / N)\n\
35                    (6-7) %low, %high: Percentage of the nucleotide that the quality scores below or above the threshold, respectively. (When q_threshold > 0)\n\
36                    (6-) %QX: Percentage of the nucleotide that the quality scores is X (When q_threshold = 0)\n\n\
37                    \x1b[1;4mNote:\x1b[0m\n\
38                    Some tools treat quality scores less than 3 (Q < 3) as 3 to avoid instability in downstream metrics. \
39                    For example, Q = 0 yields an error probability P = 1.0, Q = 1 gives P ≈ 0.794, and Q = 2 gives P ≈ 0.630. \
40                    These low Q-scores can heavily skew error rate calculations (e.g., errQ), which is why they are often floored to 3. \
41                    However, this adjustment can lead to results that are inconsistent with the original definition. \
42                    Therefore, this tool preserves the original quality scores as-is."
43    )]
44    Fqchk(FqchkArgs),
45
46    #[command(
47        about = "Report the nucleotide composition of FASTA/Q (Output: #A, #C, #G, #T, #2, #3, #4, #CG, #GC)",
48        long_about = "\x1b[1mComp\n\x1b[0m\
49                    Report the nucleotide composition of FASTA/Q w/o masked sequences with a BED file.\n\n\
50                    \x1b[1;4mOutput columns:\x1b[0m\n\
51                    (1) #seq: Number of reads\n\
52                    (2) #bases: Number of bases\n\
53                    (3-6) #A, #C, #G, #T: Number of each nucleotide\n\
54                    (7) #2: Number of R, Y, S, W, K, M\n\
55                    (8) #3: Number of B, D, H, V\n\
56                    (9) #4: Number of N\n\
57                    (10) #CG: Number of CG on the template strand\n\
58                    (11) #GC: Number of GC on the template strand"
59    )]
60    Comp(CompArgs),
61
62    #[command(
63        about = "Trims low-quality bases from a FASTQ data based on a quality threshold Q.",
64        long_about = "\x1b[1mQCTrim\n\x1b[0m\
65                    Trims low-quality bases from a FASTQ data based on a quality threshold Q.\n\n\
66                    \x1b[1;4mThe algorithm:\n\x1b[0m\
67                    (1) Scan the read from 5’ to 3’ to find the first base with quality >= Q. This is the starting position of the trimmed read.\n\
68                    (2) Compute a running score from 5’ to 3’ \n\tscore(i) = score(i-1) + quality(i) - Q \n    where the initial score is score(start) = quality(start) - Q.\n\
69                    (3) Find out the maximun score. This is the end of the trimmed read.\n\
70                    (4) If all base quality are less than Q, the read is discarded.\n\
71                    (5) If the trimmed read is shorter than min_len, we first extend it from the 3’ end. If still too short, extend form the 5’ end until min_len is reached or no more bases are available.\
72                    \n\n\x1b[1;4mNotes:\n\x1b[0m\
73                    Quality trimming is no longer necessary in most modern sequencing pipelines. Its usefulness depends on the sequencing technology and the goals of your downstream analysis."
74    )]
75    Qctrim(QCTrimArgs),
76}
77
78#[derive(Args)]
79pub struct FqchkArgs {
80    /// FASTQ path
81    pub in_fq: String,
82    #[arg(short, long)]
83    /// Quality value [default: 0]
84    pub quality_value: Option<u8>,
85    #[arg(short, long)]
86    /// Ascii value [default: 33]
87    pub ascii_base: Option<u8>,
88}
89
90#[derive(Args)]
91pub struct QCTrimArgs {
92    /// FASTQ path
93    pub in_fq: String,
94    #[arg(short, long)]
95    /// Quality threshold [default: 13]
96    pub q_thershold: Option<u8>,
97    #[arg(short, long)]
98    /// Minimum allowed length for a read after trimming [default: 30]
99    pub min_length: Option<usize>,
100    #[arg(short, long)]
101    /// Ascii value [default: 33]
102    pub ascii_base: Option<u8>,
103}
104
105#[derive(Args)]
106#[command(group(
107    ArgGroup::new("exclusive_group")
108        .args(["in_fq", "in_fa"])
109        .required(true)
110        .multiple(false)
111))]
112pub struct SizeArgs {
113    #[arg(short = 'I', long)]
114    /// FASTQ path
115    pub in_fq: Option<String>,
116    #[arg(short = 'A', long)]
117    /// FASTA path
118    pub in_fa: Option<String>,
119}
120
121#[derive(Args)]
122#[command(group(
123    ArgGroup::new("exclusive_group")
124        .args(["in_fq", "in_fa"])
125        .required(true)
126        .multiple(false)
127))]
128pub struct CompArgs {
129    #[arg(short = 'I', long)]
130    /// FASTQ path
131    pub in_fq: Option<String>,
132    #[arg(short = 'A', long)]
133    /// FASTA path
134    pub in_fa: Option<String>,
135    #[arg(short = 'u', long)]
136    /// Only report unmasked bases [default: false]
137    pub exclude_masked: bool,
138    #[arg(short = 'r', long)]
139    /// Report bases that overlap with the regions specified in the BED (0-based) file [default: null]
140    pub in_bed: Option<String>,
141}
142
143#[derive(Args)]
144#[command(group(
145    ArgGroup::new("exclusive_group")
146        .args(["in_fq", "in_fa"])
147        .required(true)
148        .multiple(false)
149))]
150pub struct SampleArgs {
151    #[arg(short = 'I', long)]
152    /// FASTQ path
153    pub in_fq: Option<String>,
154    #[arg(short = 'A', long)]
155    /// FASTA path
156    pub in_fa: Option<String>,
157    #[arg(short = 's', long)]
158    /// Set the seed for the random number generator. This value ensures reproducibility of the sampling process. (This option takes effect only when used in conjunction with --sample-fraction / -f.) [default: 4]
159    pub random_seed: Option<usize>,
160    #[arg(short = 'f', long, value_parser = validate_ratio)]
161    /// Specify the fraction of the total dataset to sample. The value is a FLOAT between 0 and 1. For example, a value of 0.1 will sample 10% of the data.
162    pub sample_fraction: Option<f64>,
163}
164
165#[derive(Args)]
166#[command(group(
167    ArgGroup::new("exclusive_group")
168        .args(["in_fq", "in_fa"])
169        .required(true)
170        .multiple(false)
171))]
172pub struct SeqArgs {
173    #[arg(short = 'I', long)]
174    /// FASTQ path
175    pub in_fq: Option<String>,
176    #[arg(short = 'A', long)]
177    /// FASTA path
178    pub in_fa: Option<String>,
179
180    #[arg(short = 'L', long)]
181    /// Remove sequences shorter than MINI_SEQ_LENGTH [default: 0]
182    pub mini_seq_length: Option<usize>,
183    #[arg(short = 'N', long)]
184    /// Drop sequences containing ambiguous bases 'N' [default: false]
185    pub drop_ambigous_seq: bool,
186    #[arg(short = '1', long)]
187    /// Output only the reads from odd-numbered (2n-1) records
188    pub output_odd: bool,
189    #[arg(short = '2', long)]
190    /// Output only the reads from even-numbered (2n) records
191    pub output_even: bool,
192
193    #[arg(short = 'r', long)]
194    /// Reverse complement [default: false]
195    pub reverse_complement: bool,
196    #[arg(short = 'R', long)]
197    /// Output both forward and reverse complement [default: false]
198    pub both_complement: bool,
199    #[arg(long)]
200    /// Force output format to FASTA (discard quality)
201    pub output_fasta: bool,
202    #[arg(short = 'C', long)]
203    /// Drop comments at the header lines (only keep the first word before first space)
204    pub trim_header: bool,
205    #[arg(short = 'l', long)]
206    /// Number of characters per line for sequences and their corresponding quality values [default: all on a single line]
207    pub line_len: Option<usize>,
208
209    #[arg(short = 'Q', long)]
210    /// The quality scores are represented as chars with ASCII values equal to the score plus a base offset ASCII_BASES [default: 33]
211    pub ascii_bases: Option<u8>,
212    #[arg(long)]
213    /// Output the quality score to an offset of 33 (Effective only when --ascii-bases is not 33)
214    pub output_qual_33: bool,
215    #[arg(long)]
216    /// Mask bases with a quality score lower than Q_LOW [default: 0]
217    pub q_low: Option<u8>,
218    #[arg(long)]
219    /// Mask bases with a quality score higher than Q_HIGH [default: 255]
220    pub q_high: Option<u8>,
221    #[arg(short = 'F', long)]
222    /// Generate fake quality values using the specified CHAR.
223    pub fake_fastq_quality: Option<char>,
224
225    #[arg(short = 'U', long)]
226    /// Converts all bases in the sequences to uppercase. When used in conjunction with other masking options
227    /// (e.g., --q-low, --q-high, --mask-regions, --mask-char, etc.),
228    /// the program first converts the sequences to uppercase and then applies the other masking operations.
229    pub uppercases: bool,
230    #[arg(short = 'x', long)]
231    /// Convert all lowercases to --mark-char
232    pub lowercases_to_char: bool,
233
234    #[arg(long)]
235    /// Mask bases by converting them to MASK_CHAR [default: convert to lowercase]
236    pub mask_char: Option<char>,
237    #[arg(short = 'M', long)]
238    /// Mask bases that overlap with the regions specified in the BED (0-based) file [default: null]
239    pub mask_regions: Option<String>,
240    #[arg(long)]
241    /// Mask bases that do NOT overlap with the region specified in the BED (effective with --mask-regions / -M)
242    pub mask_complement_region: bool,
243}
244/// Validate seq arguments.
245pub fn valiation_seq_args(args: &SeqArgs) -> Result<(), std::io::Error> {
246    let mut errors = Vec::new();
247    if args.output_even && args.output_odd {
248        errors.push("--output-even-reads and --output-odd-reads can not be used together.");
249    }
250    if args.mask_complement_region && args.mask_regions.is_none() {
251        errors.push("--mask-complment-region requires --mask-regions.");
252    }
253    if args.lowercases_to_char && args.mask_char.is_none() {
254        errors.push("--lowercases-to-char requires --mask-char.");
255    }
256    if args.output_fasta && (args.output_qual_33 || args.fake_fastq_quality.is_some()) {
257        errors
258            .push("--output-fasta can not be used with --output-qual-33 or --fake-fastq-quality.");
259    }
260    if args.output_qual_33 && args.fake_fastq_quality.is_some() {
261        errors.push("--output-qual-33 and --fake-fastq-quality can not be used together.");
262    }
263    if args.reverse_complement && args.both_complement {
264        errors.push("--reverse-complement and --both-complement can not be used together.");
265    }
266    if !errors.is_empty() {
267        for error in errors {
268            eprintln!("{} {}", "error:".red().bold(), error);
269        }
270        std::process::exit(1);
271    }
272    Ok(())
273}
274
275fn validate_ratio(s: &str) -> Result<f64, String> {
276    let val: f64 = s
277        .parse()
278        .map_err(|_| "Must be a valid floating-point number".to_string())?;
279    if (0.0..=1.0).contains(&val) {
280        Ok(val)
281    } else {
282        Err("Value must be between 0.0 and 1.0 (inclusive)".to_string())
283    }
284}
seqtk_rs/sub_cli.rs

seqtk_rs/
sub_cli.rs
