Swiss-Prot 與 ExPASy

解析 Swiss-Prot 檔案

Swiss-Prot (https://web.expasy.org/docs/swiss-prot_guideline.html) 是一個經過人工整理的蛋白質序列資料庫。Biopython 可以解析「純文字」Swiss-Prot 檔案格式，該格式仍用於結合 Swiss-Prot、TrEMBL 和 PIR-PSD 的 UniProt Knowledgebase。

儘管以下內容我們著重於較舊的人類可讀純文字格式，Bio.SeqIO 可以讀取此格式和用於帶註解的蛋白質序列的較新 UniProt XML 檔案格式。

解析 Swiss-Prot 紀錄

在第從網路解析 SwissProt 序列節中，我們描述了如何將 Swiss-Prot 記錄的序列提取為 SeqRecord 物件。或者，您可以將 Swiss-Prot 紀錄儲存在 Bio.SwissProt.Record 物件中，它實際上儲存了 Swiss-Prot 紀錄中包含的完整資訊。在本節中，我們描述如何從 Swiss-Prot 檔案中提取 Bio.SwissProt.Record 物件。

要解析 Swiss-Prot 紀錄，我們首先取得 Swiss-Prot 紀錄的控制代碼。有幾種方法可以做到這一點，具體取決於 Swiss-Prot 紀錄的儲存位置和方式

在本地開啟 Swiss-Prot 檔案

>>> handle = open("SwissProt/F2CXE6.txt")

開啟 gzipped 的 Swiss-Prot 檔案

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

透過網際網路開啟 Swiss-Prot 檔案

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

在呼叫 read 之前開啟儲存在網際網路上的檔案。

透過網際網路從 ExPASy 資料庫開啟 Swiss-Prot 檔案（請參閱擷取 Swiss-Prot 紀錄一節）
```
>>> from Bio import ExPASy
>>> handle = ExPASy.get_sprot_raw("F2CXE6")
```

關鍵是對於解析器來說，控制代碼是如何建立的並不重要，只要它指向 Swiss-Prot 格式的資料即可。如果控制代碼以二進制模式開啟，則解析器將自動將資料解碼為 ASCII（Swiss-Prot 使用的編碼）。

我們可以按照第從網路解析 SwissProt 序列節中的描述，使用 Bio.SeqIO 來取得檔案格式獨立的 SeqRecord 物件。或者，我們可以使用 Bio.SwissProt 取得 Bio.SwissProt.Record 物件，這些物件與底層檔案格式更為匹配。

要從控制代碼中讀取一個 Swiss-Prot 紀錄，我們使用函數 read()

>>> from Bio import SwissProt
>>> record = SwissProt.read(handle)

如果控制代碼指向剛好一個 Swiss-Prot 紀錄，則應使用此函數。如果未找到 Swiss-Prot 紀錄，以及找到多個紀錄，則會引發 ValueError。

我們現在可以列印出有關此紀錄的一些資訊

>>> print(record.description)
SubName: Full=Plasma membrane intrinsic protein {ECO:0000313|EMBL:BAN04711.1}; SubName: Full=Predicted protein {ECO:0000313|EMBL:BAJ87517.1};
>>> for ref in record.references:
...     print("authors:", ref.authors)
...     print("title:", ref.title)
...
authors: Matsumoto T., Tanaka T., Sakai H., Amano N., Kanamori H., Kurita K., Kikuta A., Kamiya K., Yamamoto M., Ikawa H., Fujii N., Hori K., Itoh T., Sato K.
title: Comprehensive sequence analysis of 24,783 barley full-length cDNAs derived from 12 clone libraries.
authors: Shibasaka M., Sasano S., Utsugi S., Katsuhara M.
title: Functional characterization of a novel plasma membrane intrinsic protein2 in barley.
authors: Shibasaka M., Katsuhara M., Sasano S.
title:
>>> print(record.organism_classification)
['Eukaryota', 'Viridiplantae', 'Streptophyta', 'Embryophyta', 'Tracheophyta', 'Spermatophyta', 'Magnoliophyta', 'Liliopsida', 'Poales', 'Poaceae', 'BEP clade', 'Pooideae', 'Triticeae', 'Hordeum']

要解析包含多個 Swiss-Prot 紀錄的檔案，我們改為使用 parse 函數。此函數允許我們迭代檔案中的紀錄。

例如，讓我們解析完整的 Swiss-Prot 資料庫並收集所有描述。您可以從 ExPASy FTP 站點下載此內容，作為單一 gzipped 檔案 uniprot_sprot.dat.gz（約 300MB）。這是一個壓縮檔案，其中包含單一檔案 uniprot_sprot.dat（超過 1.5GB）。

如本節開頭所述，您可以使用 Python 程式庫 gzip 開啟和解壓縮 .gz 檔案，如下所示

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

但是，解壓縮大型檔案需要時間，而且每次您以這種方式開啟檔案進行讀取時，都必須即時解壓縮。因此，如果您可以空出磁碟空間，那麼如果您先將檔案解壓縮到磁碟，以取得裡面的 uniprot_sprot.dat 檔案，從長遠來看，您將節省時間。然後您可以像平常一樣開啟檔案進行讀取

>>> handle = open("uniprot_sprot.dat")

截至 2009 年 6 月，從 ExPASy 下載的完整 Swiss-Prot 資料庫包含 468851 個 Swiss-Prot 紀錄。使用列表推導式建立紀錄描述列表的一種簡潔方式是

>>> from Bio import SwissProt
>>> handle = open("uniprot_sprot.dat")
>>> descriptions = [record.description for record in SwissProt.parse(handle)]
>>> len(descriptions)
468851
>>> descriptions[:5]
['RecName: Full=Protein MGF 100-1R;',
 'RecName: Full=Protein MGF 100-1R;',
 'RecName: Full=Protein MGF 100-1R;',
 'RecName: Full=Protein MGF 100-1R;',
 'RecName: Full=Protein MGF 100-2L;']

或者，使用 for 迴圈遍歷紀錄迭代器

>>> from Bio import SwissProt
>>> descriptions = []
>>> handle = open("uniprot_sprot.dat")
>>> for record in SwissProt.parse(handle):
...     descriptions.append(record.description)
...
>>> len(descriptions)
468851

由於這是一個如此大型的輸入檔案，無論哪種方式，在我的新桌上型電腦上都需要大約 11 分鐘（使用未壓縮的 uniprot_sprot.dat 檔案作為輸入）。

從 Swiss-Prot 紀錄中提取您想要的任何類型的資訊都同樣容易。若要查看 Swiss-Prot 紀錄的成員，請使用

>>> dir(record)
['__doc__', '__init__', '__module__', 'accessions', 'annotation_update',
'comments', 'created', 'cross_references', 'data_class', 'description',
'entry_name', 'features', 'gene_name', 'host_organism', 'keywords',
'molecule_type', 'organelle', 'organism', 'organism_classification',
'references', 'seqinfo', 'sequence', 'sequence_length',
'sequence_update', 'taxonomy_id']

解析 Prosite 紀錄

Prosite 是一個包含蛋白質結構域、蛋白質家族、功能位點以及用於識別它們的模式和設定檔的資料庫。Prosite 與 Swiss-Prot 並行開發。在 Biopython 中，Prosite 紀錄由 Bio.ExPASy.Prosite.Record 類別表示，其成員對應於 Prosite 紀錄中的不同欄位。

一般來說，一個 Prosite 檔案可以包含多個 Prosite 紀錄。例如，整組 Prosite 紀錄可以從 ExPASy FTP 站點下載為單一檔案 (prosite.dat)，其中包含 2073 個紀錄（2007 年 12 月 4 日發布的 20.24 版）。要解析此類檔案，我們再次使用迭代器

>>> from Bio.ExPASy import Prosite
>>> handle = open("myprositefile.dat")
>>> records = Prosite.parse(handle)

我們現在可以一次取得一個紀錄並列印出一些資訊。例如，使用包含完整 Prosite 資料庫的檔案，我們會發現

>>> from Bio.ExPASy import Prosite
>>> handle = open("prosite.dat")
>>> records = Prosite.parse(handle)
>>> record = next(records)
>>> record.accession
'PS00001'
>>> record.name
'ASN_GLYCOSYLATION'
>>> record.pdoc
'PDOC00001'
>>> record = next(records)
>>> record.accession
'PS00004'
>>> record.name
'CAMP_PHOSPHO_SITE'
>>> record.pdoc
'PDOC00004'
>>> record = next(records)
>>> record.accession
'PS00005'
>>> record.name
'PKC_PHOSPHO_SITE'
>>> record.pdoc
'PDOC00005'

等等。如果您有興趣知道有多少 Prosite 紀錄，可以使用

>>> from Bio.ExPASy import Prosite
>>> handle = open("prosite.dat")
>>> records = Prosite.parse(handle)
>>> n = 0
>>> for record in records:
...     n += 1
...
>>> n
2073

要從控制代碼中讀取剛好一個 Prosite，您可以使用 read 函數

>>> from Bio.ExPASy import Prosite
>>> handle = open("mysingleprositerecord.dat")
>>> record = Prosite.read(handle)

如果未找到 Prosite 紀錄，以及找到多個 Prosite 紀錄，此函數會引發 ValueError。

解析 Prosite 文件紀錄

在上面的 Prosite 範例中，record.pdoc 存取號碼 'PDOC00001'、'PDOC00004'、'PDOC00005' 等等是指 Prosite 文件。Prosite 文件紀錄可以從 ExPASy 取得，分為個別檔案和一個包含所有 Prosite 文件紀錄的檔案 (prosite.doc)。

我們使用 Bio.ExPASy.Prodoc 中的解析器來解析 Prosite 文件紀錄。例如，若要建立 Prosite 文件紀錄的所有存取號碼的列表，可以使用

>>> from Bio.ExPASy import Prodoc
>>> handle = open("prosite.doc")
>>> records = Prodoc.parse(handle)
>>> accessions = [record.accession for record in records]

同樣，提供了一個 read() 函數，用於從控制代碼中讀取剛好一個 Prosite 文件紀錄。

解析 Enzyme 紀錄

ExPASy 的 Enzyme 資料庫是關於酵素命名法資訊的儲存庫。典型的 Enzyme 紀錄如下所示

ID   3.1.1.34
DE   Lipoprotein lipase.
AN   Clearing factor lipase.
AN   Diacylglycerol lipase.
AN   Diglyceride lipase.
CA   Triacylglycerol + H(2)O = diacylglycerol + a carboxylate.
CC   -!- Hydrolyzes triacylglycerols in chylomicrons and very low-density
CC       lipoproteins (VLDL).
CC   -!- Also hydrolyzes diacylglycerol.
PR   PROSITE; PDOC00110;
DR   P11151, LIPL_BOVIN ;  P11153, LIPL_CAVPO ;  P11602, LIPL_CHICK ;
DR   P55031, LIPL_FELCA ;  P06858, LIPL_HUMAN ;  P11152, LIPL_MOUSE ;
DR   O46647, LIPL_MUSVI ;  P49060, LIPL_PAPAN ;  P49923, LIPL_PIG   ;
DR   Q06000, LIPL_RAT   ;  Q29524, LIPL_SHEEP ;
//

在此範例中，第一行顯示脂蛋白脂肪酶的 EC（酵素委員會）編號（第二行）。脂蛋白脂肪酶的其他名稱為「清除因子脂肪酶」、「二酰基甘油脂肪酶」和「雙甘油脂肪酶」（第 3 行至第 5 行）。以「CA」開頭的行顯示此酵素的催化活性。「CC」開頭的是註解行。「PR」行顯示對 Prosite 文件紀錄的參照，「DR」行顯示對 Swiss-Prot 紀錄的參照。這些項目不一定都出現在 Enzyme 紀錄中。

在 Biopython 中，Enzyme 紀錄由 Bio.ExPASy.Enzyme.Record 類別表示。此紀錄衍生自 Python 字典，並且具有對應於 Enzyme 檔案中使用的雙字母代碼的鍵。若要讀取包含一個 Enzyme 紀錄的 Enzyme 檔案，請使用 Bio.ExPASy.Enzyme 中的 read 函數

>>> from Bio.ExPASy import Enzyme
>>> with open("lipoprotein.txt") as handle:
...     record = Enzyme.read(handle)
...
>>> record["ID"]
'3.1.1.34'
>>> record["DE"]
'Lipoprotein lipase.'
>>> record["AN"]
['Clearing factor lipase.', 'Diacylglycerol lipase.', 'Diglyceride lipase.']
>>> record["CA"]
'Triacylglycerol + H(2)O = diacylglycerol + a carboxylate.'
>>> record["PR"]
['PDOC00110']

>>> record["CC"]
['Hydrolyzes triacylglycerols in chylomicrons and very low-density lipoproteins
(VLDL).', 'Also hydrolyzes diacylglycerol.']
>>> record["DR"]
[['P11151', 'LIPL_BOVIN'], ['P11153', 'LIPL_CAVPO'], ['P11602', 'LIPL_CHICK'],
['P55031', 'LIPL_FELCA'], ['P06858', 'LIPL_HUMAN'], ['P11152', 'LIPL_MOUSE'],
['O46647', 'LIPL_MUSVI'], ['P49060', 'LIPL_PAPAN'], ['P49923', 'LIPL_PIG'],
['Q06000', 'LIPL_RAT'], ['Q29524', 'LIPL_SHEEP']]

如果未找到 Enzyme 紀錄，以及找到多個 Enzyme 紀錄，read 函數會引發 ValueError。

整組 Enzyme 紀錄可以從 ExPASy FTP 站點下載為單一檔案 (enzyme.dat)，其中包含 4877 個紀錄（2009 年 3 月 3 日發布）。若要解析包含多個 Enzyme 紀錄的此類檔案，請使用 Bio.ExPASy.Enzyme 中的 parse 函數來取得迭代器

>>> from Bio.ExPASy import Enzyme
>>> handle = open("enzyme.dat")
>>> records = Enzyme.parse(handle)

我們現在可以一次迭代這些紀錄。例如，我們可以建立所有有 Enzyme 紀錄可用的 EC 編號的列表

>>> ecnumbers = [record["ID"] for record in records]

存取 ExPASy 伺服器

Swiss-Prot、Prosite 和 Prosite 文件紀錄可以從 ExPASy 網站 https://www.expasy.org 下載。ExPASy 提供四種類型的查詢

get_prodoc_entry: 以下載 HTML 格式的 Prosite 文件記錄
get_prosite_entry: 以下載 HTML 格式的 Prosite 記錄
get_prosite_raw: 以下載原始格式的 Prosite 或 Prosite 文件記錄
get_sprot_raw: 以下載原始格式的 Swiss-Prot 記錄

要從 Python 腳本存取此 Web 伺服器，我們使用 Bio.ExPASy 模組。

檢索 Swiss-Prot 記錄

假設我們正在研究蘭花的查爾酮合成酶（請參閱使用範例章節，了解研究蘭花有趣之處的一些理由）。查爾酮合成酶參與植物中的類黃酮生物合成，而類黃酮會產生許多很棒的東西，例如色素顏色和紫外線保護劑。

如果您在 Swiss-Prot 上進行搜尋，您可以找到三種蘭花的查爾酮合成酶蛋白質，ID 編號分別是 O23729、O23730、O23731。現在，讓我們編寫一個腳本來抓取這些蛋白質，並解析出一些有趣的資訊。

首先，我們使用 Bio.ExPASy 的 get_sprot_raw() 函式來抓取記錄。這個函式非常棒，因為您可以提供一個 ID，然後取得一個原始文字記錄的控制代碼（沒有 HTML 來干擾！）。然後，我們可以使用 Bio.SwissProt.read 來提取 Swiss-Prot 記錄，或是使用 Bio.SeqIO.read 來取得 SeqRecord。以下程式碼完成了我剛才寫的步驟

>>> 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)
...

如果您提供給 ExPASy.get_sprot_raw 的登錄編號不存在，則 SwissProt.read(handle) 將會引發 ValueError。您可以捕獲 ValueException 例外來偵測無效的登錄編號

>>> 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)
...

