一种具有歧义的热编码(核苷酸序列)
One hot encoding with ambiguitiy (nucleotide sequences)
核苷酸序列(或 DNA 序列)通常由 4 个碱基组成:ATGC,这是一种非常好的、简单且有效的编码方式,用于机器学习目的。
sequence = AAATGCC
ohe_sequence = [[1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1], [0, 0, 0, 1]]
但是,当您考虑到 RNA 序列和测序仪中有时会出现的错误时,会添加字母 UYRWSKMDVHBXN...当对其进行单热编码时,您最终会得到一个包含 17 行的矩阵,其中最后一行13行一般都是0。
这是非常低效的,并且没有赋予这些额外的(不明确的)字母所具有的生物学意义。
例如:
- T 和 U 可以互换
- Y代表有C或T
- N和X表示存在4个碱基中的任意一个(ATGC)
所以我制作了一本代表这种生物学意义的字典
nucleotide_dict = {'A': 'A', 'T':'T', 'U':'T', 'G':'G', 'C':'C', 'Y':['C', 'T'], 'R':['A', 'G'], 'W':['A', 'T'], 'S':['G', 'C'], 'K':['T', 'G'], 'M':['C', 'A'], 'D':['A', 'T', 'G'], 'V':['A', 'G', 'C'], 'H':['A', 'T', 'C'], 'B':['T', 'G', 'C'], 'X':['A', 'T', 'G', 'C'], 'N':['A', 'T', 'G', 'C']}
但我似乎无法弄清楚如何制作一个有效的单热编码脚本(或者是否有一种方法可以使用 scikit 学习模块来做到这一点),它利用这本字典来获得结果是这样的:
sequence = ANTUYCC
ohe_sequence = [[1, 0, 0, 0], [1, 1, 1, 1], [0, 1, 0, 0], [0, 1, 0, 0], [0, 1, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1]]
# or even better:
ohe_sequence = [[1, 0, 0, 0], [0.25, 0.25, 0.25, 0.25], [0, 1, 0, 0], [0, 1, 0, 0], [0, 0.5, 0, 0.5], [0, 0, 0, 1], [0, 0, 0, 1]]
我喜欢后一种方法,因为它更接近真实的含义: Y
并不表示 C
和 T
,而是 C
或 T
。如果没有进一步的信息可用,假设概率相等(即权重)似乎是合理的。当然,这种 non-standard 编码需要通过损失函数的选择来反映。
回答您的问题:您可以预先计算从字母到编码的映射,然后创建一个 encode
函数,该函数将 sequence
作为输入并 returns 编码序列作为 (len(sequence), 4)
形 np.array
如下:
import numpy as np
nucleotide_dict = {'A':'A', 'T':'T', 'U':'T', 'G':'G', 'C':'C', 'Y':['C', 'T'], 'R':['A', 'G'], 'W':['A', 'T'], 'S':['G', 'C'], 'K':['T', 'G'], 'M':['C', 'A'], 'D':['A', 'T', 'G'], 'V':['A', 'G', 'C'], 'H':['A', 'T', 'C'], 'B':['T', 'G', 'C'], 'X':['A', 'T', 'G', 'C'], 'N':['A', 'T', 'G', 'C']}
index_mapper = {'A': 0, 'T': 1, 'G': 2, 'C': 3}
mapper_dict = dict()
for k, v in nucleotide_dict.items():
encoding = np.zeros(4)
p = 1 / len(v)
encoding[[index_mapper[i] for i in v]] = p
mapper_dict[k] = encoding
def encode(sequence):
return np.array([mapper_dict[s] for s in sequence])
这似乎产生了预期的结果,并且可能有点高效。
一个例子:
print(encode('AYSDX'))
打印以下内容:
array([[1. , 0. , 0. , 0. ],
[0. , 0.5 , 0. , 0.5 ],
[0. , 0. , 0.5 , 0.5 ],
[0.33333333, 0.33333333, 0.33333333, 0. ],
[0.25 , 0.25 , 0.25 , 0.25 ]])
这很有趣!我认为您可以使用具有适当值的字典来做到这一点。我添加了 scikit-learn class 因为你提到你正在使用它。请参阅下面的 transform
:
from sklearn.base import BaseEstimator, TransformerMixin
import numpy as np
nucleotide_dict = {
"A": [1, 0, 0, 0],
"G": [0, 1, 0, 0],
"C": [0, 0, 1, 0],
"T": [0, 0, 0, 1],
"U": [0, 0, 0, 1],
"Y": [0, 0, 1, 1],
"R": [1, 1, 0, 0],
"W": [1, 0, 0, 1],
"S": [0, 1, 1, 0],
"K": [0, 1, 0, 1],
"M": [1, 0, 1, 0],
"D": [1, 1, 0, 1],
"V": [1, 1, 1, 0],
"H": [1, 0, 1, 1],
"B": [0, 1, 1, 1],
"X": [1, 1, 1, 1],
"N": [1, 1, 1, 1],
"-": [0, 0, 0, 0],
}
class NucleotideEncoder(BaseEstimator, TransformerMixin):
def __init__(self, norm=True):
self.norm = norm
def fit(self, X, y=None):
return self
def transform(self, X, y=None):
f1 = lambda a: list(a)
f2 = lambda g: nucleotide_dict[g]
f3 = lambda c: list(map(f2, f1(c[0])))
f4 = lambda t: np.array(f3(t)) / np.sum(np.array(f3(t)), axis=1)[:, np.newaxis]
f = f3
if self.norm:
f = f4
return np.apply_along_axis(f, 1, X)
samples = np.array([["AAATGCC"], ["ANTUYCC"]])
print(NucleotideEncoder().fit_transform(samples))
输出:
[[[1. 0. 0. 0. ]
[1. 0. 0. 0. ]
[1. 0. 0. 0. ]
[0. 0. 0. 1. ]
[0. 1. 0. 0. ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]
[[1. 0. 0. 0. ]
[0.25 0.25 0.25 0.25]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]]
另一个允许不同长度序列的版本:
import numpy as np
nucleotide_dict = {
"A": [1, 0, 0, 0],
"G": [0, 1, 0, 0],
"C": [0, 0, 1, 0],
"T": [0, 0, 0, 1],
"U": [0, 0, 0, 1],
"Y": [0, 0, 1, 1],
"R": [1, 1, 0, 0],
"W": [1, 0, 0, 1],
"S": [0, 1, 1, 0],
"K": [0, 1, 0, 1],
"M": [1, 0, 1, 0],
"D": [1, 1, 0, 1],
"V": [1, 1, 1, 0],
"H": [1, 0, 1, 1],
"B": [0, 1, 1, 1],
"X": [1, 1, 1, 1],
"N": [1, 1, 1, 1],
"-": [0, 0, 0, 0],
}
norm = True
samples = np.array(["AAATGCC", "ANTUYCC", "".join(list(nucleotide_dict.keys()))[:-1]])
def nucleotide_encode(samples, norm=True):
m = map(list, samples)
f1 = lambda x: np.array(list(map(nucleotide_dict.get, x)))
f = f1
if norm:
f = lambda x: np.nan_to_num(f1(x) / np.sum(f1(x), axis=1)[:, np.newaxis])
return list(map(f, m))
for i, j in zip(samples, nucleotide_encode(samples, norm=norm)):
print(i)
print(j)
产量:
AAATGCC
[[1. 0. 0. 0.]
[1. 0. 0. 0.]
[1. 0. 0. 0.]
[0. 0. 0. 1.]
[0. 1. 0. 0.]
[0. 0. 1. 0.]
[0. 0. 1. 0.]]
ANTUYCC
[[1. 0. 0. 0. ]
[0.25 0.25 0.25 0.25]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]
AGCTUYRWSKMDVHBXN
[[1. 0. 0. 0. ]
[0. 1. 0. 0. ]
[0. 0. 1. 0. ]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0.5 0.5 0. 0. ]
[0.5 0. 0. 0.5 ]
[0. 0.5 0.5 0. ]
[0. 0.5 0. 0.5 ]
[0.5 0. 0.5 0. ]
[0.33333333 0.33333333 0. 0.33333333]
[0.33333333 0.33333333 0.33333333 0. ]
[0.33333333 0. 0.33333333 0.33333333]
[0. 0.33333333 0.33333333 0.33333333]
[0.25 0.25 0.25 0.25 ]
[0.25 0.25 0.25 0.25 ]]
核苷酸序列(或 DNA 序列)通常由 4 个碱基组成:ATGC,这是一种非常好的、简单且有效的编码方式,用于机器学习目的。
sequence = AAATGCC
ohe_sequence = [[1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1], [0, 0, 0, 1]]
但是,当您考虑到 RNA 序列和测序仪中有时会出现的错误时,会添加字母 UYRWSKMDVHBXN...当对其进行单热编码时,您最终会得到一个包含 17 行的矩阵,其中最后一行13行一般都是0。
这是非常低效的,并且没有赋予这些额外的(不明确的)字母所具有的生物学意义。
例如:
- T 和 U 可以互换
- Y代表有C或T
- N和X表示存在4个碱基中的任意一个(ATGC)
所以我制作了一本代表这种生物学意义的字典
nucleotide_dict = {'A': 'A', 'T':'T', 'U':'T', 'G':'G', 'C':'C', 'Y':['C', 'T'], 'R':['A', 'G'], 'W':['A', 'T'], 'S':['G', 'C'], 'K':['T', 'G'], 'M':['C', 'A'], 'D':['A', 'T', 'G'], 'V':['A', 'G', 'C'], 'H':['A', 'T', 'C'], 'B':['T', 'G', 'C'], 'X':['A', 'T', 'G', 'C'], 'N':['A', 'T', 'G', 'C']}
但我似乎无法弄清楚如何制作一个有效的单热编码脚本(或者是否有一种方法可以使用 scikit 学习模块来做到这一点),它利用这本字典来获得结果是这样的:
sequence = ANTUYCC
ohe_sequence = [[1, 0, 0, 0], [1, 1, 1, 1], [0, 1, 0, 0], [0, 1, 0, 0], [0, 1, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1]]
# or even better:
ohe_sequence = [[1, 0, 0, 0], [0.25, 0.25, 0.25, 0.25], [0, 1, 0, 0], [0, 1, 0, 0], [0, 0.5, 0, 0.5], [0, 0, 0, 1], [0, 0, 0, 1]]
我喜欢后一种方法,因为它更接近真实的含义: Y
并不表示 C
和 T
,而是 C
或 T
。如果没有进一步的信息可用,假设概率相等(即权重)似乎是合理的。当然,这种 non-standard 编码需要通过损失函数的选择来反映。
回答您的问题:您可以预先计算从字母到编码的映射,然后创建一个 encode
函数,该函数将 sequence
作为输入并 returns 编码序列作为 (len(sequence), 4)
形 np.array
如下:
import numpy as np
nucleotide_dict = {'A':'A', 'T':'T', 'U':'T', 'G':'G', 'C':'C', 'Y':['C', 'T'], 'R':['A', 'G'], 'W':['A', 'T'], 'S':['G', 'C'], 'K':['T', 'G'], 'M':['C', 'A'], 'D':['A', 'T', 'G'], 'V':['A', 'G', 'C'], 'H':['A', 'T', 'C'], 'B':['T', 'G', 'C'], 'X':['A', 'T', 'G', 'C'], 'N':['A', 'T', 'G', 'C']}
index_mapper = {'A': 0, 'T': 1, 'G': 2, 'C': 3}
mapper_dict = dict()
for k, v in nucleotide_dict.items():
encoding = np.zeros(4)
p = 1 / len(v)
encoding[[index_mapper[i] for i in v]] = p
mapper_dict[k] = encoding
def encode(sequence):
return np.array([mapper_dict[s] for s in sequence])
这似乎产生了预期的结果,并且可能有点高效。
一个例子:
print(encode('AYSDX'))
打印以下内容:
array([[1. , 0. , 0. , 0. ],
[0. , 0.5 , 0. , 0.5 ],
[0. , 0. , 0.5 , 0.5 ],
[0.33333333, 0.33333333, 0.33333333, 0. ],
[0.25 , 0.25 , 0.25 , 0.25 ]])
这很有趣!我认为您可以使用具有适当值的字典来做到这一点。我添加了 scikit-learn class 因为你提到你正在使用它。请参阅下面的 transform
:
from sklearn.base import BaseEstimator, TransformerMixin
import numpy as np
nucleotide_dict = {
"A": [1, 0, 0, 0],
"G": [0, 1, 0, 0],
"C": [0, 0, 1, 0],
"T": [0, 0, 0, 1],
"U": [0, 0, 0, 1],
"Y": [0, 0, 1, 1],
"R": [1, 1, 0, 0],
"W": [1, 0, 0, 1],
"S": [0, 1, 1, 0],
"K": [0, 1, 0, 1],
"M": [1, 0, 1, 0],
"D": [1, 1, 0, 1],
"V": [1, 1, 1, 0],
"H": [1, 0, 1, 1],
"B": [0, 1, 1, 1],
"X": [1, 1, 1, 1],
"N": [1, 1, 1, 1],
"-": [0, 0, 0, 0],
}
class NucleotideEncoder(BaseEstimator, TransformerMixin):
def __init__(self, norm=True):
self.norm = norm
def fit(self, X, y=None):
return self
def transform(self, X, y=None):
f1 = lambda a: list(a)
f2 = lambda g: nucleotide_dict[g]
f3 = lambda c: list(map(f2, f1(c[0])))
f4 = lambda t: np.array(f3(t)) / np.sum(np.array(f3(t)), axis=1)[:, np.newaxis]
f = f3
if self.norm:
f = f4
return np.apply_along_axis(f, 1, X)
samples = np.array([["AAATGCC"], ["ANTUYCC"]])
print(NucleotideEncoder().fit_transform(samples))
输出:
[[[1. 0. 0. 0. ]
[1. 0. 0. 0. ]
[1. 0. 0. 0. ]
[0. 0. 0. 1. ]
[0. 1. 0. 0. ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]
[[1. 0. 0. 0. ]
[0.25 0.25 0.25 0.25]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]]
另一个允许不同长度序列的版本:
import numpy as np
nucleotide_dict = {
"A": [1, 0, 0, 0],
"G": [0, 1, 0, 0],
"C": [0, 0, 1, 0],
"T": [0, 0, 0, 1],
"U": [0, 0, 0, 1],
"Y": [0, 0, 1, 1],
"R": [1, 1, 0, 0],
"W": [1, 0, 0, 1],
"S": [0, 1, 1, 0],
"K": [0, 1, 0, 1],
"M": [1, 0, 1, 0],
"D": [1, 1, 0, 1],
"V": [1, 1, 1, 0],
"H": [1, 0, 1, 1],
"B": [0, 1, 1, 1],
"X": [1, 1, 1, 1],
"N": [1, 1, 1, 1],
"-": [0, 0, 0, 0],
}
norm = True
samples = np.array(["AAATGCC", "ANTUYCC", "".join(list(nucleotide_dict.keys()))[:-1]])
def nucleotide_encode(samples, norm=True):
m = map(list, samples)
f1 = lambda x: np.array(list(map(nucleotide_dict.get, x)))
f = f1
if norm:
f = lambda x: np.nan_to_num(f1(x) / np.sum(f1(x), axis=1)[:, np.newaxis])
return list(map(f, m))
for i, j in zip(samples, nucleotide_encode(samples, norm=norm)):
print(i)
print(j)
产量:
AAATGCC
[[1. 0. 0. 0.]
[1. 0. 0. 0.]
[1. 0. 0. 0.]
[0. 0. 0. 1.]
[0. 1. 0. 0.]
[0. 0. 1. 0.]
[0. 0. 1. 0.]]
ANTUYCC
[[1. 0. 0. 0. ]
[0.25 0.25 0.25 0.25]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0. 0. 1. 0. ]
[0. 0. 1. 0. ]]
AGCTUYRWSKMDVHBXN
[[1. 0. 0. 0. ]
[0. 1. 0. 0. ]
[0. 0. 1. 0. ]
[0. 0. 0. 1. ]
[0. 0. 0. 1. ]
[0. 0. 0.5 0.5 ]
[0.5 0.5 0. 0. ]
[0.5 0. 0. 0.5 ]
[0. 0.5 0.5 0. ]
[0. 0.5 0. 0.5 ]
[0.5 0. 0.5 0. ]
[0.33333333 0.33333333 0. 0.33333333]
[0.33333333 0.33333333 0.33333333 0. ]
[0.33333333 0. 0.33333333 0.33333333]
[0. 0.33333333 0.33333333 0.33333333]
[0.25 0.25 0.25 0.25 ]
[0.25 0.25 0.25 0.25 ]]