python 根据染色体起始终止点坐标来获取碱基序列

时间 2019-12-01

标签 python 根据染色体起始终止坐标获取碱基序列栏目 Python 繁體版

原文原文链接

在工做当中，有时候咱们知到染色体编号以及染色体起始终止坐标，咱们想知道这段序列是什么样的碱基。
其一，咱们通常用去UCSC的genome browser里面去查询 ，其实也能够从UCSC的接口去解析网页，而后在提取序列信息
好比chr17:7676091,7676196 ，那么我只须要构造下面一个网页地址 http://genome.ucsc.edu/cgi-bin/das/hg38/dna?segment=chr17:7676091,7676196
而后 hg38能够更换成hg19，dna?segment= 后面能够按照标准格式更换，就能够返回咱们想要的序列了。如今对网页返回 一个xml格式的信息，用python解析一下


  1 import requests

 2 import re  3 from bs4 import BeautifulSoup  4 import xlwt  5 import xlrd  6 from xlutils.copy import copy  7 import os ,sys  8 #print(sys.path)
 9 cwd=os.getcwd()  10 
 11 
 12 
 13 def getHTMLText(url):  14     print("111111")  15     try:  16         header = {'user-agent': 'Mozilla/5.0'}  17         r = requests.get(url,headers = header,timeout = 30 )  18  r.raise_for_status()  19         r.encoding =r.apparent_encoding  20         print("get 222222222222222")  21         return r.text  22     except:  23         return ""
 24 
 25 def fillDNAList(dnalist, html):  26     # 使用正则表达式获取dna 序列的头文件
 27     match = re.search('SEQUENCE([\s\S]*?version="1.00")', html)  28     print("match ok")  29     if match:  30         dna_header = re.search('SEQUENCE([\s\S]*?version="1.00")', html)  31         #print('10====>', dna_header.group())
 32         #dna_header 存到列表
 33  dnalist.append(dna_header.group())  34 
 35     match = re.search('DNA.*?length="(\d*)"', html)  36     if match:  37         length_header= re.search('DNA.*?length="(\d*)"', html)  38         #print('11=====>', length_header.group())
 39  dnalist.append(length_header.group())  40 
 41     #使用 BeautifulSoup
 42     print("BeautifulSoup tag属性 获取dna标签属性的字符串部分")  43     soup = BeautifulSoup(html, 'html.parser')  44     tag = soup.dna  45     seq = soup.dna.string  46     seq = seq.replace('\n','').upper()  #换行符删除掉,转换成大写
 47     # seq 存到列表
 48  dnalist.append(seq)  49     print("final======>",dnalist)  50     return dnalist  51 
 52 def write_excell(dnalist,chrnum,pos):  53     head = '>hg19' + ' ' + dnalist[0] + dnalist[1]  54     f = xlwt.Workbook(encoding='utf-8', style_compression=0) #建立新的Excel（新的workbook）
 55     sheet = f.add_sheet('test8', cell_overwrite_ok=True) #建立新的表单
 56     # 先写第一行的头文件
 57  sheet.write(0,0,head)  58     #再从第二行开始写，每行写入50 个字符串
 59     dna = dnalist[2]  60     print('=====',dna,type(dna))  61     for i in range(0,len(dna),50):  62         sheet.write((int((i+1)/50))+1,0,dna[i:i+50])  63 
 64     out_file = 'chrmosome%s_%s.xls'% (chrnum,pos)  65  f.save(out_file)  66     out_file_dir = os.path.join(cwd, out_file)  67     return out_file_dir  68 
 69 def modify_excell(out_file_dir,chrnum,pos):  70     '''
 71  改Excel表（xlutils模块）  72  :return:  73     '''
 74     rb = xlrd.open_workbook(out_file_dir)  # 打开out_file.xls文件,建立工做簿实例对象
 75     sheet = rb.sheet_by_index(0)  76     nrow11 = sheet.cell_value(10, 0) #修改第11行第一列，索引是10，0
 77     # 根据须要截取原单元格里面的内容与须要添加的内容进行拼接
 78     new_nrow11 = '[' + nrow11  79     # 同理操做nrow12
 80     nrow12 = sheet.cell_value(11, 0)  81     new_nrow12 = nrow12 + ']'
 82 
 83     wb = copy(rb)  84     ws = wb.get_sheet(0)  85     # 往单元格中写入拼接后的新字符串内容
 86     ws.write(10,0,new_nrow11)  87     ws.write(11, 0, new_nrow12)  88     modify_file = 'new_chrmosome%s_%s.xls' % (chrnum,pos)  89  wb.save(modify_file)  90 
 91 def main():  92     hg19 = "hg19"
 93     chrnum = 17
 94     pos = 7676091
 95     start = pos - 500
 96     end = pos + 500
 97     position_DNA_list = []
 99     url = f"http://genome.ucsc.edu/cgi-bin/das/{hg19}/dna?segment=chr{chrnum}:{start},{end}"