使用 UniProt 搜尋

現在，您可能會說我事先知道這些記錄的登錄編號。確實，get_sprot_raw() 需要條目名稱或登錄編號。當您手邊沒有這些資訊時，可以使用 UniProt 來搜尋它們。您也可以使用 UniProt 套件以程式化的方式搜尋蛋白質。

例如，讓我們搜尋來自特定生物（生物 ID：2697049）且已審查的蛋白質。我們可以使用以下程式碼來執行此操作

>>> from Bio import UniProt

>>> query = "(organism_id:2697049) AND (reviewed:true)"
>>> results = list(UniProt.search(query))

UniProt.search 方法會傳回搜尋結果的迭代器。該迭代器會一次傳回一個結果，並在需要時從 UniProt 中提取更多結果，直到傳回所有結果。我們可以從這個迭代器有效率地建立此特定查詢的清單，因為此查詢只會傳回幾個結果（撰寫時為 17 個）。

讓我們嘗試一個會傳回更多結果的搜尋。在撰寫本文時，查詢「Insulin AND (reviewed:true)」共有 5,147 個結果。我們可以使用切片來取得前 50 個結果的清單。

>>> from Bio import UniProt
>>> from itertools import islice

>>> query = "Insulin AND (reviewed:true)"
>>> results = UniProt.search(query, batch_size=50)[:50]

