C + + 中的简单字典

将一些代码从 Python 迁移到 C + + 。

BASEPAIRS = { "T": "A", "A": "T", "G": "C", "C": "G" }

认为地图可能是过度杀伤? 你会用什么?

281827 次浏览

字符数组中的表:

char map[256] = { 0 };
map['T'] = 'A';
map['A'] = 'T';
map['C'] = 'G';
map['G'] = 'C';
/* .... */

在我真正关心性能之前,我会使用一个函数,它接受一个基并返回它的匹配项:

char base_pair(char base)
{
switch(base) {
case 'T': return 'A';
... etc
default: // handle error
}
}

如果考虑到性能,我会将基值定义为四分之一字节。0代表 A,1代表 G,2代表 C,3代表 T。然后我把4个碱基打包成一个字节,为了得到它们的对,我只需要取补数。

这是地图的解决方案:

#include <iostream>
#include <map>


typedef std::map<char, char> BasePairMap;


int main()
{
BasePairMap m;
m['A'] = 'T';
m['T'] = 'A';
m['C'] = 'G';
m['G'] = 'C';


std::cout << "A:" << m['A'] << std::endl;
std::cout << "T:" << m['T'] << std::endl;
std::cout << "C:" << m['C'] << std::endl;
std::cout << "G:" << m['G'] << std::endl;


return 0;
}

BASEPAIRS = {“ T”: “ A”,“ A”: “ T”,“ G”: “ C”,“ C”: “ G”} 你会用什么?

也许:

static const char basepairs[] = "ATAGCG";
// lookup:
if (const char* p = strchr(basepairs, c))
// use p[1]

;-)

虽然使用 std::map是可以的,或者使用一个256大小的字符表也是可以的,但是通过简单地使用 enum,您可以为自己节省大量的空间。如果你有 C + + 11的特性,你可以使用 enum class强类型:

// First, we define base-pairs. Because regular enums
// Pollute the global namespace, I'm using "enum class".
enum class BasePair {
A,
T,
C,
G
};


// Let's cut out the nonsense and make this easy:
// A is 0, T is 1, C is 2, G is 3.
// These are indices into our table
// Now, everything can be so much easier
BasePair Complimentary[4] = {
T, // Compliment of A
A, // Compliment of T
G, // Compliment of C
C, // Compliment of G
};

用法变得简单:

int main (int argc, char* argv[] ) {
BasePair bp = BasePair::A;
BasePair complimentbp = Complimentary[(int)bp];
}

如果这对你来说太难了,你可以定义一些助手来获得人类可读的 ASCII 字符,也可以获得基对的赞美,这样你就不用一直进行 (int)转换了:

BasePair Compliment ( BasePair bp ) {
return Complimentary[(int)bp]; // Move the pain here
}


// Define a conversion table somewhere in your program
char BasePairToChar[4] = { 'A', 'T', 'C', 'G' };
char ToCharacter ( BasePair bp ) {
return BasePairToChar[ (int)bp ];
}

干净,简单,高效。

现在,突然之间,您就没有一个256字节的表了。您也没有存储字符(每个字节1字节) ,因此,如果您将其写入一个文件,您可以每个基对写入2位,而不是每个基对写入1字节(8位)。我必须使用生物信息学文件,每个文件存储一个字符的数据。好处是它是人类可读的。缺点是本来应该是250MB 的文件最终占用了1GB 的空间。移动,储存和使用是一场噩梦。当然,250MB 是 慷慨,甚至占蠕虫 DNA。无论如何,没有人会读取1GB 的碱基对。

您可以使用以下语法:

#include <map>


std::map<char, char> my_map = {
{ 'A', '1' },
{ 'B', '2' },
{ 'C', '3' }
};

如果你想进行优化,并假设输入总是四个字符中的一个,下面的函数可能值得一试,作为地图的替代品:

char map(const char in)
{ return ((in & 2) ? '\x8a' - in : '\x95' - in); }

它的工作原理是基于您正在处理两个对称对的事实。条件句的作用是区分 A/T 对和 G/C 对(‘ G’和‘ C’碰巧有第二个最低有效位的共同点)。其余的算法执行对称映射。它是基于 a = (a + b)-b 对于任何 a,b 都是真的这一事实。

这是我能想到的最快,最简单,最小的空间解决方案。一个好的编译器最佳化甚至可以降低访问对和名称数组的成本。这个解决方案在 C 语言中同样适用。

#include <iostream>


enum Base_enum { A, C, T, G };
typedef enum Base_enum Base;
static const Base pair[4] = { T, G, A, C };
static const char name[4] = { 'A', 'C', 'T', 'G' };
static const Base base[85] =
{ -1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1, -1,  A, -1,  C, -1, -1,
-1,  G, -1, -1, -1, -1, -1, -1, -1, -1,
-1, -1, -1, -1,  T };


const Base
base2 (const char b)
{
switch (b)
{
case 'A': return A;
case 'C': return C;
case 'T': return T;
case 'G': return G;
default: abort ();
}
}


int
main (int argc, char *args)
{
for (Base b = A; b <= G; b++)
{
std::cout << name[b] << ":"
<< name[pair[b]] << std::endl;
}
for (Base b = A; b <= G; b++)
{
std::cout << name[base[name[b]]] << ":"
<< name[pair[base[name[b]]]] << std::endl;
}
for (Base b = A; b <= G; b++)
{
std::cout << name[base2(name[b])] << ":"
<< name[pair[base2(name[b])]] << std::endl;
}
};

Base []是一个到 Base 的快速 ascii 字符(即 int 在0到3之间)查找,有点难看。一个好的编译器最佳化应该能够处理 base2() ,但我不确定是否有人能够做到这一点。