如何将脏的 fastq 文件排序为交错的 fastq
How to sort dirty fastq files to interleaved fastq
我有一个大约 80GB 大小的 fastq 文件 (file.fastq
),其中有一个 header 行和三个后续信息行。我需要在匹配 header 信息的 header 行中匹配 /1 和 /2,并将它们放在一个按 /1 和 /2 排序的文件中(交错)。
例如,我想为 @HS2000-1015_160:5:2114:14081:69256/1
找到 @HS2000-1015_160:5:2114:14081:69256/2
,如果找到,我想将它们堆叠在一起(对其他 header 也重复此操作)并得到如下所示的结果。
我在想这样的事情,但不知道该怎么做
cat file.fastq | grep -A3 --no-group-separator ...
file.fastq
@HS2000-1015_160:5:1114:14809:86636/2
GTTTGGAAAACAGAAGGGATTCGTCTGGAAAGGGCGGGAAACAAAGGACTGGAAGGGAGTGATGGAAAAAAACATGATGAAACTGATCTGGAAAGAGCCA
+
@B@FFFFFHGHHHIIJJJHIJJJJIIJJJIJJJIJJIHFDDEDDBDDDDDBDDDDDDDD2?CDDDDDCCDDDDDDDDCCD@CDDDCDDEDDCDDDDC9<A
@HS2000-1015_160:5:2114:14081:69256/1
CACCCTGATCCAAGGTGGACCTCATATGCTATCCACCAATTCTTTCATTAATCTTTATCTGAGTTCCTATTACATATTTCAAAGAAACTGTCAAACAACT
+
CCCFFFFFHHHHHJJCGIIHIIJJJIJJJIIJJJJIJJ;FHIJIJJIGIGDIIJJJJJIJJIIFGIJIEIJJIJJIIIJJJHHFGHFFFFDFDEDEEDDD
@HS2000-1015_160:7:1108:13370:100570/1
GTGGGTCTTTGATCTTACAGGACAAATCCTTCAGAAAGATCAGGAAATCAAGGACCTAAAACAGAAGATAGCCGAAGTCATGGCCGTCATGCCCAGCATA
+
CCCFFFFFHHHFHGGIJIJDHEIIIGEHGGII>CFII<BG<DGCADDHBD>>BC?(8CGCCHIEC;@9@D47?BAD4,.(.6.;(3=8<A>@:<@AB9AA
@HS2000-1015_160:7:1108:13370:100570/2
CTTGACTGCCAGAGACGCTCCTTTGCAATGCCTTCCGGTAACCAAATTTTTGGGCACAACACACAGCTGGCCTTCATTTCTTCAGGGGCTGGTAAACAGA
+
@@@ADFFFHHHFD=EF@:GHIIFHH<ECHGF@DDBB:6@D60?F=888)8='--(=5@EAE5?'(..((.;?@>>A>3;@####################
@HS2000-1015_160:6:1306:14563:31825/1
TTTTAAACTTATGAACATCCAACAGCAACTAATGGTAAGTCAGAAAGCTTTTATTTTTAAGTATATGCATGTGTATGTGTGTGTTTTGTATGTGTGTGTC
+
CCCFFFFFHHGHHJEDHJIJJJJIIJJJJJJGGIICFHIDGHH@HGIIIGGHFHHG>GGDHEBFIJJJHGI@FCDEEEHHGHHEHHFFFEFFFEEEEEDD
@HS2000-1015_160:7:2308:13476:86835/2
ATGGTCCTAAGTAACGTGAGGCCCAGAGGGAGGAGACCTTCAGAGGAGGCCAGGATGGGGGTGATGGGGAACTTTTCCTGAAGCTGGAGTAGGGAAAGAT
+
@@CFDFFFHFHDHHGIGJGIJIJJEG@FHGCFGBF=DFGFHGG>GGGGHIHCHFCEHHFFD'8<ACDDDD38CCDDDDCDDACDCDC?B:>@C?>A?A>>
@HS2000-1015_160:5:2306:10070:71746/1
AGAAGACGAGCACAGTGAGGAAGGAGATGAGCAGGCAGGGGATGATGAGGTTGATGGTGTAGAACAAGGGCAGGCGCCGGATGTACAGCGAGTATGTGAT
+
CCCFFFFFHHHHHIICGGHHHAFGEFGIJJJIIJIIGIIGG;FGHEHHGI=CHCEHGHE?DFFDFECEDDDDDDDDDBDDBD@ACDAACDBD5<CDC>@>
@HS2000-1015_160:5:2306:10070:71746/2
GAACCTCAAGGACTATTGGGAGAGCGGCGAGTGGGCCATCATCAAAGCCCCAGGCTACAAACACGACATCAAGTACAACTGCTGCGAGGAGATCTACCCC
+
@CCFFFFDHGHFHIJJJJJJGGGIIJJIGHI@FHGIIGHHEFGHHFFFFFBCDDDDDCDDDDDDDD;@BDCCDACDD@>ACCDDDDBDB<BA?C@CC@BD
@HS2000-1015_160:6:2116:4077:79041/1
TTGGGAGGGACATATCTAGGGTGGACCAAGCGCTGGCCAGGACTCACCTGGGCTGAAGAGGTGAGCCCAGAGGAGCTGGTGGTCTTGAAGCAGCTCCGCC
+
CCCFFFFFHHHHHIJIIIJFGCHIIIJJJIJJIJIJJJJJJJJJJJJJJJGHHHFFDDFC=ACEDDDB=?BD@ABDBDDCCBDCCCCCCCCCDDCDDDD?
@HS2000-1015_160:6:2116:4077:79041/2
GGTCCCCGCCTACGCCCACTGGGTTGGTGCACCTGGTGGTGGTGGCCGCCAAGAAGCTGGTGAACCGCCTCCAAGTGGCTCCCAAGACGCAGCTGGATGA
+
CCCFFFFFHHFHHJHJJJJJJJJGGHJHIGAGIIJIFHJ;@F;CHHFHFDDDDDCDDCDD9CCCDDBDDBBDDCDACDD8@BD3>?BCDBDDDACCDC@>
@HS2000-1015_160:5:1110:16577:55703/2
AGAAGCTGAAGTACGGTGTAGGAATGTTTCTACCACGCATTTTCCCTGCACCCTCTCACATCATTCATGGGCCATTTGTAGACAGATGGGAAAGAGATGT
+
CC@FFFFFHGDDHGIG:C<FBHEFHIGGGIIJIIIJIFJJJJIIHIJIHHIIDHIIIJJJJJJJIHHGEHHGEFFFCEDECDDEE>@CDDCCDBCDCCDC
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
@HS2000-1015_160:5:1216:14609:76938/2
TCAAAAGGCCTCTGGAGAGATTTCTGGGTGTGGGAGGTGCTATTTTGGGGCTTGGCCAGTCATATGGAGATAAGCCTACAAGGTTGGGGACCTGGCAGAT
+
CCCFFFFFFHHGHJJHHCCAFHIGHHJJ:E?FHHHDG?FGGIIJJIJIIJIJJIJHGCEECC@CBEFECEC@C@@CCDDDDDD>CCBDDD@??CBDC72:
@HS2000-1015_160:6:2209:18284:44195/1
GCAGCCAGCGGACGTCCAGGAACCGGGATGCCTCCAGCAGTGAGGCGGTCAGCCTGCAGCATGGGATGGCTGTGGATCTTTGGGGCAGCCCTGGGGCAGT
+
CCCFFFFFHHHHHIJHIIJJJJJJIJJEHIJJGIJJJJGHCHEHGFDDDDDCDDDDDDCDDCDDDCCBDCDDBCDDDDDDDDDBBDBDDDDDBDDDDBDC
@HS2000-1015_160:6:2209:18284:44195/2
TAAAATGTCACAAAGCTGGAAACTCTTCCCTATCACAAACCAAAACTTAAAAGGACGTTACCTGGCTGGGTCTAAACTCCACATAACTCGCTTGCAGTTG
+
CCCFFFFEHHHGHJIIIJJIJJHIIJEHJJHIJJJIIJJIJIJJIJIIHJJIJGGHGHGIIHHIIIIHFH@DFFFDEEEECDDDCDDDDBDDBBDCDACC
我想要的结果:
@HS2000-1015_160:7:1108:13370:100570/1
GTGGGTCTTTGATCTTACAGGACAAATCCTTCAGAAAGATCAGGAAATCAAGGACCTAAAACAGAAGATAGCCGAAGTCATGGCCGTCATGCCCAGCATA
+
CCCFFFFFHHHFHGGIJIJDHEIIIGEHGGII>CFII<BG<DGCADDHBD>>BC?(8CGCCHIEC;@9@D47?BAD4,.(.6.;(3=8<A>@:<@AB9AA
@HS2000-1015_160:7:1108:13370:100570/2
CTTGACTGCCAGAGACGCTCCTTTGCAATGCCTTCCGGTAACCAAATTTTTGGGCACAACACACAGCTGGCCTTCATTTCTTCAGGGGCTGGTAAACAGA
+
@@@ADFFFHHHFD=EF@:GHIIFHH<ECHGF@DDBB:6@D60?F=888)8='--(=5@EAE5?'(..((.;?@>>A>3;@####################
+
@HS2000-1015_160:5:2306:10070:71746/1
AGAAGACGAGCACAGTGAGGAAGGAGATGAGCAGGCAGGGGATGATGAGGTTGATGGTGTAGAACAAGGGCAGGCGCCGGATGTACAGCGAGTATGTGAT
+
CCCFFFFFHHHHHIICGGHHHAFGEFGIJJJIIJIIGIIGG;FGHEHHGI=CHCEHGHE?DFFDFECEDDDDDDDDDBDDBD@ACDAACDBD5<CDC>@>
@HS2000-1015_160:5:2306:10070:71746/2
GAACCTCAAGGACTATTGGGAGAGCGGCGAGTGGGCCATCATCAAAGCCCCAGGCTACAAACACGACATCAAGTACAACTGCTGCGAGGAGATCTACCCC
+
@CCFFFFDHGHFHIJJJJJJGGGIIJJIGHI@FHGIIGHHEFGHHFFFFFBCDDDDDCDDDDDDDD;@BDCCDACDD@>ACCDDDDBDB<BA?C@CC@BD
@HS2000-1015_160:6:2116:4077:79041/1
TTGGGAGGGACATATCTAGGGTGGACCAAGCGCTGGCCAGGACTCACCTGGGCTGAAGAGGTGAGCCCAGAGGAGCTGGTGGTCTTGAAGCAGCTCCGCC
+
CCCFFFFFHHHHHIJIIIJFGCHIIIJJJIJJIJIJJJJJJJJJJJJJJJGHHHFFDDFC=ACEDDDB=?BD@ABDBDDCCBDCCCCCCCCCDDCDDDD?
@HS2000-1015_160:6:2116:4077:79041/2
GGTCCCCGCCTACGCCCACTGGGTTGGTGCACCTGGTGGTGGTGGCCGCCAAGAAGCTGGTGAACCGCCTCCAAGTGGCTCCCAAGACGCAGCTGGATGA
+
CCCFFFFFHHFHHJHJJJJJJJJGGHJHIGAGIIJIFHJ;@F;CHHFHFDDDDDCDDCDD9CCCDDBDDBBDDCDACDD8@BD3>?BCDBDDDACCDC@>
+
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
+
@HS2000-1015_160:6:2209:18284:44195/1
GCAGCCAGCGGACGTCCAGGAACCGGGATGCCTCCAGCAGTGAGGCGGTCAGCCTGCAGCATGGGATGGCTGTGGATCTTTGGGGCAGCCCTGGGGCAGT
+
CCCFFFFFHHHHHIJHIIJJJJJJIJJEHIJJGIJJJJGHCHEHGFDDDDDCDDDDDDCDDCDDDCCBDCDDBCDDDDDDDDDBBDBDDDDDBDDDDBDC
@HS2000-1015_160:6:2209:18284:44195/2
TAAAATGTCACAAAGCTGGAAACTCTTCCCTATCACAAACCAAAACTTAAAAGGACGTTACCTGGCTGGGTCTAAACTCCACATAACTCGCTTGCAGTTG
+
CCCFFFFEHHHGHJIIIJJIJJHIIJEHJJHIJJJIIJJIJIJJIJIIHJJIJGGHGHGIIHHIIIIHFH@DFFFDEEEECDDDCDDDDBDDBBDCDACC
Assumptions/assertions:
- 所有组都有4行(即没有3行或5行组)
- 文件的大小 (80 GB) 无法将整个文件加载到内存中(例如,
awk
数组)用于 grouping/sorting 目的
- 文件中不存在
~
字符
- 每个唯一的 header 最多匹配文件中的 2 行,并且所述匹配具有后缀
/1
和 /2
一个awk/sort/awk
想法:
- 使用
~
字符将每组 4 行附加在一起(有效地将 \n
替换为 ~
)
- 使用通用
sort
对附加行进行排序
- 使用
awk
仅去除那些匹配 /1
和 /2
条目的行
- 将
~
字符替换为换行符 (\n
)
注意:fastq.dat
包含问题中提供的 'before' 数据的副本
awk '{ (FNR % 4 == 0) ? sfx="\n" : sfx="~"; printf "%s%s", [=10=] , sfx }' fastq.dat |
sort |
awk -F '/' 'prevhdr != { prevhdr= ; prevline=[=10=] ; next } { gsub(/~/, "\n", prevline) ; print prevline ; gsub(/~/, "\n", [=10=]) ; print }' > fastq.final
注释:
- 显然需要足够的 space 用于另一个 80 GB 的文件(可能在不同的文件系统上,以便最大限度地减少 2 个进程 reading/writing 80 GB 在同一底层设备上的 FS 抖动?)
- 假设此时 OS 可以处理飞行中的内存需求,否则可能需要在写入各个输出文件时分别执行每个步骤 (???)
OP 评论了内存使用和生成输出速度慢的问题。可能有兴趣查看将进程分解为单独的命令是否有帮助,尽管以需要更多磁盘为代价 space,例如:
awk '{ (FNR % 4 == 0) ? sfx="\n" : sfx="~"; printf "%s%s", [=11=] , sfx }' fastq.dat > fastq.appended
sort fastq.appended > fastq.sorted
awk -F '/' 'prevhdr != { prevhdr= ; prevline=[=11=] ; next } { gsub(/~/, "\n", prevline) ; print prevline ; gsub(/~/, "\n", [=11=]) ; print }' fastq.sorted > fastq.final
'rm' -rf fastq.appended fastq.sorted 2>/dev/null
假设sort
操作可能是大猪(内存和慢),这个answer about speeding up sort可能是有趣的。大意是增加sort
使用的内存,允许并行线程,eg:
sort -S 50% --parallel=4
显示已维护4行群:
$ head -8 fastq.final
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
仅拉出 @HS
行以显示数据已按需要排序:
$ grep '@HS' fastq.final
@HS2000-1015_160:5:2113:11446:94436/1
@HS2000-1015_160:5:2113:11446:94436/2
@HS2000-1015_160:5:2306:10070:71746/1
@HS2000-1015_160:5:2306:10070:71746/2
@HS2000-1015_160:6:2116:4077:79041/1
@HS2000-1015_160:6:2116:4077:79041/2
@HS2000-1015_160:6:2209:18284:44195/1
@HS2000-1015_160:6:2209:18284:44195/2
@HS2000-1015_160:7:1108:13370:100570/1
@HS2000-1015_160:7:1108:13370:100570/2
注意:手动添加的空白行仅用于显示目的
评论有点长...想知道以下任何一个是否更适合处理这些大型数据集。
注意:我不使用 fastq
数据,以下是基于约 15 分钟的谷歌搜索,所以 fwiw ...
1 - 以允许更快处理的方式生成 fastq 文件
fastq-dump 的文档提到使用 --split-3
标志会将 fastq 数据转储到 3x 文件中......一个 'left' 文件,一个 'right' 文件和一个 'singleton' 文件。
'left' 和 'right' 文件会分别引用您的 /1
和 /2
后缀吗?如果 'yes' 那么我想知道您是否可以生成 3x 文件(然后丢弃 'singleton' 文件)?
也有兴趣查看 'left' 和 'right' 文件中的数据是否以相同的顺序排列 headers(在这种情况下,某种 'merge' 在 2x 文件之间加入可能比手动 pre-sorting 然后匹配更快)。
假设可以生成 3x 文件,这就引出了第二个想法...
2 - 是否有设计用于匹配 left/right (/1 & /2
) 数据集的工具?
fastq-pair的简单介绍似乎表明可以匹配多个文件(如fastq-dump
生成的'left'和'right'文件 - 上面??)比许多蛮力解决方案更快(例如,参考的 awk
解决方案与我的第一个答案有点相似)。
我有一个大约 80GB 大小的 fastq 文件 (file.fastq
),其中有一个 header 行和三个后续信息行。我需要在匹配 header 信息的 header 行中匹配 /1 和 /2,并将它们放在一个按 /1 和 /2 排序的文件中(交错)。
例如,我想为 @HS2000-1015_160:5:2114:14081:69256/1
找到 @HS2000-1015_160:5:2114:14081:69256/2
,如果找到,我想将它们堆叠在一起(对其他 header 也重复此操作)并得到如下所示的结果。
我在想这样的事情,但不知道该怎么做
cat file.fastq | grep -A3 --no-group-separator ...
file.fastq
@HS2000-1015_160:5:1114:14809:86636/2
GTTTGGAAAACAGAAGGGATTCGTCTGGAAAGGGCGGGAAACAAAGGACTGGAAGGGAGTGATGGAAAAAAACATGATGAAACTGATCTGGAAAGAGCCA
+
@B@FFFFFHGHHHIIJJJHIJJJJIIJJJIJJJIJJIHFDDEDDBDDDDDBDDDDDDDD2?CDDDDDCCDDDDDDDDCCD@CDDDCDDEDDCDDDDC9<A
@HS2000-1015_160:5:2114:14081:69256/1
CACCCTGATCCAAGGTGGACCTCATATGCTATCCACCAATTCTTTCATTAATCTTTATCTGAGTTCCTATTACATATTTCAAAGAAACTGTCAAACAACT
+
CCCFFFFFHHHHHJJCGIIHIIJJJIJJJIIJJJJIJJ;FHIJIJJIGIGDIIJJJJJIJJIIFGIJIEIJJIJJIIIJJJHHFGHFFFFDFDEDEEDDD
@HS2000-1015_160:7:1108:13370:100570/1
GTGGGTCTTTGATCTTACAGGACAAATCCTTCAGAAAGATCAGGAAATCAAGGACCTAAAACAGAAGATAGCCGAAGTCATGGCCGTCATGCCCAGCATA
+
CCCFFFFFHHHFHGGIJIJDHEIIIGEHGGII>CFII<BG<DGCADDHBD>>BC?(8CGCCHIEC;@9@D47?BAD4,.(.6.;(3=8<A>@:<@AB9AA
@HS2000-1015_160:7:1108:13370:100570/2
CTTGACTGCCAGAGACGCTCCTTTGCAATGCCTTCCGGTAACCAAATTTTTGGGCACAACACACAGCTGGCCTTCATTTCTTCAGGGGCTGGTAAACAGA
+
@@@ADFFFHHHFD=EF@:GHIIFHH<ECHGF@DDBB:6@D60?F=888)8='--(=5@EAE5?'(..((.;?@>>A>3;@####################
@HS2000-1015_160:6:1306:14563:31825/1
TTTTAAACTTATGAACATCCAACAGCAACTAATGGTAAGTCAGAAAGCTTTTATTTTTAAGTATATGCATGTGTATGTGTGTGTTTTGTATGTGTGTGTC
+
CCCFFFFFHHGHHJEDHJIJJJJIIJJJJJJGGIICFHIDGHH@HGIIIGGHFHHG>GGDHEBFIJJJHGI@FCDEEEHHGHHEHHFFFEFFFEEEEEDD
@HS2000-1015_160:7:2308:13476:86835/2
ATGGTCCTAAGTAACGTGAGGCCCAGAGGGAGGAGACCTTCAGAGGAGGCCAGGATGGGGGTGATGGGGAACTTTTCCTGAAGCTGGAGTAGGGAAAGAT
+
@@CFDFFFHFHDHHGIGJGIJIJJEG@FHGCFGBF=DFGFHGG>GGGGHIHCHFCEHHFFD'8<ACDDDD38CCDDDDCDDACDCDC?B:>@C?>A?A>>
@HS2000-1015_160:5:2306:10070:71746/1
AGAAGACGAGCACAGTGAGGAAGGAGATGAGCAGGCAGGGGATGATGAGGTTGATGGTGTAGAACAAGGGCAGGCGCCGGATGTACAGCGAGTATGTGAT
+
CCCFFFFFHHHHHIICGGHHHAFGEFGIJJJIIJIIGIIGG;FGHEHHGI=CHCEHGHE?DFFDFECEDDDDDDDDDBDDBD@ACDAACDBD5<CDC>@>
@HS2000-1015_160:5:2306:10070:71746/2
GAACCTCAAGGACTATTGGGAGAGCGGCGAGTGGGCCATCATCAAAGCCCCAGGCTACAAACACGACATCAAGTACAACTGCTGCGAGGAGATCTACCCC
+
@CCFFFFDHGHFHIJJJJJJGGGIIJJIGHI@FHGIIGHHEFGHHFFFFFBCDDDDDCDDDDDDDD;@BDCCDACDD@>ACCDDDDBDB<BA?C@CC@BD
@HS2000-1015_160:6:2116:4077:79041/1
TTGGGAGGGACATATCTAGGGTGGACCAAGCGCTGGCCAGGACTCACCTGGGCTGAAGAGGTGAGCCCAGAGGAGCTGGTGGTCTTGAAGCAGCTCCGCC
+
CCCFFFFFHHHHHIJIIIJFGCHIIIJJJIJJIJIJJJJJJJJJJJJJJJGHHHFFDDFC=ACEDDDB=?BD@ABDBDDCCBDCCCCCCCCCDDCDDDD?
@HS2000-1015_160:6:2116:4077:79041/2
GGTCCCCGCCTACGCCCACTGGGTTGGTGCACCTGGTGGTGGTGGCCGCCAAGAAGCTGGTGAACCGCCTCCAAGTGGCTCCCAAGACGCAGCTGGATGA
+
CCCFFFFFHHFHHJHJJJJJJJJGGHJHIGAGIIJIFHJ;@F;CHHFHFDDDDDCDDCDD9CCCDDBDDBBDDCDACDD8@BD3>?BCDBDDDACCDC@>
@HS2000-1015_160:5:1110:16577:55703/2
AGAAGCTGAAGTACGGTGTAGGAATGTTTCTACCACGCATTTTCCCTGCACCCTCTCACATCATTCATGGGCCATTTGTAGACAGATGGGAAAGAGATGT
+
CC@FFFFFHGDDHGIG:C<FBHEFHIGGGIIJIIIJIFJJJJIIHIJIHHIIDHIIIJJJJJJJIHHGEHHGEFFFCEDECDDEE>@CDDCCDBCDCCDC
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
@HS2000-1015_160:5:1216:14609:76938/2
TCAAAAGGCCTCTGGAGAGATTTCTGGGTGTGGGAGGTGCTATTTTGGGGCTTGGCCAGTCATATGGAGATAAGCCTACAAGGTTGGGGACCTGGCAGAT
+
CCCFFFFFFHHGHJJHHCCAFHIGHHJJ:E?FHHHDG?FGGIIJJIJIIJIJJIJHGCEECC@CBEFECEC@C@@CCDDDDDD>CCBDDD@??CBDC72:
@HS2000-1015_160:6:2209:18284:44195/1
GCAGCCAGCGGACGTCCAGGAACCGGGATGCCTCCAGCAGTGAGGCGGTCAGCCTGCAGCATGGGATGGCTGTGGATCTTTGGGGCAGCCCTGGGGCAGT
+
CCCFFFFFHHHHHIJHIIJJJJJJIJJEHIJJGIJJJJGHCHEHGFDDDDDCDDDDDDCDDCDDDCCBDCDDBCDDDDDDDDDBBDBDDDDDBDDDDBDC
@HS2000-1015_160:6:2209:18284:44195/2
TAAAATGTCACAAAGCTGGAAACTCTTCCCTATCACAAACCAAAACTTAAAAGGACGTTACCTGGCTGGGTCTAAACTCCACATAACTCGCTTGCAGTTG
+
CCCFFFFEHHHGHJIIIJJIJJHIIJEHJJHIJJJIIJJIJIJJIJIIHJJIJGGHGHGIIHHIIIIHFH@DFFFDEEEECDDDCDDDDBDDBBDCDACC
我想要的结果:
@HS2000-1015_160:7:1108:13370:100570/1
GTGGGTCTTTGATCTTACAGGACAAATCCTTCAGAAAGATCAGGAAATCAAGGACCTAAAACAGAAGATAGCCGAAGTCATGGCCGTCATGCCCAGCATA
+
CCCFFFFFHHHFHGGIJIJDHEIIIGEHGGII>CFII<BG<DGCADDHBD>>BC?(8CGCCHIEC;@9@D47?BAD4,.(.6.;(3=8<A>@:<@AB9AA
@HS2000-1015_160:7:1108:13370:100570/2
CTTGACTGCCAGAGACGCTCCTTTGCAATGCCTTCCGGTAACCAAATTTTTGGGCACAACACACAGCTGGCCTTCATTTCTTCAGGGGCTGGTAAACAGA
+
@@@ADFFFHHHFD=EF@:GHIIFHH<ECHGF@DDBB:6@D60?F=888)8='--(=5@EAE5?'(..((.;?@>>A>3;@####################
+
@HS2000-1015_160:5:2306:10070:71746/1
AGAAGACGAGCACAGTGAGGAAGGAGATGAGCAGGCAGGGGATGATGAGGTTGATGGTGTAGAACAAGGGCAGGCGCCGGATGTACAGCGAGTATGTGAT
+
CCCFFFFFHHHHHIICGGHHHAFGEFGIJJJIIJIIGIIGG;FGHEHHGI=CHCEHGHE?DFFDFECEDDDDDDDDDBDDBD@ACDAACDBD5<CDC>@>
@HS2000-1015_160:5:2306:10070:71746/2
GAACCTCAAGGACTATTGGGAGAGCGGCGAGTGGGCCATCATCAAAGCCCCAGGCTACAAACACGACATCAAGTACAACTGCTGCGAGGAGATCTACCCC
+
@CCFFFFDHGHFHIJJJJJJGGGIIJJIGHI@FHGIIGHHEFGHHFFFFFBCDDDDDCDDDDDDDD;@BDCCDACDD@>ACCDDDDBDB<BA?C@CC@BD
@HS2000-1015_160:6:2116:4077:79041/1
TTGGGAGGGACATATCTAGGGTGGACCAAGCGCTGGCCAGGACTCACCTGGGCTGAAGAGGTGAGCCCAGAGGAGCTGGTGGTCTTGAAGCAGCTCCGCC
+
CCCFFFFFHHHHHIJIIIJFGCHIIIJJJIJJIJIJJJJJJJJJJJJJJJGHHHFFDDFC=ACEDDDB=?BD@ABDBDDCCBDCCCCCCCCCDDCDDDD?
@HS2000-1015_160:6:2116:4077:79041/2
GGTCCCCGCCTACGCCCACTGGGTTGGTGCACCTGGTGGTGGTGGCCGCCAAGAAGCTGGTGAACCGCCTCCAAGTGGCTCCCAAGACGCAGCTGGATGA
+
CCCFFFFFHHFHHJHJJJJJJJJGGHJHIGAGIIJIFHJ;@F;CHHFHFDDDDDCDDCDD9CCCDDBDDBBDDCDACDD8@BD3>?BCDBDDDACCDC@>
+
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
+
@HS2000-1015_160:6:2209:18284:44195/1
GCAGCCAGCGGACGTCCAGGAACCGGGATGCCTCCAGCAGTGAGGCGGTCAGCCTGCAGCATGGGATGGCTGTGGATCTTTGGGGCAGCCCTGGGGCAGT
+
CCCFFFFFHHHHHIJHIIJJJJJJIJJEHIJJGIJJJJGHCHEHGFDDDDDCDDDDDDCDDCDDDCCBDCDDBCDDDDDDDDDBBDBDDDDDBDDDDBDC
@HS2000-1015_160:6:2209:18284:44195/2
TAAAATGTCACAAAGCTGGAAACTCTTCCCTATCACAAACCAAAACTTAAAAGGACGTTACCTGGCTGGGTCTAAACTCCACATAACTCGCTTGCAGTTG
+
CCCFFFFEHHHGHJIIIJJIJJHIIJEHJJHIJJJIIJJIJIJJIJIIHJJIJGGHGHGIIHHIIIIHFH@DFFFDEEEECDDDCDDDDBDDBBDCDACC
Assumptions/assertions:
- 所有组都有4行(即没有3行或5行组)
- 文件的大小 (80 GB) 无法将整个文件加载到内存中(例如,
awk
数组)用于 grouping/sorting 目的 - 文件中不存在
~
字符 - 每个唯一的 header 最多匹配文件中的 2 行,并且所述匹配具有后缀
/1
和/2
一个awk/sort/awk
想法:
- 使用
~
字符将每组 4 行附加在一起(有效地将\n
替换为~
) - 使用通用
sort
对附加行进行排序 - 使用
awk
仅去除那些匹配/1
和/2
条目的行 - 将
~
字符替换为换行符 (\n
)
注意:fastq.dat
包含问题中提供的 'before' 数据的副本
awk '{ (FNR % 4 == 0) ? sfx="\n" : sfx="~"; printf "%s%s", [=10=] , sfx }' fastq.dat |
sort |
awk -F '/' 'prevhdr != { prevhdr= ; prevline=[=10=] ; next } { gsub(/~/, "\n", prevline) ; print prevline ; gsub(/~/, "\n", [=10=]) ; print }' > fastq.final
注释:
- 显然需要足够的 space 用于另一个 80 GB 的文件(可能在不同的文件系统上,以便最大限度地减少 2 个进程 reading/writing 80 GB 在同一底层设备上的 FS 抖动?)
- 假设此时 OS 可以处理飞行中的内存需求,否则可能需要在写入各个输出文件时分别执行每个步骤 (???)
OP 评论了内存使用和生成输出速度慢的问题。可能有兴趣查看将进程分解为单独的命令是否有帮助,尽管以需要更多磁盘为代价 space,例如:
awk '{ (FNR % 4 == 0) ? sfx="\n" : sfx="~"; printf "%s%s", [=11=] , sfx }' fastq.dat > fastq.appended
sort fastq.appended > fastq.sorted
awk -F '/' 'prevhdr != { prevhdr= ; prevline=[=11=] ; next } { gsub(/~/, "\n", prevline) ; print prevline ; gsub(/~/, "\n", [=11=]) ; print }' fastq.sorted > fastq.final
'rm' -rf fastq.appended fastq.sorted 2>/dev/null
假设sort
操作可能是大猪(内存和慢),这个answer about speeding up sort可能是有趣的。大意是增加sort
使用的内存,允许并行线程,eg:
sort -S 50% --parallel=4
显示已维护4行群:
$ head -8 fastq.final
@HS2000-1015_160:5:2113:11446:94436/1
AATTGCTTGAACCCGGGAGGCAGAGGTTGCAGTGAGCCAAGATCGCACCACTGCAGTCCAGCCTGGGAAACCAGAGCAAAACTTGGTCTCAAAAAAAAAA
+
CCCFFFFFHHHHGHGFGFHJIBF1DC<@;)9D*099?F;8=(BFF5(;'.=CEHFBDC7.>;?C@BBB35>(<?<<A>AC?B>4+(+>A:(:@(2)5@B#
@HS2000-1015_160:5:2113:11446:94436/2
CGTCAGGGCCAACCCCGCCCCACCCTGACCCTACCTGGCACCCCTCACCTGTGGCCTGCCAGCACAGCCTCGCCCCTGCTGGCCAATGTGTCCCCCGTCA
+
?@@DA@DDFHH?DHI)<@@FHDBGGCHCBDH;DFA<)6.=7D;@CBCHD)).7@=>;?==AABC95<(5(5309@D########################
仅拉出 @HS
行以显示数据已按需要排序:
$ grep '@HS' fastq.final
@HS2000-1015_160:5:2113:11446:94436/1
@HS2000-1015_160:5:2113:11446:94436/2
@HS2000-1015_160:5:2306:10070:71746/1
@HS2000-1015_160:5:2306:10070:71746/2
@HS2000-1015_160:6:2116:4077:79041/1
@HS2000-1015_160:6:2116:4077:79041/2
@HS2000-1015_160:6:2209:18284:44195/1
@HS2000-1015_160:6:2209:18284:44195/2
@HS2000-1015_160:7:1108:13370:100570/1
@HS2000-1015_160:7:1108:13370:100570/2
注意:手动添加的空白行仅用于显示目的
评论有点长...想知道以下任何一个是否更适合处理这些大型数据集。
注意:我不使用 fastq
数据,以下是基于约 15 分钟的谷歌搜索,所以 fwiw ...
1 - 以允许更快处理的方式生成 fastq 文件
fastq-dump 的文档提到使用 --split-3
标志会将 fastq 数据转储到 3x 文件中......一个 'left' 文件,一个 'right' 文件和一个 'singleton' 文件。
'left' 和 'right' 文件会分别引用您的 /1
和 /2
后缀吗?如果 'yes' 那么我想知道您是否可以生成 3x 文件(然后丢弃 'singleton' 文件)?
也有兴趣查看 'left' 和 'right' 文件中的数据是否以相同的顺序排列 headers(在这种情况下,某种 'merge' 在 2x 文件之间加入可能比手动 pre-sorting 然后匹配更快)。
假设可以生成 3x 文件,这就引出了第二个想法...
2 - 是否有设计用于匹配 left/right (/1 & /2
) 数据集的工具?
fastq-pair的简单介绍似乎表明可以匹配多个文件(如fastq-dump
生成的'left'和'right'文件 - 上面??)比许多蛮力解决方案更快(例如,参考的 awk
解决方案与我的第一个答案有点相似)。