您可以使用 len 方法取得搜尋結果的總數（不論批次大小為何）

>>> from Bio import UniProt

>>> query = "Insulin AND (reviewed:true)"
>>> result_iterator = UniProt.search(query, batch_size=0)
>>> len(result_iterator)
5147

檢索 Prosite 和 Prosite 文件記錄

Prosite 和 Prosite 文件記錄可以 HTML 格式或原始格式檢索。若要使用 Biopython 解析 Prosite 和 Prosite 文件記錄，您應該以原始格式檢索記錄。但是，對於其他用途，您可能對 HTML 格式的這些記錄感興趣。

若要以原始格式檢索 Prosite 或 Prosite 文件記錄，請使用 get_prosite_raw()。例如，以下載 Prosite 記錄並以原始文字格式列印，請使用

>>> from Bio import ExPASy
>>> handle = ExPASy.get_prosite_raw("PS00001")
>>> text = handle.read()
>>> print(text)

若要檢索 Prosite 記錄並將其解析為 Bio.Prosite.Record 物件，請使用

>>> from Bio import ExPASy
>>> from Bio import Prosite
>>> handle = ExPASy.get_prosite_raw("PS00001")
>>> record = Prosite.read(handle)

相同的函式可用於檢索 Prosite 文件記錄，並將其解析為 Bio.ExPASy.Prodoc.Record 物件

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

