我正在尝试压缩基因组序列。这些是字母“A”、“C”、“G”和“T”的字符串。在文本文件形式中,它们存储为字符。由于它们只有 4 个,因此在一个文件中可以表示为两位。
换句话说:ACTG -> 00 01 10 11,而不是 8 位字符
此数据将被写回一个文件,其中每个字节代表 4 个字符。在 bash 脚本或 C 程序中执行此操作的最有效方法是什么?
谢谢!
最佳答案
这是一个首先对序列进行最低有效位编码的过滤器:
#include <stdio.h>
int main(void) {
unsigned i = 0;
int c, d = 0;
while ((c = getchar()) != EOF) {
switch (c) {
case 'A': d |= 0 << (2 * (i & 3)); break;
case 'C': d |= 1 << (2 * (i & 3)); break;
case 'T': d |= 2 << (2 * (i & 3)); break;
case 'G': d |= 3 << (2 * (i & 3)); break;
default: continue; // ignore all other characters
}
if ((++i & 3) == 0) {
putchar(d);
d = 0;
}
}
if (i & 3) {
putchar(d);
}
return 0;
}
这里是最高有效位在前(又名像素顺序)
#include <stdio.h>
int main(void) {
unsigned i = 0;
int c, d = 0;
while ((c = getchar()) != EOF) {
switch (c) {
case 'A': d = (d << 2) | 0; break;
case 'C': d = (d << 2) | 1; break;
case 'T': d = (d << 2) | 2; break;
case 'G': d = (d << 2) | 3; break;
default: continue; // ignore all other characters
}
if ((++i & 3) == 0) {
putchar(d);
d = 0;
}
}
if (i & 3) {
putchar(d << (2 * (3 - (i & 3))));
}
return 0;
}
注意事项:
序列被隐式填充
A
,最多为 4 碱基的倍数。将二进制数据写入
stdout
可能会在stdout
默认处于文本模式且语义与二进制模式不同的系统上产生不正确的输出(例如 Windows,不同于 OS/X 或 Unix)。
关于有效地将每个字 rune 件的 4 个字符序列压缩为 2 位,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/40559042/