[笔记]python对FASTA文件的处理
这学期选了生信的选修课—perl/python在生物信息学中的应用
把结课作业的代码整理出来主要是python对FASTA文件的读取和数据处理
FASTA文件数据处理
**
**
FASTA文件读取：
只含一个基因序列
将FASTA文件的基因序列读取到一个列表中，列表中的每个元素为每一行基因序列构成的字符串
f=open('/home/miaoyr/perl_practice/test1_file/DTNBP1.fasta')
ls=[]
for line in f:
if not line.startswith('>'):
ls.append(line.replace('\n','')) #去掉行尾的换行符真的很重要！
f.close()
python
含有多个序列信息
将数据保存在字典中，键为基因名，值为基因序列
f=open('/home/miaoyr/perl_practice/test4_file/3.fasta_seq.txt')
seq={}
for line in f:
if line.startswith('>'):
name=line.replace('>','').split()[0]
seq[name]=''
else:
seq[name]+=line.replace('\n','').strip()
f.close()
python
碱基数目和序列长度统计：
上一步将基因序列写入列表，为了使用str.count()函数（使用对象为字符串）直接计算，将全部基因序列合并为一个字符串_str
A=C=T=G=0
_str=''
for i in range(0,len(ls)):
_str+=ls[i]
print _str
A+=_str.count('A')
C+=_str.count('C')
T+=_str.count('T')
G+=_str.count('G')
print 'A:',A,'\n','C:',C,'\n','T:',T,'\n','G:',G,'\n'
python
输出互补序列：
html`
`
这个程序好像没有解决反向输出的问题。。
参考：seq[name][:-1] 没有试过不知道行不行ORZ
输出基因信息表中特定的基因信息
这里以本题为例：在基因的信息A文件中，找出B文件中存在的基因名的信息
文件B中只包含基因名，这里讲文件B按行读取，保存在一个列表中，列表中的每一个元素对应一个基因名
import os
get_name=[]
for _name in open('/home/miaoyr/perl_practice/test2_file/2.1name.txt'):
get_name.append(_name.replace('\n',''))
print get_name
f=open('2-1-output.txt','w')
for name in get_name:
m_name="'$3~/^"+name+"/ {print $0}'"
cmd="'awk' %s '/home/miaoyr/perl_practice/test2_file/2.1Predicted_Targets_Info.txt'" % m_name
output=os.popen(cmd).readlines()
f.writelines(output)
f.close()
python
由于当时脑子抽了，想用shell语言（程序都是在linux上编写的）python的交互来写。。。
用了os.popen()函数实现和shell语言的awk命令交互，多此一举了。。
变量传递的时候遇到了很大的问题，最后用字符串解决的。
实现shell和python交互可用的函数主要由 os.popen()和os.system()，二者有一些区别。
python比较简单的方法思路：
将文件A按行读取，用split()函数分割。
遍历B文件中的基因名，通过pattern正则表达式，在文件A读取的行中查找符合条件的基因信息并输出。
问题：程序速度较慢，我也想不出更好的办法了QAQ
以本题为例：同样是在基因的信息A文件中，找出B文件中存在的基因名的信息，这里用python直接处理。
get_id=[]
for _id in open('/home/miaoyr/perl_practice/test2_file/2.2gene_id'):
get_id.append(_id.replace('\n','')) #文件B中只包含基因名，直接将基因名按行读取，保存在一个列表中，每一个元素为一个基因名
info={} #定义一个字典，键为基因名，值用于储存基因信息
for l in open('/home/miaoyr/perl_practice/test2_file/2.2gene_DE_info'):
ls=l.replace('\n','').split('\t') #文件A中含有多种信息，用split()分开储存
name=ls[0]
info[name]=ls[1]
for _id in get_id: #遍历文件B中的基因名找到匹配的键
if _id==name:
print name,'\n',info[name]
python
简单的数据筛选可直接用linux命令实现
常用的包括 awk wc grep等
awk真的。。超级好用
以上题目对python的应用主要表现在
文件的读取和写入
输出格式化
一些库和函数的简单应用
我可能
还没有入门ORZ