對於不存在的登錄編號，ExPASy.get_prosite_raw 會傳回一個空字串的控制代碼。當遇到空字串時，Prosite.read 和 Prodoc.read 將會引發 ValueError。您可以捕獲這些例外來偵測無效的登錄編號。

函式 get_prosite_entry() 和 get_prodoc_entry() 用於下載 HTML 格式的 Prosite 和 Prosite 文件記錄。若要建立顯示一則 Prosite 記錄的網頁，您可以使用

>>> from Bio import ExPASy
>>> handle = ExPASy.get_prosite_entry("PS00001")
>>> html = handle.read()
>>> with open("myprositerecord.html", "w") as out_handle:
...     out_handle.write(html)
...

Prosite 文件記錄也是如此

>>> from Bio import ExPASy
>>> handle = ExPASy.get_prodoc_entry("PDOC00001")
>>> html = handle.read()
>>> with open("myprodocrecord.html", "w") as out_handle:
...     out_handle.write(html)
...

對於這些函式，無效的登錄編號會傳回 HTML 格式的錯誤訊息。

掃描 Prosite 資料庫

ScanProsite 可讓您透過提供 UniProt 或 PDB 序列識別碼或序列本身，線上掃描蛋白質序列以對照 Prosite 資料庫。如需有關 ScanProsite 的更多資訊，請參閱 ScanProsite 文件以及以程式化方式存取 ScanProsite 的文件。