100    
101     print(url) 102     html = getHTMLText(url) 103     dnalist = fillDNAList(position_DNA_list,html) 104     out_file_dir = write_excell(dnalist,chrnum,pos) 105  modify_excell(out_file_dir,chrnum,pos) 106 
107 main()

结果以下：html

http://genome.ucsc.edu/cgi-bin/das/hg19/dna?segment=chr17:7675591,7676591

match ok
BeautifulSoup tag属性 获取dna标签属性的字符串部分
final======> ['SEQUENCE id="chr17" start="7675591" stop="7676591" version="1.00"', 'DNA length="1001"', 'CCCAAGAGCCTTCAGTATACACATCAATAAAAATAATTTTAATTATTCTGATAAAAGATAAACATGAAAAGTTATGGTATGCAAAGTTGAATGACAACAACTGATACTATTTGAAATAATTGACAGAATTATATTCCGTAACAATTTATAAGCAAAGCCAAAAAAACAATGATCCCTTTGTTGAATGCACAGAACAAATCCATCTTGTCCACGGCTACTGAGCATGCCTGTGATCTCCAGGGGTCACTCAGGTTTGACTCAAAGGATCCAACAGCCTGTAGACCCTGTGCTTGAAGGCATGAGGGTCACCTCTGAGTTCACACTCACTAGTGTCCCTCCTTTCTTCAGAAAGCTAGGAACTGGGAAGACAAGGGGAAAATCAATCAAGGCCTGAGGTATGGGGCTGTAGGCTGGGAGGAAACTAACATTATTGAGAAGCTACTGATGTGAATACATTTCAATTACTACTCACATTGGTTTTTTGTTTGTTTGTTTGTTTGTTTGTTTGTTTGTTTTTTAAGACGGAGTTTTGCTCTCGTTGCCCAGGCTGGAGTGCAATGGAATGATCTAAGGTCACCACAACCTCCACCTCCCGGTTCAAGCAATTCTCCTGCCTCAGCCTCCCAAGTAGCTGGGACTACAGGCGTGTGCCACCACACCCAGCTAAGTTTGTATTTTTTTAGTAGAGACGGTGTTTCACCATGTTGGTCAGGCTGGTCTCGAACTCCTGACCTCAAGTGATCCACCCACCTCGGCCTCCCAAAGTGCTGGGATTATAGGCATGAGCCACCACACCCAGCCTCACGTTGGTTTTTGAGATGGATTTTATTGCCATTTTGTACACAAAAAGGTCAAAACTCAGTGAGGTGAATTGACATGACAGTAAGTGAAAGAACTACTATCTGATTGGGGGTCTTCTGCCGCCTGCTCTGGGACTCTTTCTGCTATGACATGAAGGACATTGGCAACCCCAGTCCTTGCAGATTTCTTTCACTGTGTGC']
===== CCCAAGAGCCTTCAGTATACACATCAATAAAAATAATTTTAATTATTCTGATAAAAGATAAACATGAAAAGTTATGGTATGCAAAGTTGAATGACAACAACTGATACTATTTGAAATAATTGACAGAATTATATTCCGTAACAATTTATAAGCAAAGCCAAAAAAACAATGATCCCTTTGTTGAATGCACAGAACAAATCCATCTTGTCCACGGCTACTGAGCATGCCTGTGATCTCCAGGGGTCACTCAGGTTTGACTCAAAGGATCCAACAGCCTGTAGACCCTGTGCTTGAAGGCATGAGGGTCACCTCTGAGTTCACACTCACTAGTGTCCCTCCTTTCTTCAGAAAGCTAGGAACTGGGAAGACAAGGGGAAAATCAATCAAGGCCTGAGGTATGGGGCTGTAGGCTGGGAGGAAACTAACATTATTGAGAAGCTACTGATGTGAATACATTTCAATTACTACTCACATTGGTTTTTTGTTTGTTTGTTTGTTTGTTTGTTTGTTTGTTTTTTAAGACGGAGTTTTGCTCTCGTTGCCCAGGCTGGAGTGCAATGGAATGATCTAAGGTCACCACAACCTCCACCTCCCGGTTCAAGCAATTCTCCTGCCTCAGCCTCCCAAGTAGCTGGGACTACAGGCGTGTGCCACCACACCCAGCTAAGTTTGTATTTTTTTAGTAGAGACGGTGTTTCACCATGTTGGTCAGGCTGGTCTCGAACTCCTGACCTCAAGTGATCCACCCACCTCGGCCTCCCAAAGTGCTGGGATTATAGGCATGAGCCACCACACCCAGCCTCACGTTGGTTTTTGAGATGGATTTTATTGCCATTTTGTACACAAAAAGGTCAAAACTCAGTGAGGTGAATTGACATGACAGTAAGTGAAAGAACTACTATCTGATTGGGGGTCTTCTGCCGCCTGCTCTGGGACTCTTTCTGCTATGACATGAAGGACATTGGCAACCCCAGTCCTTGCAGATTTCTTTCACTGTGTGC <class 'str'>