python基因
Ⅰ python-gffutils 修改gff3文件基因名称
基因结构注释文件一般为gff3的格式,一共是9列,依次为基因组序列id,注释来源,类型,起始位置,终止位置,得分,正负链,相位,属性。
基因结构注释文件中,基因包含mRNA,mRNA包含exon, CDS, UTR等信息,同时在注释文件中除基因行外,其他行在第9列会通过Parent指明该行从属的上一级ID,也就是一个基因的gene行、mRNA行、CDS行、exon行都会通过Parent层层关联在一起。
EVM软件整合出的基因结构注释文件一般如下所示,目前有一个需求,就是修改基因名称,修改为类似AT01G000001这种形式
gffutils是一个用于解析gff3/gtf格式的python包,读取gff3文件非常方便,利用这个包实现对基因结构注释文件重新命名python脚本如下
帮助文档,用法
我输入的change.bed文件内容如下
最终得到更换名称后的注释文件
Ⅱ python小程序--实现DNA到RNA的转换
在处理基因组信息的时候,我们时不时需要对DNA序列进行处理,比如将DNA转换成RNA,对于较少批量的数据来说,可以人工处理,但是在面对大的fastq文件的时候,会显得苍白无力,下面分别介绍两种处理方式
1.单个序列
2.文件处理
以上小程序可以进行转换。
Ⅲ 如何用PYTHON提取基因序列里的小写字符串(内含子)
s="序列"
printre.findall("[a-z]+",序列)
Ⅳ Python 基于基因表达量绘制热图
:
Ⅳ python A, C, T, 和 G构成的字符串来构建一个基因组,,在 TAG, TAA, 或者 TGA之前结束
#-*-coding:utf-8-*-
sstr=raw_input('PleaseinputtheGeneString: ')
endsplit=['TAG','TAA','TGA']
if'ATG'insstr:
foriinsstr.split('ATG'):
forjinendsplit:
ifjini:
printi.split(j)[0]
else:
print'error'
>>>
Please input the Gene String:
TTATGTTTTAAGGATGGGGCGTTAGTT
TTT
GGGCGT
Please input the Gene String:
TGTGTGTATAT
error
Ⅵ 利用gff提取某个基因的最长转录本(Python实现)
from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
import gffutils
import re
#######################
#从cds获得每条转录本的长度,从gff获得每个基因包含哪些转录本
fout = 'output/longest.fa'
cds = 'data/rna.fna'
gff = 'data/genomic.gff'
fout = 'output/longest.fa'
seq_index = SeqIO.index(cds ,'fasta')
gffdb = gffutils.create_db(gff,'gff.db'.force=True,
merge_strategy='creat_unique')
#以基因为key,转录本为value的字典
gene2longestcds = {}
#对基因进行迭代
for g in gffdb.all_features(featuretype='gene'):
g_id= g.id
for m in gffdb.children(g, featuretype='mRNA'):##只关注mRNA信息
m_id = m.id.replace('rna-' , '')
m_len = len(seq_index[m_id].seq)
#一个基因中最长转录本的名字
if g_id not in gene2longestcds:
gene2longestcds[g_id] = m_id
elif m_len >len(seq_index[gene2longestcds[g_id]].seq):
gene2longestcds[g_id] = m_id
sr_list = [v for k,v in seq_index.items() if k in gene2longestcds.values()]
##输出
SeqIO.write(sr_list, fout, 'fasta')
Ⅶ python有没有简单的遗传算法库
首先遗传算法是一种优化算法,通过模拟基因的优胜劣汰,进行计算(具体的算法思路什么的就不赘述了)。大致过程分为初始化编码、个体评价、选择,交叉,变异。
以目标式子 y = 10 * sin(5x) + 7 * cos(4x)为例,计算其最大值
首先是初始化,包括具体要计算的式子、种群数量、染色体长度、交配概率、变异概率等。并且要对基因序列进行初始化
[python]view plain
pop_size=500#种群数量
max_value=10#基因中允许出现的最大值
chrom_length=10#染色体长度
pc=0.6#交配概率
pm=0.01#变异概率
results=[[]]#存储每一代的最优解,N个二元组
fit_value=[]#个体适应度
fit_mean=[]#平均适应度
pop=geneEncoding(pop_size,chrom_length)
defgeneEncoding(pop_size,chrom_length):
pop=[[]]
foriinrange(pop_size):
temp=[]
forjinrange(chrom_length):
temp.append(random.randint(0,1))
pop.append(temp)
returnpop[1:]
#0.0coding:utf-80.0
#解码并计算值
importmath
defdecodechrom(pop,chrom_length):
temp=[]
foriinrange(len(pop)):
t=0
forjinrange(chrom_length):
t+=pop[i][j]*(math.pow(2,j))
temp.append(t)
returntemp
defcalobjValue(pop,chrom_length,max_value):
temp1=[]
obj_value=[]
temp1=decodechrom(pop,chrom_length)
foriinrange(len(temp1)):
x=temp1[i]*max_value/(math.pow(2,chrom_length)-1)
obj_value.append(10*math.sin(5*x)+7*math.cos(4*x))
returnobj_value
#0.0coding:utf-80.0
#淘汰(去除负值)
defcalfitValue(obj_value):
fit_value=[]
c_min=0
foriinrange(len(obj_value)):
if(obj_value[i]+c_min>0):
temp=c_min+obj_value[i]
else:
temp=0.0
fit_value.append(temp)
returnfit_value
#0.0coding:utf-80.0
#选择
importrandom
defsum(fit_value):
total=0
foriinrange(len(fit_value)):
total+=fit_value[i]
returntotal
defcumsum(fit_value):
foriinrange(len(fit_value)-2,-1,-1):
t=0
j=0
while(j<=i):
t+=fit_value[j]
j+=1
fit_value[i]=t
fit_value[len(fit_value)-1]=1
defselection(pop,fit_value):
newfit_value=[]
#适应度总和
total_fit=sum(fit_value)
foriinrange(len(fit_value)):
newfit_value.append(fit_value[i]/total_fit)
#计算累计概率
cumsum(newfit_value)
ms=[]
pop_len=len(pop)
foriinrange(pop_len):
ms.append(random.random())
ms.sort()
fitin=0
newin=0
newpop=pop
#转轮盘选择法
whilenewin<pop_len:
if(ms[newin]<newfit_value[fitin]):
newpop[newin]=pop[fitin]
newin=newin+1
else:
fitin=fitin+1
pop=newpop
- 以上代码主要进行了3个操作,首先是计算个体适应度总和,然后在计算各自的累积适应度。这两步都好理解,主要是第三步,转轮盘选择法。这一步首先是生成基因总数个0-1的小数,然后分别和各个基因的累积个体适应度进行比较。如果累积个体适应度大于随机数则进行保留,否则就淘汰。这一块的核心思想在于:一个基因的个体适应度越高,他所占据的累计适应度空隙就越大,也就是说他越容易被保留下来。
#0.0coding:utf-80.0
#交配
importrandom
defcrossover(pop,pc):
pop_len=len(pop)
foriinrange(pop_len-1):
if(random.random()<pc):
cpoint=random.randint(0,len(pop[0]))
temp1=[]
temp2=[]
temp1.extend(pop[i][0:cpoint])
temp1.extend(pop[i+1][cpoint:len(pop[i])])
temp2.extend(pop[i+1][0:cpoint])
temp2.extend(pop[i][cpoint:len(pop[i])])
pop[i]=temp1
pop[i+1]=temp2
- 变异:
#0.0coding:utf-80.0
#基因突变
importrandom
defmutation(pop,pm):
px=len(pop)
py=len(pop[0])
foriinrange(px):
if(random.random()<pm):
mpoint=random.randint(0,py-1)
if(pop[i][mpoint]==1):
pop[i][mpoint]=0
else:
pop[i][mpoint]=1
- 整个遗传算法的实现完成了,总的调用入口代码如下
#0.0coding:utf-80.0
importmatplotlib.pyplotasplt
importmath
fromselectionimportselection
fromcrossoverimportcrossover
frommutationimportmutation
frombestimportbest
print'y=10*math.sin(5*x)+7*math.cos(4*x)'
#计算2进制序列代表的数值
defb2d(b,max_value,chrom_length):
t=0
forjinrange(len(b)):
t+=b[j]*(math.pow(2,j))
t=t*max_value/(math.pow(2,chrom_length)-1)
returnt
pop_size=500#种群数量
max_value=10#基因中允许出现的最大值
chrom_length=10#染色体长度
pc=0.6#交配概率
pm=0.01#变异概率
results=[[]]#存储每一代的最优解,N个二元组
fit_value=[]#个体适应度
fit_mean=[]#平均适应度
#pop=[[0,1,0,1,0,1,0,1,0,1]foriinrange(pop_size)]
pop=geneEncoding(pop_size,chrom_length)
foriinrange(pop_size):
obj_value=calobjValue(pop,chrom_length,max_value)#个体评价
fit_value=calfitValue(obj_value)#淘汰
best_indivial,best_fit=best(pop,fit_value)#第一个存储最优的解,第二个存储最优基因
results.append([best_fit,b2d(best_indivial,max_value,chrom_length)])
selection(pop,fit_value)#新种群复制
crossover(pop,pc)#交配
mutation(pop,pm)#变异
results=results[1:]
results.sort()
X=[]
Y=[]
foriinrange(500):
X.append(i)
t=results[i][0]
Y.append(t)
plt.plot(X,Y)
plt.show()
- 最后调用了一下matplotlib包,把500代最优解的变化趋势表现出来。
其中genEncodeing是自定义的一个简单随机生成序列的函数,具体实现如下
[python]view plain
编码完成之后就是要进行个体评价,个体评价主要是计算各个编码出来的list的值以及对应带入目标式子的值。其实编码出来的就是一堆2进制list。这些2进制list每个都代表了一个数。其值的计算方式为转换为10进制,然后除以2的序列长度次方减一,也就是全一list的十进制减一。根据这个规则就能计算出所有list的值和带入要计算式子中的值,代码如下
[python]view plain
有了具体的值和对应的基因序列,然后进行一次淘汰,目的是淘汰掉一些不可能的坏值。这里由于是计算最大值,于是就淘汰负值就好了
[python]view plain
然后就是进行选择,这是整个遗传算法最核心的部分。选择实际上模拟生物遗传进化的优胜劣汰,让优秀的个体尽可能存活,让差的个体尽可能的淘汰。个体的好坏是取决于个体适应度。个体适应度越高,越容易被留下,个体适应度越低越容易被淘汰。具体的代码如下
[python]view plain
选择完后就是进行交配和变异,这个两个步骤很好理解。就是对基因序列进行改变,只不过改变的方式不一样
交配:
[python]view plain
[python]view plain
[python]view plain
完整代码可以在github查看
欢迎访问我的个人博客
阅读全文
Ⅷ Python 脚本之统计基因组文件中染色体长度及N碱基数目
re模块是Python中的正则表达式调用模块,在python中,通过将正则表达式内嵌集成re模块,程序员们可以直接调用来实现正则匹配。
正则表达式的大致匹配过程是:
re模块所支持的方法有如下:
其中,pattern为匹配模式,由re.compile生成,例如: pattern = re.compile(r'hello')
参数flag是匹配模式,取值可以使用按位或运算符’|’表示同时生效,如 re.I | re.M。可选值有:
注:以上七个方法中的flags同样是代表匹配模式的意思,如果在pattern生成时已经指明了flags,那么在下面的方法中就不需要传入这个参数了。
1. re.match(pattern, string[, flags])
2. re.search(pattern, string[, flags])
3. re.split(pattern, string[, maxsplit])
4. re.findall(pattern, string[, flags])
5. re.finditer(pattern, string[, flags])
6. re.sub(pattern, repl, string[, count])
7. re.subn(pattern, repl, string[, count])
描述
注意:该方法只能删除开头或是结尾的字符,不能删除中间部分的字符。
语法
参数
返回值
实例
以上实例输出结果如下:
参考
https://www.cnblogs.com/chengege/p/11190782.html
https://www.runoob.com/python/att-string-strip.html