python编程实战手记_生信编程实战第9题（python）

image.png

这个题目不难，但是我想说明的是

大的数据集和小的数据集的脚本很多时候是不一样的

比如这道题，如果使用的是小的数据集>chr1

ATCGTCGaaAATGAANccNNttGTA

AGGTCTNAAccAAttGggG

>chr2

ATCGAATGATCGANNNGccTA

AGGTCTNAAAAGG

>chr3

ATCGTCGANNNGTAATggGA

AGGTCTNAAAAGG

>chr4

ATCGTCaaaGANNAATGANGgggTA

我可以用构建字典的方式，将想要查找的染色体和所有该染色体所有碱基的坐标以及对应碱基构建字典import sys

args=sys.argv

filename=args[1]import collections

chr_num=args[2]

base_num=args[3]

aDict=collections.OrderedDict()

num=0with open(filename) as fh:

chrome=">"+chr_num  for line in fh:

line=line.strip()     if line == chrome :

aDict[chrome]=collections.OrderedDict()

num=1

continue

if num:        if line.startswith(">"):           break

for i in line:

aDict[chrome][num]=i

num+=1for k,v in aDict[chrome].items():    if k == int(base_num):

print(aDict[chrome][k])

这样，如果想得到chr2的坐标为8，9，10的碱基python3 chr_base.py test.fa chr2 8Gpython3 chr_base.py test.fa chr2 9Apython3 chr_base.py test.fa chr2 10T

但是，如果对于hg38这种大的数据集，这种方法显然是不合适的

构建字典的过程会消耗大量的内存

所以尽量避免使用字典import sys

args=sys.argv

filename=args[1]import collections

chr_num=args[2]

base_num=args[3]

up_stream=int(base_num)-1down_stream=int(base_num)+1num=0stop=0with open(filename) as fh:

chrome=">"+chr_num  for line in fh:

line=line.strip()     if line == chrome :

num=1

continue

if num:        if line.startswith(">"):           break

for i in line:           if num == up_stream:

print(i)

stop=stop+1

if num == int(base_num):

print(i)

stop=stop+1

if num == down_stream:

print(i)

stop=stop+1

num+=1

if stop == 3:         break

这里需要注意的是输入的坐标位置一定要转成int

还有就是，我设置stop这个值是为了提高速度，也就是得到了结果后，后面的就不要遍历了。time python3 chr_base_bigdata.py hg38.fa chr5 8397384GACreal    0m21.503suser    0m21.008ssys 0m0.488s

如果不加stopimport sys

args=sys.argv

filename=args[1]import collections

chr_num=args[2]

base_num=args[3]

up_stream=int(base_num)-1down_stream=int(base_num)+1num=0#stop=0with open(filename) as fh:

chrome=">"+chr_num  for line in fh:

line=line.strip()     if line == chrome :

num=1

continue

if num:        if line.startswith(">"):           break

for i in line:           if num == up_stream:

print(i)#              stop=stop+1

if num == int(base_num):

print(i)#              stop=stop+1

if num == down_stream:

print(i)#              stop=stop+1

num+=1#     if stop == 3:#         breaktime python3 chr_base_bigdata.py hg38.fa chr5 8397384GACreal    1m47.082suser    1m46.152ssys 0m0.692s

很明显，一个21s,另外一个1min46s

作者：天秤座的机器狗

链接：https://www.jianshu.com/p/e88190b1e1ed