您可以使用 Biopython 的 Bio.ExPASy.ScanProsite 模組從 Python 掃描 Prosite 資料庫。此模組可協助您以程式化的方式存取 ScanProsite，並解析 ScanProsite 傳回的結果。若要在以下蛋白質序列中掃描 Prosite 模式

MEHKEVVLLLLLFLKSGQGEPLDDYVNTQGASLFSVTKKQLGAGSIEECAAKCEEDEEFT
CRAFQYHSKEQQCVIMAENRKSSIIIRMRDVVLFEKKVYLSECKTGNGKNYRGTMSKTKN

您可以使用以下程式碼

>>> sequence = (
...     "MEHKEVVLLLLLFLKSGQGEPLDDYVNTQGASLFSVTKKQLGAGSIEECAAKCEEDEEFT"
...     "CRAFQYHSKEQQCVIMAENRKSSIIIRMRDVVLFEKKVYLSECKTGNGKNYRGTMSKTKN"
... )
>>> from Bio.ExPASy import ScanProsite
>>> handle = ScanProsite.scan(seq=sequence)

透過執行 handle.read()，您可以取得原始 XML 格式的搜尋結果。讓我們改用 Bio.ExPASy.ScanProsite.read 將原始 XML 解析為 Python 物件

>>> result = ScanProsite.read(handle)
>>> type(result)
<class 'Bio.ExPASy.ScanProsite.Record'>

Bio.ExPASy.ScanProsite.Record 物件衍生自清單，清單中的每個元素都儲存一個 ScanProsite 命中。此物件也會儲存命中次數，以及 ScanProsite 傳回的搜尋序列數量。此 ScanProsite 搜尋產生了六次命中

>>> result.n_seq
1
>>> result.n_match
1
>>> len(result)
1
>>> result[0]
{'sequence_ac': 'USERSEQ1', 'start': 16, 'stop': 98, 'signature_ac': 'PS50948', 'score': '8.873', 'level': '0'}

其他 ScanProsite 參數可以關鍵字引數的形式傳遞；如需更多資訊，請參閱以程式化方式存取 ScanProsite 的文件。例如，傳遞 lowscore=1 以包含低分比對，讓我們找到一個額外的命中

>>> handle = ScanProsite.scan(seq=sequence, lowscore=1)
>>> result = ScanProsite.read(handle)
>>> result.n_match
2