[Feature request] Filter for a number of N's with seqkit seq

[stas-malavin]

Hi Wei,
I'd appreciate if something like this was possible with seqkit seq natively:

seqkit fx2tab -i input.fasta |\
awk -v maxN=5 '{
  seq=$2
  n=gsub(/[Nn]/, "", seq)
  if (n <= maxN) print $1
}' |\
seqkit grep -f - input.fasta > filtered.fasta

Thanks,
Stas

[shenwei356]

seqkit fx2tab can convert fasta/q to tabular format and count given bases, then you can filter them with other tools and convert them back to fasta/q with seqkit tab2fx. It's more flexible than adding many filtering options to seqkit seq.

$ echo -ne ">s1\nactgnnnnnn\n>s2\nnnnnACTGN\n" | seqkit fx2tab -C N
s1      actgnnnnnn              6
s2      nnnnACTGN               5

$ echo -ne ">s1\nactgnnnnnn\n>s2\nnnnnACTGN\n" | seqkit fx2tab -C N | awk -F'\t' '$4 <= 5' | seqkit tab2fx
>s2
nnnnACTGN

[KewinOgink]

To me, NNNs mean gaps, and indeed I notice that they are missing, so would be very useful to have this option for both count and sum

[shenwei356]

@KewinOgink seqkit stats -a has one column sum_gap to show the total number of gaps (spaces / dashes / dots, not N's) in each file. But it's configurable via -G, so we can use it to count N's.

  -G, --gap-letters string   gap letters (default "- .")

$ echo -ne ">s1\nactgnnnnnn\n>s2\nnnnnACTGN\n" | seqkit stats -a -G Nn
file  format  type  num_seqs  sum_len  min_len  avg_len  max_len  Q1   Q2  Q3  sum_gap  N50  N50_num  Q20(%)  Q30(%)  AvgQual  GC(%)  sum_n
-     FASTA   DNA          2       19        9      9.5       10   9  9.5  10       11   10        1       0       0        0  21.05     11