k-mer字符串的生成

Posted huanping

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了k-mer字符串的生成相关的知识,希望对你有一定的参考价值。

原文

简要介绍

在生信中,k-mer指生物序列中长度为k的子序列。$k$-mer包含着生物序列的两个基本特征:
1. 单体组分信息
2. 序列顺序信息

通过两个信息可以基本确定一个序列。在许多生物相关领域都有着广泛应用,如用于序列组装(构建De Bruijn图),生物序列特征提取(机器学习的输入特征),序列分析等领域。

同时在不属于生物相关的其他领域,如计算语言学中,它被称为n-gram。

一条长度为$L$的序列,有$L-k+1$个k-mer。比如序列: "ATGCA"有5个1-mer(A,T,G,C,A),4个2-mer(AT,TG,GC,CA),2个4-mer(ATGC,TGCA)。对于一个由n种单体组成的序列,理论上可能存在$n^k$种k-mer, 比如对于DNA序列来说可能存在$4^4$种不同的4-mer。对于蛋白质一级结构序列来说,可能存在$20^3$种不同3-mer。

计算k-mer,当前已有优化过的软件来计算k-mer。这里代码k-mer,仅供学习使用。

代码

从这里看到一份有点技巧性的的代码

from itertools import islice

seq = "ATGCGGCAA"
n = 3
def n_grams(a, n):
    z = (islice(a, i, None) for i in range(n))
    return zip(*z)

for i in n_grams(seq, n):
    print(i)
#('A', 'T', 'G')
#('T', 'G', 'C')
#('G', 'C', 'G')
#('C', 'G', 'G')
#('G', 'G', 'C')
#('G', 'C', 'A')
#('C', 'A', 'A')

下文分解代码~

islice(iterable, start, stop[, step])

islice创建一个生成器,返回从 iterable 里选中的元素。start,stop为下标索引值,start表示从什么位置开始,stop 表示到什么位置结束,step为步伐,默认为1。

from itertools import islice
seq = "ATGCGGCAA"
k1 = islice(seq, 0, None)  # 从索引0开始,返回一个生成器
for i in k1: 
    print(i, end=" ")  # 默认是换行符,改成空格,一行输出
#  A T G C G G C A A

k2 = islice(seq, 1, None)  # 从索引1开始
for i in k2: 
    print(i, end=" ")  
#  T G C G G C A A

k3 = islice(seq, 2, None)  # 从索引2开始,
for i in k3: 
    print(i, end=" ")  
#  G C G G C A A

z = (islice(seq, i, None) for i in range(n))

z 也是一个生成器,可生成 islice(seq, 0, None) ,islice(seq, 1, None), islice(seq, 2, None)

zip(*z)

z = (islice(seq, i, None) for i in range(n))
z1 = zip(*z) ## *拆包,会将z迭代器拆成三个元素,每个元素是一个生成
器
### (islice(seq, 0, None) ,islice(seq, 1, None), islice(seq, 2, None) = *z
"""
kmer1 = [i for i in z1]

z2 = zip(('A', 'T', 'G', 'C', 'G', 'G', 'C', 'A', 'A'),
         ('T', 'G', 'C', 'G', 'G', 'C', 'A', 'A'),
         ('G', 'C', 'G', 'G', 'C', 'A', 'A'), )
kmer2 = [i for i in z2]

kmer1 等于kmer2。

zip 的作用是: 同时迭代zip() 里的每个可迭代对象,
并将其组合成一个新元组, 一旦其中某个对象迭代完成就停止。
如:
A T G
T G C
...
C A A
"""

for i in z1: 
    print(i) 
#('A', 'T', 'G')
#('T', 'G', 'C')
#('G', 'C', 'G')
#('C', 'G', 'G')
#('G', 'G', 'C')
#('G', 'C', 'A')
#('C', 'A', 'A')
技巧
假如我要生成3-ker的字符串, 
开始生成这样的三个字符串,利用islice的特性
ATGCGGCAA  # 即从0开始
TGCGGCAA  ## 去掉第1个字符,从1开始
GCGGCAA   ## 去掉第2个字符, 从2开始
之后同时迭代形成新的元组,利用zip的特性
(A, T, G)
(T, G, C)
(G, C, G)
...
这样就可保持这样的取值形式(seq[i], seq[i+1], seq[i+2])

其它

然而。。。我以为弄那花里胡哨的的,会提升代码运行速度。实际并没有。还不如下面这种形式直接点。

seq = "ATGCGGCAA"
for i in range(len(seq)-3+1):
    print(seq[i:i+3])
#('A', 'T', 'G')
#('T', 'G', 'C')
#('G', 'C', 'G')
#('C', 'G', 'G')
#('G', 'G', 'C')
#('G', 'C', 'A')
#('C', 'A', 'A')

参考

Sliding windows (n -grams) using zip and iterators
k-mer

以上是关于k-mer字符串的生成的主要内容,如果未能解决你的问题,请参考以下文章

关于K-mer

在滑动窗口中查找 k-mers

reads k-mer scaffold 知乎

c_cpp 从fermi2指数获得k-mer计数

基因组survey——K-mer频谱

用k-mer分析进行基因组调查:(一)基本原理