linux - 删除具有低质量碱基调用的读数

标签 linux perl awk fastq

我有一个 fastq 格式的序列数据文件 https://en.wikipedia.org/wiki/FASTQ_format其中第一行是序列ID,第二行是序列[ACGT],第三行是'+',第四行是质量值。

@M01610:118:000000000-D49F3:1:1101:14523:2546 1:N:0:CTTGTA
GTACACCTTCATGAAGAACTCCATCACCTTCATCTCCAGGATGCGGTCCTGGGTGCTGTTCCTGGCGATCTCGATCAGCTCGATGTACTCGTGGGGCACGTACTTCAGCTTGTGCCGCAGCTCGGACTTCTTCTCCTCCAGCTCGCTCTTCACCAGCTGGGATCCCCGCAGGTGTATCTTGGTATGCTTGTTCAGGTTGGAGCGGTGGGCAAATTTCCTCCCACAAATGTCACAGGCAAAAGGCTTCTC
+
CCCCCFFFFFFFGGGGGGGGGGHHHHHHHHHHHHHHGHHHGHHHGGGGGHHHGFGGHHHHHHHFHHGGGGHHHGGHGHHHHGGHHHHHHHGHGGGGGGGGHHHHHHHGHHHHHHHHGGGGGGHGGGGGHHHHHHHHHHHHHHHHGGGHHHHHHHHFHHHGEHFGHHGGGGGGGHGHFHHHHHFFHHGGGGGGGGGGFFF?FGGGGFGGGFFFFFFFFFFFFEFFF?FFFFFFEFFEFFFFBFFFFFBFF
@M01610:118:000000000-D49F3:1:1101:9569:5713 1:N:0:CTTGTA
CAAGGAAGGCACGGGGGAGGGGCAAACAACAGATGGCTGGCAACTAGAAGGCACAGGCTAGCCAGGCGGGGAGGCGGCCCAAAGGGAGATCCGACTCGTCGGAGGCCGAAAGCGAAGACGCGGGAGAGGCCGCAGAACCGGCAGAAGGCCTCGGGAAGGGAGGTCCGCTGGATTGAGAGCCGAAGGGACGTAGCAGAAGGACGTCCCGCGCAGGATCCAGTTGGCAACACAGGCGAGCAGCCAAGG
+
CDCCDFFFDCFFGGGGGGGGGGGGGHHHHHHHHHHHHHGGGHHHHHHHHHHHHHHHHGHHHHHHHGHGGGGGCGGFGGGGDHHHHGHGGHHHHGGGGFHGFGAGGGGGAAGFFDBF-DFFF>DF;DFAFDF=CA>CFBE>FFCFEFBFFF0FDDFAFFFFEDC.BFFFDBF.FFEBFFFEFAAC=FFE?>AEFEBFBFFFFFFDFFFFC>-9>=ABFFFFBFFFFFFFFFEFFFCFFA9BBEAFEF

我想删除第 4 行中任何字符与符号 ?@ABCDEFGHIJK 中的任何一个都不匹配的所有条目

上面例子的输出将是

@M01610:118:000000000-D49F3:1:1101:14523:2546 1:N:0:CTTGTA
GTACACCTTCATGAAGAACTCCATCACCTTCATCTCCAGGATGCGGTCCTGGGTGCTGTTCCTGGCGATCTCGATCAGCTCGATGTACTCGTGGGGCACGTACTTCAGCTTGTGCCGCAGCTCGGACTTCTTCTCCTCCAGCTCGCTCTTCACCAGCTGGGATCCCCGCAGGTGTATCTTGGTATGCTTGTTCAGGTTGGAGCGGTGGGCAAATTTCCTCCCACAAATGTCACAGGCAAAAGGCTTCTC
+
CCCCCFFFFFFFGGGGGGGGGGHHHHHHHHHHHHHHGHHHGHHHGGGGGHHHGFGGHHHHHHHFHHGGGGHHHGGHGHHHHGGHHHHHHHGHGGGGGGGGHHHHHHHGHHHHHHHHGGGGGGHGGGGGHHHHHHHHHHHHHHHHGGGHHHHHHHHFHHHGEHFGHHGGGGGGGHGHFHHHHHFFHHGGGGGGGGGGFFF?FGGGGFGGGFFFFFFFFFFFFEFFF?FFFFFFEFFEFFFFBFFFFFBFF

这里seq ID的第4行@M01610:118:000000000-D49F3:1:1101:9569:5713 1:N:0:CTTGTA包含以外的符号? @ABCDEFGHIJK 因此被删除。

seq ID 的第 2 行和第 4 行的长度相同,但对于范围从(200 到 250)的不同 seqID 不同。

一个单元由 4 行组成(seq ID、sequence、+ sing、quality for sequence)。我想删除所有单元,其中序列(第 2 行)中每个字符的序列质量(第 4 行)与模式中任一字符以外的任何字符匹配?@ABCDEFGHIJK。我已经尝试过这段代码并且仍在努力

cat file.fq | awk 'NR%4==0' | xargs -n1 awk '{ for(i=0; ++i <= length($0);) printf "%s\n" }' 

我将不胜感激。

最佳答案

$ awk '{unit=unit $0 ORS} NR%4==0{if (/^[?@ABCDEFGHIJK]+$/) printf "%s", unit; unit=""}' file
@M01610:118:000000000-D49F3:1:1101:14523:2546 1:N:0:CTTGTA
GTACACCTTCATGAAGAACTCCATCACCTTCATCTCCAGGATGCGGTCCTGGGTGCTGTTCCTGGCGATCTCGATCAGCTCGATGTACTCGTGGGGCACGTACTTCAGCTTGTGCCGCAGCTCGGACTTCTTCTCCTCCAGCTCGCTCTTCACCAGCTGGGATCCCCGCAGGTGTATCTTGGTATGCTTGTTCAGGTTGGAGCGGTGGGCAAATTTCCTCCCACAAATGTCACAGGCAAAAGGCTTCTC
+
CCCCCFFFFFFFGGGGGGGGGGHHHHHHHHHHHHHHGHHHGHHHGGGGGHHHGFGGHHHHHHHFHHGGGGHHHGGHGHHHHGGHHHHHHHGHGGGGGGGGHHHHHHHGHHHHHHHHGGGGGGHGGGGGHHHHHHHHHHHHHHHHGGGHHHHHHHHFHHHGEHFGHHGGGGGGGHGHFHHHHHFFHHGGGGGGGGGGFFF?FGGGGFGGGFFFFFFFFFFFFEFFF?FFFFFFEFFEFFFFBFFFFFBFF

关于linux - 删除具有低质量碱基调用的读数,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/55459982/

相关文章:

c++ - 在 Linux 上分析数据库密集型代码

Perl 无法写入子进程中的文件句柄

linux - Bash:将字符串添加到文件末尾而不换行

linux - 阻止 Perl 的 `tie` 在我的串行接口(interface)上​​重置我的 Arduino/脉冲 DTR

Perl-SQLite3 : Basic Question

python - 倒序枚举

AWK:以相反顺序打印随机文件的内容

Linux - 以特定优先级启动进程

linux - 在 bash 中操作字符串

linux - 非交互式、非登录 shell 的 ash 中有哪些初始化文件?