Biopython으로 Swiss-prot과 ExPASy 데이터베이스 탐방하기

근데 스위스프롯은 긁어와서 저장 안됨? 

---

첫빠따는 파싱이 국룰이지

파싱 방법이 네 가지가 있는데 gZIP은 생략. gZIP 파일을 못 구했다. 

 

handle=open('/home/koreanraichu/Q63HQ2.txt')
print(handle)

import gzip
handle = gzip.open("myswissprotfile.dat.gz", "rt")

from urllib.request import urlopen
url = "https://raw.githubusercontent.com/biopython/biopython/master/Tests/SwissProt/F2CXE6.txt"
handle = urlopen(url)
print(handle)

from Bio import ExPASy
handle = ExPASy.get_sprot_raw("Q63HQ2")
print(handle)

위에서부터 순서대로 

1. 갖고 있는 파일 열기 
2. gZIP 핸들링해서 열기(...)
3. URL 입력해서 갖고오기
4. ExPASy에서 가져오기

이상이다. 

네? 출력이요? 그냥 파싱에 포문 빠져서 그런갑다 하자... 

 

실제로 가져와봤습니다

from Bio import ExPASy
from Bio import SwissProt
# 일단 여기서는 4번, ExPASY에 접속해서 가져와보자.
handle = ExPASy.get_sprot_raw("Q63HQ2")
record = SwissProt.read(handle)
print(record.description)

RecName: Full=Pikachurin; AltName: Full=Agrin-like protein; AltName: Full=EGF-like, fibronectin type-III and laminin G-like domain-containing protein; Flags: Precursor;

참고로 위 방법 중 4번을 선택했다. 네? 결과에 피카츄요? 내가 뭘 잘못 본 거 아니냐고요? 정상인데요? 단백질 이름이 피카츄린이니까 피카츄가 들어가지. 

 

for ref in record.references:
    print("authors:", ref.authors)
    print("title:", ref.title)
    # 아마도 관련 참고문헌 정보?

이런것도 된다. 

 

Prosite 레코드 파싱하기

Prosite는 단백질 도메인, 단백질에서 기능하는 부분과 인식 프로파일 정보, functional family에 대한 정보... 가 들어가 있는데... 이거 CATH 아니냐? CATH와는 다르다 CATH와는! 

 

ftp://ftp.expasy.org/databases/prosite/prosite.dat 

dat파일은 여기서 받자. 리눅스 유저는 터미널에서 wget -i 들어가면 되고... 윈도우는...... 여기다 올려주려고 했는데 10MB 넘어가서 안됨... 리눅스로 갈아타자 뭔소리여 들어가면 텍스트 뜨는데 걍 메모장에 복붙하고 prosite.dat으로 이름 지어주면 된다. 깃헙에요? 저게 올라감? ㄷㄷ 

 

from Bio.ExPASy import Prosite
handle = open("/home/koreanraichu/prosite.dat")
records = Prosite.parse(handle)
print(records)

이렇게 가져오면 된다는데 결과가 또 뭔가 괴랄하다. (dat파일 받는 링크는 위에 올려드림) 

from Bio.ExPASy import Prosite
from urllib.request import urlopen
handle=open("/home/koreanraichu/prosite.dat")
records = Prosite.parse(handle)
record = next(records)
print(record.accession)

이것까지 주면 원하는 정보만 볼 수는 있는데... 어? next? 

 

from Bio.ExPASy import Prosite
from urllib.request import urlopen
handle=open("/home/koreanraichu/prosite.dat")
records = Prosite.parse(handle)
for i in range(0,5):
    record = next(records)
    print(record.accession)

For문 마려우셨죠? 됩니다. 

from Bio.ExPASy import Prosite
from urllib.request import urlopen
handle=open("/home/koreanraichu/prosite.dat")
records = Prosite.parse(handle)
i=0
while i < 5:
    record = next(records)
    print(record.accession)
    i=i+1

아 For문 되면 while도 됨 

저 코드 각개로 해보는 걸 추천한다. 연달아서 하면 1-(2~6)-(7~11) 이렇게 나온다. 

 

n=0
for record in records: n+=1
print(n)

이 코드로 레코드 숫자도 볼 수 있다. 

 

효소 레코드 파싱하기

from Bio.ExPASy import Enzyme
with open("/home/koreanraichu/RuBisCO.txt") as handle:
    record = Enzyme.read(handle)
    print(record['ID'])

4.1.1.39

위 번호는 EC번호로, 효소가 반응을 '어떤 방식으로' 촉매하는가에 따라 번호가 붙는다. 

 

from Bio.ExPASy import Enzyme
with open("/home/koreanraichu/RuBisCO.txt") as handle:
    record = Enzyme.read(handle)
    print(record['ID']) # EC no.
    print(record['DE']) # description
    print(record['AN']) # 대충 synonyms같은건가? 뭐 얘 이렇게도 불러요 이런거
    print(record["CA"]) # 촉매하는 반응(오 이거 식으로 나온다)
    print(record["PR"]) # 이건 모르겠다... 데이터베이스 번호인가... 
    print(record["CC"]) # 아마도 뭐 하는 효소인가에 대한 설명인 듯 
    print(record['DR']) # 뭔진 모르겠지만 일단 잘못했어요... 뭐가 되게 많이떴는데 넘파이 마려웠음

이거 다 띄웠더니 뭐가 많아서 결과는 생략. 작동 하기는 한다. 

 

