K核苷酸频率(KNF,k-nucleotide frequencies)或K-mer频率
KNF描述了序列中存在k个核苷酸的所有可能的多核苷酸的频率。如果k=2,则计算的为双核苷酸频率(即AA、AT、AG、AC、……TT),共42=16种;如果k=3,则计算的为双核苷酸频率(即AAA、AAT、AAG、AAC、……TTT),共43=64种;以此类推。
K-mer频率亦如此。
方法一:
#提取核苷酸类型(排列组合)
from itertools import product
def nucleotide_type(k):
z = []
for i in product('ACGT', repeat = k): #笛卡尔积(有放回抽样排列)
z.append(''.join(i)) #把('A,A,A')转变成(AAA)形式
return z
# 碱基对数量统计
def char_count(sequence,num,k):
n = 0
char = nucleotide_type(k) #调用提取核苷酸类型模块
for i in range(len(sequence)): #统计相应字符出现的数量
if sequence[i:i+k] == char[num]:
n += 1
return n/(len(sequence)-k+1) #返回频率(出现的次数/总次数)总次数=序列长度-取几个碱基+1
def feature_336d(seq,k):
list = []
for i in range(4**k): #根据核苷酸类型数量来取值(二、三、四核苷酸分别循环16、64、256次)
list.append(char_count(seq,i,k))
return (list)
# 逐行调用特征编码
def Sequence_replacement(sequ,k):
sequen = [None]*len(sequ)
for i in range(len(sequ)):
s = sequ[i]
sequen[i] = feature_336d(s,k)
return sequen
#用具体数据进行调用
feature_knf = Sequence_replacement(data,k) #data为具体数据,k的值根据需要自己设定
方法二:
#首先把数据划分成K-mer形式
def Kmers_funct(seq,x):
X = [None]*len(seq) #若数据只有一个序列,可不用此定义
for i in range(len(seq)): #若数据只有一个序列,可不用此循环
a = seq[i]
t=0
l=[]
for index in range(len(a)):
t=a[index:index+x]
if (len(t))==x:
l.append(t)
X[i] = l
return X #具体看返回需要,也可直接:return X
#提取核苷酸类型(排列组合)
from itertools import product
def nucleotide_type(k):
z = []
for i in product('ACGU', repeat = k): #笛卡尔积(有放回抽样排列)
z.append(''.join(i)) #把('A,A,A')转变成(AAA)形式
return z
#定义K-mer频率模块
def Kmers_frequency(seq,x):
X = []
char = nucleotide_type(x) #调用提取核苷酸类型(排列组合)代码
for i in range(len(seq)):
s = seq[i]
frequence = []
for a in char:
number = s.count(a) #依次统计字符数量
char_frequence = number/len(s) #计算频率
frequence.append(char_frequence)
X.append(frequence)
return X
#调用K-mer模块代码
feature_kmer = Kmers_frequency(data,k)
#用k-mer生成的数据调用K-mer频率模块
feature_kmer_frequency = Kmers_frequency(feature_kmer,k)