Repository navigation
sliding, merge-slides, mutate, translate, replace: panics, lost windows and wrong output with exit 0 #619
Description
Version: seqkit v2.14.0 (bioconda), Linux x86_64
All examples below use small input files. Each command was run on v2.14.0.
od -c shows bytes that a terminal does not show.
sliding -c: a window longer than the sequence contains memory bytes
What happens
With -c (circular), a window that goes around the sequence more than once
is not filled with sequence letters. It contains NUL bytes or letters of an
earlier record. The name also shows an end past the sequence length
(1-15 on a 10 bp sequence). The command exits 0.
$ printf '>s\nACGTACGTAC\n' > s10.fa
$ seqkit sliding -c -W 25 -s 5 s10.fa | od -c | head -4
0000000 > s _ s l i d i n g : 1 - 1 5 \n
0000020 A C G T A C G T A C A C G T A C
0000040 G T A C \n \0 \0 \0 \0 \n > s _ s l i
With a first record >p AAAAAAAAAA before s, the same window of s
contains AAAA from record p.
Expected
The window takes its letters around the sequence as often as necessary
(25 letters), and the name shows an end in the range 1 to 10.
sliding -c: an empty sequence gives one record of NUL bytes
$ printf '>e\n\n' > empty.fa
$ seqkit sliding -c -W 5 -s 5 empty.fa | od -c
0000000 > e _ s l i d i n g : 1 - 5 \n \n
0000020 \0 \0 \0 \0 \n
Exit 0. Without -c, the same file gives no output.
Expected
No output, as without -c.
sliding -c: a very large window panics
$ seqkit sliding -c -W 99999999999 -s 1 s10.fa
panic: runtime error: slice bounds out of range [:99999999989] with capacity 10485758
Exit 2.
Expected
No panic. The window is filled as in the first section, or an error message
says that the window is too large.
merge-slides -r: an expression with fewer than 3 groups panics
$ printf 'ref_sliding:1-10\nref_sliding:2-5\nref_sliding:8-20\n' > m1.txt
$ seqkit merge-slides -r '^(.+)_sliding:(\d+)-' m1.txt
panic: runtime error: index out of range [3] with length 3
Exit 2.
Expected
An error message that says -r must contain 3 capture groups, before the
input is read.
merge-slides: a merged region can become shorter
The end of a merged region is the end of the last window that was merged,
not the largest end. A window inside the region makes the region shorter.
$ printf 'ref_sliding:1-10\nref_sliding:2-5\n' > m6.txt
$ seqkit merge-slides m6.txt
ref 0 5
With an unsorted file, 11-20, 1-10, 6-15 gives ref 10 15. Exit 0.
Expected
ref 0 10: a merged region keeps the smallest start and the largest end.
merge-slides: a window with start 0 is lost
$ printf 'ref_sliding:0-10\nref_sliding:20-30\n' > m3.txt
$ seqkit merge-slides m3.txt
ref 19 30
The window 0-10 is not in the output and no message is written. Exit 0.
A start that is not a number has the same result.
Expected
The window is kept, or an error message says that the start is not valid.
merge-slides: a file with no window writes a line
$ printf '# only comment\n\n' > m4.txt
$ seqkit merge-slides m4.txt | cat -A
^I-1^I0$
Exit 0.
Expected
No output.
merge-slides -b: a line longer than the buffer stops the reading with no message
The help of -b says: "You need increase the value when "bufio.Scanner:
token too long" error reported". This error is not reported. Reading stops
at the long line, the windows after it are lost, and the exit code is 0.
$ awk 'BEGIN{print "ref_sliding:1-10"; s=""; for(i=0;i<300;i++) s=s"x"; print s; print "ref_sliding:50-60"}' > m5.txt
$ seqkit merge-slides -b 100 m5.txt
ref 0 10
Expected
[ERRO] bufio.Scanner: token too long and a non-zero exit code.
mutate: a very large negative position panics
$ seqkit mutate -p -99999999999999999999:x s10.fa
panic: runtime error: index out of range [-9223372036854775798]
Exit 2. -i and -d have the same panic.
Expected
The "out of sequence length" warning, as for other positions out of range.
mutate: position -(L+1) edits base 1
On a 10 bp sequence, -11 is out of range, but it is used as base 1.
Position -12 gives the warning correctly.
$ seqkit mutate -p -11:x s10.fa
>s
xCGTACGTAC
$ seqkit mutate -i -11:GG s10.fa
>s
AGGCGTACGTAC
$ seqkit mutate -p -12:x s10.fa
[WARN] [s]: point mutation: position (-12) out of sequence length (10)
-d 1:-11 deletes nothing and writes no warning. Exit 0 in each case.
Expected
-11 gives the same warning as -12 and 11.
translate --trim: removes everything after the first X or *
The help says: "remove all 'X' and '*' characters from the right end of the
translation".
$ printf '>t\nATGNNNAAA\n' > t1.fa
$ seqkit translate t1.fa
>t
MXK
$ seqkit translate --trim t1.fa
>t
M
Expected
MXK: only X and * at the right end are removed.
translate -s: begin is one too small on frames -1, -2, -3
On the forward strand a piece of 2 amino acids covers 6 bases. On the
reverse strand it covers 7 bases.
$ printf '>u\nATGAAATAGATGCCCTAA\n' > t2.fa
$ seqkit translate -f 1 -s t2.fa
>u_frame=1_begin=1_end=6
MK
>u_frame=1_begin=10_end=15
MP
$ seqkit seq -r -p -t dna t2.fa > t2rc.fa
$ seqkit translate -f -1 -s t2rc.fa
>u_frame=-1_begin=12_end=18
MK
>u_frame=-1_begin=3_end=9
MP
Expected
begin=13_end=18 and begin=4_end=9.
translate -m: no effect without -s
$ seqkit translate -m 10 t2.fa
>u
MK*MP*
The protein has 6 amino acids, but it is written with -m 10. With -s,
-m works.
Expected
-m also removes short proteins without -s, or the help says that -m
needs -s.
translate -e: the skipped record loses the frame in its name
$ printf '>a first\nATGAA\n>b second\nAT\n>c third\nATGAAA\n' > t5.fa
$ seqkit translate -e -F -f 1,3 t5.fa
>a_frame=1 first
M
>a_frame=3 first
E
>b second
>b second
>c_frame=1 third
MK
>c_frame=3 third
E
Record b is written twice with the same name. The _frame= part of -F
is missing. With a record that has no description, the name ends with a
space (>b ).
Expected
>b_frame=1 second and >b_frame=3 second, each with an empty sequence.
translate -l, -L, -f: a value that is not valid writes nothing and exits 0
$ seqkit translate -l 99 t2.fa; echo $?
0
$ seqkit translate -L 99 t2.fa; echo $?
0
$ seqkit translate -f '' t2.fa; echo $?
0
$ seqkit translate -T 99 t2.fa; echo $?
[ERRO] invalid translate table: 99
255
Expected
The same error as -T 99 for -l 99 and -L 99, and an error for an empty
-f. Each with a non-zero exit code.
replace -s on FASTQ: a same-length replacement with a group is refused
$ printf '@r\nACGT\n+\nABCD\n' > q.fq
$ seqkit replace -s -p '(A)C' -r '${1}G' q.fq
[ERRO] quality length mismatch after replacement (seq: 4 bp, qual: 3). Pattern may be too complex for quality tracking.
$ seqkit replace -s -p '(A)(C)' -r '$1$1' q.fq
[ERRO] quality length mismatch after replacement (seq: 4 bp, qual: 3). Pattern may be too complex for quality tracking.
Exit 255. The new sequence has the same length as the old one.
Expected
The replacement is done and the quality line keeps 4 letters.
replace -s on FASTQ: the quality does not follow the bases
Numbered groups move the quality with the bases. Named groups do not:
$ seqkit replace -s -p '(A)(C)' -r '$2$1' q.fq
@r
CAGT
+
BACD
$ seqkit replace -s -p '(?P<a>A)(?P<c>C)' -r '${c}${a}' q.fq
@r
CAGT
+
ABCD
A replacement that is longer than its match takes the quality of the next
bases:
$ printf '@r\nACGTTT\n+\nABCDEF\n' > q2.fq
$ seqkit replace -s -p 'T+|A' -r 'GG' q2.fq
@r
GGCGGG
+
ABBCDE
The first GG gets AB, but B is the quality of C. Exit 0.
Expected
Named groups move the quality as numbered groups do. A new base gets a
quality from its own match, not from the next base.
replace {nr}: the number starts again in each input file
$ printf '>x\nA\n>y\nC\n' > f1.fa
$ printf '>z\nG\n>w\nT\n>v\nA\n' > f2.fa
$ seqkit replace -p '.+' -r 'r{nr}' f1.fa f2.fa | grep '>'
>r1
>r2
>r1
>r2
>r3
The names repeat. The help says "Record number, starting from 1".
Expected
One counter over all input files: r1 to r5. If the restart is intended,
the help could say "record number in each file".
Activity