ExPASy 서버 접속해보기

이제 본격적으로 서버 털러 가 봅시다. 일단 ExPASy 관련해서는 네 가지 쿼리가 있는데 

get_prodoc_entry: To download a Prosite documentation record in HTML format
get_prosite_entry: To download a Prosite record in HTML format
get_prosite_raw: To download a Prosite or Prosite documentation record in raw format
get_sprot_raw: To download a Swiss-Prot record in raw format

이렇게 네 개가 있다. 

다운로드래서 행복회로 태웠는데 저장 안되더만 swissprot은... 

 

Swissprot 검색하기

참고로 이거 Uniprot 번호 먹힌다. 

 

from Bio import ExPASy
from Bio import SwissProt
accessions = ["O23729", "O23730", "O23731"]
records = []
# 저 세개에 대해 찾고싶은갑다.
for accession in accessions:
    handle = ExPASy.get_sprot_raw(accession)
    record = SwissProt.read(handle)
    records.append(record)
    print(records) # IDLE 아니라 이거 있어야 나옴

[<Bio.SwissProt.Record object at 0x7f1abf171f70>, <Bio.SwissProt.Record object at 0x7f1abf9a3d00>, <Bio.SwissProt.Record object at 0x7f1abf117550>]

이러시는 이유가 있으실 것 아니예요. 

 

from Bio import ExPASy
from Bio import SwissProt
accessions = ["C3RT18"] # 단백질 이름은 신경쓰지 맙시다 
records = []
# 저 세개에 대해 찾고싶은갑다.
for accession in accessions:
    handle = ExPASy.get_sprot_raw(accession)
    record = SwissProt.read(handle)
    records.append(record)
print(records) #IDLE 아니라서 이거 써야 나온다

단식 검색도 된다. 

 

>>> for accession in accessions:
...     handle = ExPASy.get_sprot_raw(accession)
...     try:
...         record = SwissProt.read(handle)
...     except ValueException:
...         print("WARNING: Accession %s not found" % accession)
...     records.append(record)

검색했는데 결과가 없을 때 ValueError가 여러분을 반길 것이다. 그럴때는 당황하지 말고 에외 처리를 해 보자. 

 

Prosite 검색하기

from Bio import ExPASy
handle = ExPASy.get_prosite_raw("PS00001")
text = handle.read()
print(text)
# 가장 원시적인 긁는 방법

ID   ASN_GLYCOSYLATION; PATTERN.
AC   PS00001;
DT   01-APR-1990 CREATED; 01-APR-1990 DATA UPDATE; 01-APR-1990 INFO UPDATE.
DE   N-glycosylation site.
PA   N-{P}-[ST]-{P}.
CC   /SITE=1,carbohydrate;
CC   /SKIP-FLAG=TRUE;
CC   /VERSION=1;
PR   PRU00498;
DO   PDOC00001;
//

가장 원시적인 긁어오기 신공. 아 근데 개구린 하려고 했는데 prosite에는 개구린이 없다... (개구린: 단백질 이름) 

 

from Bio import Prosite
handle = ExPASy.get_prosite_raw("PS51036")
record = Prosite.read(handle)
print(record)

이건 된다고 올려놨는데 첫 줄에서 ImportError: cannot import name 'Prosite' from 'Bio' (/home/koreanraichu/PycharmProjects/pythonProject/venv/lib/python3.8/site-packages/Bio/init.py가 나를 반긴다. 넌 또 왜그러냐. 

 

from Bio.ExPASy import Prodoc
handle = ExPASy.get_prosite_raw("PDOC00001")
record = Prodoc.read(handle)
print(record)

Prodoc을 가져올 수도 있다. 

 

handle = ExPASy.get_prosite_entry("PS00001")
html = handle.read()
with open("myprositerecord.html", "w") as out_handle:
    out_handle.write(html)
# HTML format으로 다운로드 받을 수 있다.

handle = ExPASy.get_prodoc_entry("PDOC51036")
html = handle.read()
with open("myprositerecord2.html", "w") as out_handle:
    out_handle.write(html)
# 얘는 prodoc 다운로드 하는 코드

HTML 형식으로 저장도 된다. 

이런 식으로 생성되고(경로 지정을 안 했더니 파이참 폴더에 박아버리네...)

위는 prosite, 아래는 prodoc. 이런 식으로 뜬다. 

 

Prosite database 스캔하기

스캔 아니고 시퀀스 검색 아니냐. 

 

from Bio.ExPASy import ScanProsite
seq="MEHKEVVLLLLLFLKSGQGEPLDDYVNTQGASLFSVTKKQLGAGSIEECAAKCEEDEEFTCRAFQYHSKEQQCVIMAENRKSSIIIRMRDVVLFEKKVYLSECKTGNGKNYRGTMSKTKN"
handle = ScanProsite.scan(seq=seq)
result = ScanProsite.read(handle)
print(result)

[{'sequence_ac': 'USERSEQ1', 'start': 16, 'stop': 98, 'signature_ac': 'PS50948', 'score': '8.873', 'level': '0'}]

어? 딕셔너리네? 

 

print(result.n_seq)
print(result.n_match)
print(result[0])

1
1
{'sequence_ac': 'USERSEQ1', 'start': 16, 'stop': 98, 'signature_ac': 'PS50948', 'score': '8.873', 'level': '0'}

뭐 이렇다는데... 이거 예제꺼 긁어왔는데 왜 결과가 다르냐...