序列比對
序列比對是兩個或多個已彼此對齊的序列的集合 — 通常會插入間隔，並添加前導或尾隨間隔 — 以使所有序列字串的長度相同。
比對可能延伸覆蓋每個序列的完整長度，或者可能僅限於每個序列的子部分。在 Biopython 中，所有序列比對都由 Alignment 物件表示，如比對物件節中所述。Alignment 物件可以透過解析比對軟體（例如 Clustal 或 BLAT）的輸出（如讀取和寫入比對節中所述）或使用 Biopython 的雙序列比對器（可以將兩個序列彼此對齊，如雙序列比對章中所述）來獲得。
有關較舊的 MultipleSeqAlignment 類別和 Bio.AlignIO 中解析序列比對軟體輸出的解析器（產生 MultipleSeqAlignment 物件），請參閱多序列比對物件章。
比對物件
Alignment 類別定義於 Bio.Align 中。通常，您會透過解析比對程式的輸出（讀取和寫入比對節）或執行 Biopython 的雙序列比對器（雙序列比對章）來取得 Alignment 物件。然而，為了本節的利益，我們將從頭開始建立一個 Alignment 物件。
從序列和座標建立比對物件
假設您有三個序列
>>> seqA = "CCGGTTTTT"
>>> seqB = "AGTTTAA"
>>> seqC = "AGGTTT"
>>> sequences = [seqA, seqB, seqC]
要建立 Alignment 物件，我們還需要定義序列如何彼此對齊的座標。我們使用 NumPy 陣列來做到這一點
>>> import numpy as np
>>> coordinates = np.array([[1, 3, 4, 7, 9], [0, 2, 2, 5, 5], [0, 2, 3, 6, 6]])
這些座標定義了以下序列片段的比對
SeqA[1:3]、SeqB[0:2]和SeqC[0:2]彼此對齊；SeqA[3:4]和SeqC[2:3]彼此對齊，seqB中有一個核苷酸的間隔；SeqA[4:7]、SeqB[2:5]和SeqC[3:6]彼此對齊；SeqA[7:9]沒有與seqB或seqC對齊。
請注意，比對不包括 seqA 的第一個核苷酸和 seqB 的最後兩個核苷酸。
現在我們可以建立 Alignment 物件
>>> from Bio.Align import Alignment
>>> alignment = Alignment(sequences, coordinates)
>>> alignment
<Alignment object (3 rows x 8 columns) at ...>
比對物件有一個屬性 sequences 指向此比對中包含的序列
>>> alignment.sequences
['CCGGTTTTT', 'AGTTTAA', 'AGGTTT']
以及一個具有比對座標的屬性 coordinates
>>> alignment.coordinates
array([[1, 3, 4, 7, 9],
[0, 2, 2, 5, 5],
[0, 2, 3, 6, 6]])
列印 Alignment 物件以明確顯示比對
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
其中每個序列的起始和結束座標分別顯示在比對的左側和右側。
從已比對的序列建立比對物件
如果您從已比對的序列開始，其中破折號代表間隔，則可以使用 parse_printed_alignment 類別方法來計算座標。此方法主要用於 Biopython 的比對解析器中（請參閱讀取和寫入比對節），但它可能適用於其他用途。例如，您可以從已比對的序列建構 Alignment 物件，如下所示
>>> lines = ["CGGTTTTT", "AG-TTT--", "AGGTTT--"]
>>> for line in lines:
... print(line)
...
CGGTTTTT
AG-TTT--
AGGTTT--
>>> lines = [line.encode() for line in lines] # convert to bytes
>>> lines
[b'CGGTTTTT', b'AG-TTT--', b'AGGTTT--']
>>> sequences, coordinates = Alignment.parse_printed_alignment(lines)
>>> sequences
[b'CGGTTTTT', b'AGTTT', b'AGGTTT']
>>> sequences = [sequence.decode() for sequence in sequences]
>>> sequences
['CGGTTTTT', 'AGTTT', 'AGGTTT']
>>> print(coordinates)
[[0 2 3 6 8]
[0 2 2 5 5]
[0 2 3 6 6]]
seqA 的初始 G 核苷酸和 seqB 的最後 CC 核苷酸未包含在比對中，因此這裡遺失了。但是這很容易修復
>>> from Bio.Seq import Seq
>>> sequences[0] = "C" + sequences[0]
>>> sequences[1] = sequences[1] + "AA"
>>> sequences
['CCGGTTTTT', 'AGTTTAA', 'AGGTTT']
>>> coordinates[0, :] += 1
>>> print(coordinates)
[[1 3 4 7 9]
[0 2 2 5 5]
[0 2 3 6 6]]
現在我們可以建立 Alignment 物件
>>> alignment = Alignment(sequences, coordinates)
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
它與從序列和座標建立比對物件節中建立的 Alignment 物件相同。
預設情況下，Alignment 初始化器的 coordinates 引數為 None，這表示比對中沒有間隔。無間隔比對中的所有序列必須具有相同的長度。如果 coordinates 引數為 None，則初始化器會為您填寫 Alignment 物件的 coordinates 屬性
>>> ungapped_alignment = Alignment(["ACGTACGT", "AAGTACGT", "ACGTACCT"])
>>> ungapped_alignment
<Alignment object (3 rows x 8 columns) at ...>
>>> print(ungapped_alignment.coordinates)
[[0 8]
[0 8]
[0 8]]
>>> print(ungapped_alignment)
0 ACGTACGT 8
0 AAGTACGT 8
0 ACGTACCT 8
常見的比對屬性
以下屬性在 Alignment 物件中很常見
sequences：這是彼此對齊的序列列表。根據比對的建立方式，序列可以具有以下類型普通的 Python 字串；
Seq;MutableSeq;SeqRecord;位元組;位元組陣列;資料類型為
numpy.int32的 NumPy 陣列；任何其他具有格式為
"c"、"B"、"i"或"I"的連續緩衝區的物件；在建立比對的
PairwiseAligner物件的alphabet屬性中定義的物件的列表或元組（請參閱廣義雙序列比對節）。
對於雙序列比對（表示兩個序列的比對），屬性
target和query分別是sequences[0]和sequences[1]的別名。coordinates：儲存定義序列如何彼此對齊的序列索引的整數 NumPy 陣列；score：比對分數，如比對檔案中的解析器所發現，或由PairwiseAligner所計算（請參閱基本用法節）；annotations：儲存與比對相關聯的大多數其他註解的字典；column_annotations：儲存沿著比對延伸並與比對長度相同的註解的字典，例如共識序列（有關範例，請參閱ClustalW節）。
由 Bio.Align 中的解析器建立的 Alignment 物件，可能會因為讀取對齊檔案格式的不同而有額外的屬性。
對齊的切片和索引
形式為 alignment[k, i:j] 的切片，其中 k 是整數，而 i 和 j 是整數或缺席，會返回一個字串，顯示目標（如果 k=0）或查詢（如果 k=1）的對齊序列（包括間隙），其中只包含列印的對齊中從 i 到 j 的欄位。
為了說明這一點，在以下示例中，列印的對齊有 8 個欄位
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
>>> alignment.length
8
要單獨取得對齊的序列字串，請使用
>>> alignment[0]
'CGGTTTTT'
>>> alignment[1]
'AG-TTT--'
>>> alignment[2]
'AGGTTT--'
>>> alignment[0, :]
'CGGTTTTT'
>>> alignment[1, :]
'AG-TTT--'
>>> alignment[0, 1:-1]
'GGTTTT'
>>> alignment[1, 1:-1]
'G-TTT-'
也可以使用整數的可迭代物件來選擇要包含的欄位
>>> alignment[0, (1, 2, 4)]
'GGT'
>>> alignment[1, range(0, 5, 2)]
'A-T'
要取得列印的對齊中位置 [i, j] 的字母，請使用 alignment[i, j]；如果該位置有間隙，則會傳回 "-"
>>> alignment[0, 2]
'G'
>>> alignment[2, 6]
'-'
要取得對齊中的特定欄位，請使用
>>> alignment[:, 0]
'CAA'
>>> alignment[:, 1]
'GGG'
>>> alignment[:, 2]
'G-G'
形式為 alignment[i:j:k] 的切片會傳回一個新的 Alignment 物件，其中僅包含對齊的序列 [i:j:k]
>>> alignment[1:]
<Alignment object (2 rows x 6 columns) at ...>
>>> print(alignment[1:])
target 0 AG-TTT 5
0 ||-||| 6
query 0 AGGTTT 6
形式為 alignment[:, i:j] 的切片，其中 i 和 j 是整數或缺席，會傳回一個新的 Alignment 物件，其中僅包含列印的對齊中從 i 到 j 的欄位。
提取以上範例對齊的前 4 個欄位會得到
>>> alignment[:, :4]
<Alignment object (3 rows x 4 columns) at ...>
>>> print(alignment[:, :4])
1 CGGT 5
0 AG-T 3
0 AGGT 4
同樣地，提取最後 6 個欄位會得到
>>> alignment[:, -6:]
<Alignment object (3 rows x 6 columns) at ...>
>>> print(alignment[:, -6:])
3 GTTTTT 9
2 -TTT-- 5
2 GTTT-- 6
欄位索引也可以是整數的可迭代物件
>>> print(alignment[:, (1, 3, 0)])
0 GTC 3
0 GTA 3
0 GTA 3
呼叫 alignment[:, :] 會傳回對齊的副本。
取得有關對齊的資訊
對齊形狀
對齊序列的數量由 len(alignment) 傳回
>>> len(alignment)
3
對齊長度定義為列印的對齊中的欄位數。這等於每個序列中比對、不比對數目和間隙總長度的總和
>>> alignment.length
8
shape 屬性會傳回一個由對齊長度和列印的對齊中的欄位數組成的元組
>>> alignment.shape
(3, 8)
比較對齊
如果 alignment1.sequences 和 alignment2.sequences 中的每個序列都彼此相等，並且 alignment1.coordinates 和 alignment2.coordinates 包含相同的座標，則兩個對齊相等（表示 alignment1 == alignment2 的值為 True）。如果這些條件中的任何一個未滿足，則 alignment1 == alignment2 的值為 False。兩個對齊的不等式（例如，alignment1 < alignment2）的建立方法是先比較 alignment1.sequences 和 alignment2.sequences，如果它們相等，則比較 alignment1.coordinates 和 alignment2.coordinates。
尋找對齊序列的索引
對於成對對齊，對齊的 aligned 屬性會傳回目標和查詢序列中彼此對齊的子序列的起始和結束索引。如果目標 (t) 和查詢 (q) 之間的對齊由 \(N\) 個區塊組成，則會得到兩個長度為 \(N\) 的元組
(((t_start1, t_end1), (t_start2, t_end2), ..., (t_startN, t_endN)),
((q_start1, q_end1), (q_start2, q_end2), ..., (q_startN, q_endN)))
例如，
>>> pairwise_alignment = alignment[:2, :]
>>> print(pairwise_alignment)
target 1 CGGTTTTT 9
0 .|-|||-- 8
query 0 AG-TTT-- 5
>>> print(pairwise_alignment.aligned)
[[[1 3]
[4 7]]
[[0 2]
[2 5]]]
請注意，不同的對齊可能會有相同的子序列彼此對齊。特別是，如果對齊僅在間隙位置方面彼此不同，則可能會發生這種情況
>>> pairwise_alignment1 = Alignment(["AAACAAA", "AAAGAAA"],
... np.array([[0, 3, 4, 4, 7], [0, 3, 3, 4, 7]])) # fmt: skip
...
>>> pairwise_alignment2 = Alignment(["AAACAAA", "AAAGAAA"],
... np.array([[0, 3, 3, 4, 7], [0, 3, 4, 4, 7]])) # fmt: skip
...
>>> print(pairwise_alignment1)
target 0 AAAC-AAA 7
0 |||--||| 8
query 0 AAA-GAAA 7
>>> print(pairwise_alignment2)
target 0 AAA-CAAA 7
0 |||--||| 8
query 0 AAAG-AAA 7
>>> pairwise_alignment1.aligned
array([[[0, 3],
[4, 7]],
[[0, 3],
[4, 7]]])
>>> pairwise_alignment2.aligned
array([[[0, 3],
[4, 7]],
[[0, 3],
[4, 7]]])
indices 屬性會傳回一個 2D NumPy 陣列，其中包含對齊中每個字母的序列索引，間隙以 -1 表示
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
>>> alignment.indices
array([[ 1, 2, 3, 4, 5, 6, 7, 8],
[ 0, 1, -1, 2, 3, 4, -1, -1],
[ 0, 1, 2, 3, 4, 5, -1, -1]])
inverse_indices 屬性會傳回 1D NumPy 陣列的列表，對齊的每個序列都有一個陣列，其中包含序列中每個字母在對齊中的欄位索引。未包含在對齊中的字母以 -1 表示
>>> alignment.sequences
['CCGGTTTTT', 'AGTTTAA', 'AGGTTT']
>>> alignment.inverse_indices
[array([-1, 0, 1, 2, 3, 4, 5, 6, 7]),
array([ 0, 1, 3, 4, 5, -1, -1]),
array([0, 1, 2, 3, 4, 5])]
計算比對、不比對和間隙
counts 方法會計算成對對齊的比對、不比對和間隙的數量。對於超過兩個序列的對齊，會計算並加總對齊中所有成對序列的比對、不比對和間隙的數量。這三個數字會以 AlignmentCounts 物件的形式傳回，此物件是一個具備欄位 gaps、identities 和 mismatches 的 namedtuple。此方法目前不接受引數，但未來可能會修改為接受允許自訂其行為的選用引數。
>>> print(pairwise_alignment)
target 1 CGGTTTTT 9
0 .|-|||-- 8
query 0 AG-TTT-- 5
>>> pairwise_alignment.counts()
AlignmentCounts(gaps=3, identities=4, mismatches=1)
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
>>> alignment.counts()
AlignmentCounts(gaps=8, identities=14, mismatches=2)
字母頻率
frequencies 方法會計算每個字母在對齊的每個欄位中出現的頻率
>>> alignment.frequencies
{'C': array([1., 0., 0., 0., 0., 0., 0., 0.]),
'G': array([0., 3., 2., 0., 0., 0., 0., 0.]),
'T': array([0., 0., 0., 3., 3., 3., 1., 1.]),
'A': array([2., 0., 0., 0., 0., 0., 0., 0.]),
'-': array([0., 0., 1., 0., 0., 0., 2., 2.])}
取代
使用 substitutions 方法尋找每對核苷酸之間的取代次數
>>> m = alignment.substitutions
>>> print(m)
A C G T
A 1.0 0.0 0.0 0.0
C 2.0 0.0 0.0 0.0
G 0.0 0.0 4.0 0.0
T 0.0 0.0 0.0 9.0
請注意，矩陣不是對稱的：列字母 R 和欄字母 C 的計數是序列中字母 R 被下方序列中字母 C 取代的次數。例如，對齊到稍後序列中 A 的 C 的數量為
>>> m["C", "A"]
2.0
而對齊到稍後序列中 C 的 A 的數量為
>>> m["A", "C"]
0.0
要取得對稱矩陣，請使用
>>> m += m.transpose()
>>> m /= 2.0
>>> print(m)
A C G T
A 1.0 1.0 0.0 0.0
C 1.0 0.0 0.0 0.0
G 0.0 0.0 4.0 0.0
T 0.0 0.0 0.0 9.0
>>> m["A", "C"]
1.0
>>> m["C", "A"]
1.0
對齊中 A 和 T 之間的取代總數為 1.0 + 1.0 = 2。
將對齊視為陣列
使用 NumPy，您可以將 alignment 物件變成字母陣列。特別是，這對於快速計算對齊內容可能很有用。
>>> align_array = np.array(alignment)
>>> align_array.shape
(3, 8)
>>> align_array
array([[b'C', b'G', b'G', b'T', b'T', b'T', b'T', b'T'],
[b'A', b'G', b'-', b'T', b'T', b'T', b'-', b'-'],
[b'A', b'G', b'G', b'T', b'T', b'T', b'-', b'-']], dtype='|S1')
依預設，這會為您提供一個 bytes 字元陣列（資料類型為 dtype='|S1'）。您可以使用 dtype='U' 來建立 Unicode（Python 字串）字元的陣列
>>> align_array = np.array(alignment, dtype="U")
>>> align_array
array([['C', 'G', 'G', 'T', 'T', 'T', 'T', 'T'],
['A', 'G', '-', 'T', 'T', 'T', '-', '-'],
['A', 'G', 'G', 'T', 'T', 'T', '-', '-']], dtype='<U1')
（列印的 dtype 將為「<U1」或「>U1」，具體取決於您的系統是小端還是大端）。請注意，alignment 物件和 NumPy 陣列 align_array 是記憶體中獨立的物件 - 編輯其中一個不會更新另一個！
對齊的操作
排序對齊
sort 方法會排序對齊序列。依預設，排序會根據每個序列的 id 屬性（如果有的話）或序列內容進行。
>>> print(alignment)
1 CGGTTTTT 9
0 AG-TTT-- 5
0 AGGTTT-- 6
>>> alignment.sort()
>>> print(alignment)
0 AGGTTT-- 6
0 AG-TTT-- 5
1 CGGTTTTT 9
或者，您可以提供一個 key 函式來決定排序順序。例如，您可以依據 GC 含量遞增來排序序列
>>> from Bio.SeqUtils import gc_fraction
>>> alignment.sort(key=gc_fraction)
>>> print(alignment)
0 AG-TTT-- 5
0 AGGTTT-- 6
1 CGGTTTTT 9
請注意，key 函式會套用至整個序列（包括 seqB 的初始 A 和最後的 GG 核苷酸），而不僅僅是對齊的部分。
reverse 引數可讓您反轉排序順序，以取得 GC 含量遞減的序列
>>> alignment.sort(key=gc_fraction, reverse=True)
>>> print(alignment)
1 CGGTTTTT 9
0 AGGTTT-- 6
0 AG-TTT-- 5
反向互補對齊
反向互補對齊會取得每個序列的反向互補，並重新計算座標
>>> alignment.sequences
['CCGGTTTTT', 'AGGTTT', 'AGTTTAA']
>>> rc_alignment = alignment.reverse_complement()
>>> print(rc_alignment.sequences)
['AAAAACCGG', 'AAACCT', 'TTAAACT']
>>> print(rc_alignment)
0 AAAAACCG 8
0 --AAACCT 6
2 --AAA-CT 7
>>> alignment[:, :4].sequences
['CCGGTTTTT', 'AGGTTT', 'AGTTTAA']
>>> print(alignment[:, :4])
1 CGGT 5
0 AGGT 4
0 AG-T 3
>>> rc_alignment = alignment[:, :4].reverse_complement()
>>> rc_alignment[:, :4].sequences
['AAAAACCGG', 'AAACCT', 'TTAAACT']
>>> print(rc_alignment[:, :4])
4 ACCG 8
2 ACCT 6
4 A-CT 7
反向互補對齊會保留其欄位註解（依相反順序），但會捨棄所有其他註解。
新增對齊
如果對齊具有相同的列數，則可以將對齊加在一起以形成延伸的對齊。例如，讓我們先建立兩個對齊
>>> from Bio.Seq import Seq
>>> from Bio.SeqRecord import SeqRecord
>>> a1 = SeqRecord(Seq("AAAAC"), id="Alpha")
>>> b1 = SeqRecord(Seq("AAAC"), id="Beta")
>>> c1 = SeqRecord(Seq("AAAAG"), id="Gamma")
>>> a2 = SeqRecord(Seq("GTT"), id="Alpha")
>>> b2 = SeqRecord(Seq("TT"), id="Beta")
>>> c2 = SeqRecord(Seq("GT"), id="Gamma")
>>> left = Alignment(
... [a1, b1, c1], coordinates=np.array([[0, 3, 4, 5], [0, 3, 3, 4], [0, 3, 4, 5]])
... )
>>> left.annotations = {"tool": "demo", "name": "start"}
>>> left.column_annotations = {"stats": "CCCXC"}
>>> right = Alignment(
... [a2, b2, c2], coordinates=np.array([[0, 1, 2, 3], [0, 0, 1, 2], [0, 1, 1, 2]])
... )
>>> right.annotations = {"tool": "demo", "name": "end"}
>>> right.column_annotations = {"stats": "CXC"}
現在，讓我們看看這兩個對齊
>>> print(left)
Alpha 0 AAAAC 5
Beta 0 AAA-C 4
Gamma 0 AAAAG 5
>>> print(right)
Alpha 0 GTT 3
Beta 0 -TT 2
Gamma 0 G-T 2
新增這兩個對齊會將這兩個對齊以列的方式結合
>>> combined = left + right
>>> print(combined)
Alpha 0 AAAACGTT 8
Beta 0 AAA-C-TT 6
Gamma 0 AAAAGG-T 7
為了使此操作有效，兩個對齊都必須具有相同數量的序列（這裡它們都有 3 列）
>>> len(left)
3
>>> len(right)
3
>>> len(combined)
3
這些序列是 SeqRecord 物件，可以將它們加在一起。有關如何處理註解的詳細資訊，請參閱第 序列註解物件 章。此範例是一種特殊情況，因為原始對齊都共享相同的名稱，這表示當新增列時，它們也會取得相同的名稱。
任何常見的註解都會被保留，但不同的註解則會遺失。這與 SeqRecord 註解中使用的行為相同，目的是防止不適當的值被意外傳播。
>>> combined.annotations
{'tool': 'demo'}
同樣地，任何常見的每列註解都會被合併。
>>> combined.column_annotations
{'stats': 'CCCXCCXC'}
映射成對序列比對
假設您有一個轉錄本與染色體之間的成對比對。
>>> chromosome = "AAAAAAAACCCCCCCAAAAAAAAAAAGGGGGGAAAAAAAA"
>>> transcript = "CCCCCCCGGGGGG"
>>> sequences1 = [chromosome, transcript]
>>> coordinates1 = np.array([[8, 15, 26, 32], [0, 7, 7, 13]])
>>> alignment1 = Alignment(sequences1, coordinates1)
>>> print(alignment1)
target 8 CCCCCCCAAAAAAAAAAAGGGGGG 32
0 |||||||-----------|||||| 24
query 0 CCCCCCC-----------GGGGGG 13
以及轉錄本與序列（例如，透過 RNA-seq 取得）之間的成對比對。
>>> rnaseq = "CCCCGGGG"
>>> sequences2 = [transcript, rnaseq]
>>> coordinates2 = np.array([[3, 11], [0, 8]])
>>> alignment2 = Alignment(sequences2, coordinates2)
>>> print(alignment2)
target 3 CCCCGGGG 11
0 |||||||| 8
query 0 CCCCGGGG 8
使用 alignment1 上的 map 方法，並以 alignment2 作為引數，以找出 RNA 序列與基因組的比對。
>>> alignment3 = alignment1.map(alignment2)
>>> print(alignment3)
target 11 CCCCAAAAAAAAAAAGGGG 30
0 ||||-----------|||| 19
query 0 CCCC-----------GGGG 8
>>> print(alignment3.coordinates)
[[11 15 26 30]
[ 0 4 4 8]]
>>> format(alignment3, "psl")
'8\t0\t0\t0\t0\t0\t1\t11\t+\tquery\t8\t0\t8\ttarget\t40\t11\t30\t2\t4,4,\t0,4,\t11,26,\n'
為了能夠列印序列，在此範例中，我們使用具有定義序列內容的序列建構了 alignment1 和 alignment2。然而，映射比對不依賴序列內容；只有 alignment1 和 alignment2 的坐標被用來建構 alignment3 的坐標。
map 方法也可以用來提升不同基因組組裝之間的比對。在這種情況下，self 是兩個基因組組裝之間的 DNA 比對，而引數是轉錄本針對其中一個基因組組裝的比對。
>>> from Bio import Align
>>> chain = Align.read("Blat/panTro5ToPanTro6.over.chain", "chain")
>>> chain.sequences[0].id
'chr1'
>>> len(chain.sequences[0].seq)
228573443
>>> chain.sequences[1].id
'chr1'
>>> len(chain.sequences[1].seq)
224244399
>>> import numpy as np
>>> np.set_printoptions(threshold=5) # print 5 array elements per row
>>> print(chain.coordinates)
[[122250000 122250400 122250400 ... 122909818 122909819 122909835]
[111776384 111776784 111776785 ... 112019962 112019962 112019978]]
這顯示了黑猩猩基因組組裝 panTro5 上 chr1 的範圍 122250000:122909835 與黑猩猩基因組組裝 panTro6 上 chr1 的範圍 111776384:112019978 對齊。有關 chain 檔案格式的更多資訊，請參閱 UCSC chain 檔案格式 章節。
>>> transcript = Align.read("Blat/est.panTro5.psl", "psl")
>>> transcript.sequences[0].id
'chr1'
>>> len(transcript.sequences[0].seq)
228573443
>>> transcript.sequences[1].id
'DC525629'
>>> len(transcript.sequences[1].seq)
407
>>> print(transcript.coordinates)
[[122835789 122835847 122840993 122841145 122907212 122907314]
[ 32 90 90 242 242 344]]
這顯示了表達序列標籤 DC525629 的核苷酸範圍 32:344 與黑猩猩基因組組裝 panTro5 上 chr1 的範圍 122835789:122907314 對齊。請注意，目標序列 chain.sequences[0].seq 和目標序列 transcript.sequences[0] 具有相同的長度。
>>> len(chain.sequences[0].seq) == len(transcript.sequences[0].seq)
True
我們交換 chain 的目標和查詢，使得 chain 的查詢對應於 transcript 的目標。
>>> chain = chain[::-1]
>>> chain.sequences[0].id
'chr1'
>>> len(chain.sequences[0].seq)
224244399
>>> chain.sequences[1].id
'chr1'
>>> len(chain.sequences[1].seq)
228573443
>>> print(chain.coordinates)
[[111776384 111776784 111776785 ... 112019962 112019962 112019978]
[122250000 122250400 122250400 ... 122909818 122909819 122909835]]
>>> np.set_printoptions(threshold=1000) # reset the print options
現在，我們可以透過呼叫 chain.map，並以 transcript 作為引數，取得 DC525629 針對黑猩猩基因組組裝 panTro6 的坐標。
>>> lifted_transcript = chain.map(transcript)
>>> lifted_transcript.sequences[0].id
'chr1'
>>> len(lifted_transcript.sequences[0].seq)
224244399
>>> lifted_transcript.sequences[1].id
'DC525629'
>>> len(lifted_transcript.sequences[1].seq)
407
>>> print(lifted_transcript.coordinates)
[[111982717 111982775 111987921 111988073 112009200 112009302]
[ 32 90 90 242 242 344]]
這顯示了表達序列標籤 DC525629 的核苷酸範圍 32:344 與黑猩猩基因組組裝 panTro6 上 chr1 的範圍 111982717:112009302 對齊。請注意，DC525629 在黑猩猩基因組組裝 panTro5 上的基因組跨度為 122907314 - 122835789 = 71525 bp，而在 panTro6 上，基因組跨度為 112009302 - 111982717 = 26585 bp。
映射多重序列比對
考慮黑猩猩、人類、獼猴、狨猴、小鼠和大鼠的基因組序列的多重比對。
>>> from Bio import Align
>>> path = "Blat/panTro5.maf"
>>> genome_alignment = Align.read(path, "maf")
>>> for record in genome_alignment.sequences:
... print(record.id, len(record.seq))
...
panTro5.chr1 228573443
hg19.chr1 249250621
rheMac8.chr1 225584828
calJac3.chr18 47448759
mm10.chr3 160039680
rn6.chr2 266435125
>>> print(genome_alignment.coordinates)
[[133922962 133922962 133922970 133922970 133922972 133922972 133922995
133922998 133923010]
[155784573 155784573 155784581 155784581 155784583 155784583 155784606
155784609 155784621]
[130383910 130383910 130383918 130383918 130383920 130383920 130383943
130383946 130383958]
[ 9790455 9790455 9790463 9790463 9790465 9790465 9790488
9790491 9790503]
[ 88858039 88858036 88858028 88858026 88858024 88858020 88857997
88857997 88857985]
[188162970 188162967 188162959 188162959 188162957 188162953 188162930
188162930 188162918]]
>>> print(genome_alignment)
panTro5.c 133922962 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
hg19.chr1 155784573 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
rheMac8.c 130383910 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
calJac3.c 9790455 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAGCTTAAAggct
mm10.chr3 88858039 TATAATAATTGTATATGTCACAGAAAAAAATGAATTTTCAAT---GACTTAATAGCC
rn6.chr2 188162970 TACAATAATTG--TATGTCATAGAAAAAAATGAATTTTCAAT---AACTTAATAGCC
panTro5.c 133923010
hg19.chr1 155784621
rheMac8.c 130383958
calJac3.c 9790503
mm10.chr3 88857985
rn6.chr2 188162918
假設我們想將舊版本的基因組組裝（panTro5、hg19、rheMac8、calJac3、mm10 和 rn6）替換為其當前版本（panTro6、hg38、rheMac10、calJac4、mm39 和 rn7）。為此，我們需要每個物種的新舊組裝版本之間的成對比對。這些由 UCSC 以 chain 檔案的形式提供，通常用於 UCSC 的 liftOver 工具。Biopython 原始碼發行版中 Tests/Align 子目錄中的 .chain 檔案是從 UCSC 的 .chain 檔案中提取出來的，僅包含相關的基因組區域。例如，要將 panTro5 提升到 panTro6，我們使用檔案 panTro5ToPanTro6.chain，其內容如下：
chain 1198066 chr1 228573443 + 133919957 133932620 chr1 224244399 + 130607995 130620657 1
4990 0 2
1362 3 0
6308
為了提升每個物種的基因組組裝，我們讀入對應的 .chain 檔案。
>>> paths = [
... "Blat/panTro5ToPanTro6.chain",
... "Blat/hg19ToHg38.chain",
... "Blat/rheMac8ToRheMac10.chain",
... "Blat/calJac3ToCalJac4.chain",
... "Blat/mm10ToMm39.chain",
... "Blat/rn6ToRn7.chain",
... ]
>>> liftover_alignments = [Align.read(path, "chain") for path in paths]
>>> for liftover_alignment in liftover_alignments:
... print(liftover_alignment.target.id, liftover_alignment.coordinates[0, :])
...
chr1 [133919957 133924947 133924947 133926309 133926312 133932620]
chr1 [155184381 156354347 156354348 157128497 157128497 157137496]
chr1 [130382477 130383872 130383872 130384222 130384222 130388520]
chr18 [9786631 9787941 9788508 9788508 9795062 9795065 9795737]
chr3 [66807541 74196805 74196831 94707528 94707528 94708176 94708178 94708718]
chr2 [188111581 188158351 188158351 188171225 188171225 188228261 188228261
188236997]
請注意，物種的順序在 liftover_alignments 和 genome_alignment.sequences 中相同。現在我們可以將多重序列比對提升到新的基因組組裝版本。
>>> genome_alignment = genome_alignment.mapall(liftover_alignments)
>>> for record in genome_alignment.sequences:
... print(record.id, len(record.seq))
...
chr1 224244399
chr1 248956422
chr1 223616942
chr18 47031477
chr3 159745316
chr2 249053267
>>> print(genome_alignment.coordinates)
[[130611000 130611000 130611008 130611008 130611010 130611010 130611033
130611036 130611048]
[155814782 155814782 155814790 155814790 155814792 155814792 155814815
155814818 155814830]
[ 95186253 95186253 95186245 95186245 95186243 95186243 95186220
95186217 95186205]
[ 9758318 9758318 9758326 9758326 9758328 9758328 9758351
9758354 9758366]
[ 88765346 88765343 88765335 88765333 88765331 88765327 88765304
88765304 88765292]
[174256702 174256699 174256691 174256691 174256689 174256685 174256662
174256662 174256650]]
由於 .chain 檔案不包含序列內容，因此我們無法直接列印序列比對。相反地，我們分別讀入每個物種的基因組序列（作為 .2bit 檔案，因為它允許延遲載入；請參閱 序列檔案作為字典 章節）。
>>> from Bio import SeqIO
>>> names = ("panTro6", "hg38", "rheMac10", "calJac4", "mm39", "rn7")
>>> for i, name in enumerate(names):
... filename = f"{name}.2bit"
... genome = SeqIO.parse(filename, "twobit")
... chromosome = genome_alignment.sequences[i].id
... assert len(genome_alignment.sequences[i]) == len(genome[chromosome])
... genome_alignment.sequences[i] = genome[chromosome]
... genome_alignment.sequences[i].id = f"{name}.{chromosome}"
...
>>> print(genome_alignment)
panTro6.c 130611000 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
hg38.chr1 155814782 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
rheMac10. 95186253 ---ACTAGTTA--CA----GTAACAGAAAATAAAATTTAAATAGAAACTTAAAggcc
calJac4.c 9758318 ---ACTAGTTA--CA----GTAACAGAaaataaaatttaaatagaagcttaaaggct
mm39.chr3 88765346 TATAATAATTGTATATGTCACAGAAAAAAATGAATTTTCAAT---GACTTAATAGCC
rn7.chr2 174256702 TACAATAATTG--TATGTCATAGAAAAAAATGAATTTTCAAT---AACTTAATAGCC
panTro6.c 130611048
hg38.chr1 155814830
rheMac10. 95186205
calJac4.c 9758366
mm39.chr3 88765292
rn7.chr2 174256650
mapall 方法也可以用來從相應氨基酸序列的多重序列比對中建立密碼子序列的多重比對（有關詳細資訊，請參閱 產生密碼子序列的多重序列比對 章節）。
Alignments 類別
Alignments（複數）類別繼承自 AlignmentsAbstractBaseClass 和 list，並且可以作為列表來儲存 Alignment 物件。Alignments 物件的行為與 list 物件在兩個重要方面不同。
Alignments物件是它自己的迭代器，與Bio.Align.parse（請參閱 讀取比對 章節）或成對比對器（請參閱 成對序列比對 章節）所傳回的迭代器一致。在迭代器上呼叫iter將始終傳回Alignments物件本身。相反地，在 list 物件上呼叫iter每次都會建立一個新的迭代器，讓您對給定的列表有多個獨立的迭代器。在此範例中，
alignment_iterator1和alignment_iterator2是從列表取得的，並且彼此獨立運作。>>> alignment_list = [alignment1, alignment2, alignment3] >>> alignment_iterator1 = iter(alignment_list) >>> alignment_iterator2 = iter(alignment_list) >>> next(alignment_iterator1) <Alignment object (2 rows x 24 columns) at ...> >>> next(alignment_iterator2) <Alignment object (2 rows x 24 columns) at ...> >>> next(alignment_iterator1) <Alignment object (2 rows x 8 columns) at ...> >>> next(alignment_iterator1) <Alignment object (2 rows x 19 columns) at ...> >>> next(alignment_iterator2) <Alignment object (2 rows x 8 columns) at ...> >>> next(alignment_iterator2) <Alignment object (2 rows x 19 columns) at ...>
相反地，透過在
Alignments物件上呼叫iter取得的alignment_iterator1和alignment_iterator2彼此相同。>>> from Bio.Align import Alignments >>> alignments = Alignments([alignment1, alignment2, alignment3]) >>> alignment_iterator1 = iter(alignments) >>> alignment_iterator2 = iter(alignments) >>> alignment_iterator1 is alignment_iterator2 True >>> next(alignment_iterator1) <Alignment object (2 rows x 24 columns) at ...> >>> next(alignment_iterator2) <Alignment object (2 rows x 8 columns) at ...> >>> next(alignment_iterator1) <Alignment object (2 rows x 19 columns) at ...> >>> next(alignment_iterator2) Traceback (most recent call last): File "<stdin>", line 1, in <module> StopIteration
在
Alignments物件上呼叫iter會將迭代器重設為其第一個項目，因此您可以再次迴圈。您也可以使用for迴圈多次迭代比對，這會隱式呼叫迭代器上的iter。>>> for item in alignments: ... print(repr(item)) ... <Alignment object (2 rows x 24 columns) at ...> <Alignment object (2 rows x 8 columns) at ...> <Alignment object (2 rows x 19 columns) at ...> >>> for item in alignments: ... print(repr(item)) ... <Alignment object (2 rows x 24 columns) at ...> <Alignment object (2 rows x 8 columns) at ...> <Alignment object (2 rows x 19 columns) at ...>
此行為與常規 Python 列表以及
Bio.Align.parse（請參閱 讀取比對 章節）或成對比對器（請參閱 成對序列比對 章節）所傳回的迭代器一致。中繼資料可以作為屬性儲存在
Alignments物件上，而普通的list不接受屬性。>>> alignment_list.score = 100 Traceback (most recent call last): ... AttributeError: 'list' object has no attribute 'score'... >>> alignments.score = 100 >>> alignments.score 100
讀取和寫入比對
序列比對軟體（如 Clustal）的輸出可以透過 Bio.Align.read 和 Bio.Align.parse 函式解析為 Alignment 物件。它們的用法與 Bio.SeqIO 中的 read 和 parse 函式類似（請參閱 解析或讀取序列 章節）：read 函式用於讀取包含單個比對的輸出檔案，並傳回 Alignment 物件，而 parse 函式傳回一個迭代器，以迭代儲存在包含一個或多個比對的輸出檔案中的比對。 比對檔案格式 章節說明可以在 Bio.Align 中解析的比對格式。Bio.Align 還提供了一個 write 函式，可以用這些格式中的大多數寫入比對。
讀取比對
使用 Bio.Align.parse 來解析序列比對的檔案。例如，檔案 ucsc_mm9_chr10.maf 包含 48 個 MAF（多重比對格式）格式的多重序列比對（請參閱 多重比對格式 (MAF) 章節）。
>>> from Bio import Align
>>> alignments = Align.parse("MAF/ucsc_mm9_chr10.maf", "maf")
>>> alignments
<Bio.Align.maf.AlignmentIterator object at 0x...>
其中 "maf" 是檔案格式。Bio.Align.parse 傳回的比對物件可能包含儲存在檔案中的中繼資料的屬性，例如用於建立比對的軟體版本號碼。每個檔案格式儲存的特定屬性在 比對檔案格式 章節中說明。對於 MAF 檔案，我們可以取得檔案格式版本和使用的評分方案。
>>> alignments.metadata
{'MAF Version': '1', 'Scoring': 'autoMZ.v1'}
由於比對檔案可能非常大，Align.parse 會回傳一個比對的迭代器，因此您不必同時將所有比對儲存在記憶體中。您可以迭代這些比對，並印出例如每個比對中對齊序列的數量。
>>> for a in alignments:
... print(len(a.sequences))
...
2
4
5
6
...
15
14
7
6
您也可以在比對上呼叫 len 來取得比對的數量。
>>> len(alignments)
48
根據檔案格式，比對的數量可能會明確儲存在檔案中（例如 bigBed、bigPsl 和 bigMaf 檔案的情況），否則比對的數量會透過迴圈遍歷一次來計算（並將迭代器返回到原始位置）。如果檔案很大，len 可能需要相當長的時間才能返回。但是，由於比對的數量會被快取，後續對 len 的呼叫將會快速返回。
如果比對的數量不是非常大且可以放入記憶體，您可以將比對迭代器轉換為比對列表。要做到這一點，您可以對 alignments 呼叫 list。
>>> alignment_list = list(alignments)
>>> len(alignment_list)
48
>>> alignment_list[27]
<Alignment object (3 rows x 91 columns) at 0x...>
>>> print(alignment_list[27])
mm9.chr10 3019377 CCCCAGCATTCTGGCAGACACAGTG-AAAAGAGACAGATGGTCACTAATAAAATCTGT-A
felCat3.s 46845 CCCAAGTGTTCTGATAGCTAATGTGAAAAAGAAGCATGTGCCCACCAGTAAGCTTTGTGG
canFam2.c 47545247 CCCAAGTGTTCTGATTGCCTCTGTGAAAAAGAAACATGGGCCCGCTAATAagatttgcaa
mm9.chr10 3019435 TAAATTAG-ATCTCAGAGGATGGATGGACCA 3019465
felCat3.s 46785 TGAACTAGAATCTCAGAGGATG---GGACTC 46757
canFam2.c 47545187 tgacctagaatctcagaggatg---ggactc 47545159
但這樣會遺失元數據資訊。
>>> alignment_list.metadata
Traceback (most recent call last):
...
AttributeError: 'list' object has no attribute 'metadata'
相反地，您可以要求比對的完整切片。
>>> type(alignments)
<class 'Bio.Align.maf.AlignmentIterator'>
>>> alignments = alignments[:]
>>> type(alignments)
<class 'Bio.Align.Alignments'>
這會回傳一個 Bio.Align.Alignments 物件，該物件可以像列表一樣使用，同時保留元數據資訊。
>>> len(alignments)
48
>>> print(alignments[11])
mm9.chr10 3014742 AAGTTCCCTCCATAATTCCTTCCTCCCACCCCCACA 3014778
calJac1.C 6283 AAATGTA-----TGATCTCCCCATCCTGCCCTG--- 6311
otoGar1.s 175262 AGATTTC-----TGATGCCCTCACCCCCTCCGTGCA 175231
loxAfr1.s 9317 AGGCTTA-----TG----CCACCCCCCACCCCCACA 9290
>>> alignments.metadata
{'MAF Version': '1', 'Scoring': 'autoMZ.v1'}
寫入比對
要將比對寫入檔案，請使用
>>> from Bio import Align
>>> target = "myfile.txt"
>>> Align.write(alignments, target, "clustal")
其中 alignments 是單個比對或比對列表，target 是檔案名稱或類似檔案的開啟物件，而 "clustal" 是要使用的檔案格式。由於某些檔案格式允許或要求將元數據與比對一起儲存，因此您可能需要使用 Alignments（複數）類別，而不是簡單的比對列表（請參閱第 比對類別 節），允許您將元數據字典作為 alignments 物件的屬性儲存。
>>> from Bio import Align
>>> alignments = Align.Alignments(alignments)
>>> metadata = {"Program": "Biopython", "Version": "1.81"}
>>> alignments.metadata = metadata
>>> target = "myfile.txt"
>>> Align.write(alignments, target, "clustal")
列印比對
對於文字（非二進制）格式，您可以對比對呼叫 Python 內建的 format 函數，以取得一個字串，顯示所要求格式的比對，或在格式化 (f-) 字串中使用 Alignment 物件。如果在沒有參數的情況下呼叫，format 函數會回傳比對的字串表示形式。
>>> str(alignment)
' 1 CGGTTTTT 9\n 0 AGGTTT-- 6\n 0 AG-TTT-- 5\n'
>>> format(alignment)
' 1 CGGTTTTT 9\n 0 AGGTTT-- 6\n 0 AG-TTT-- 5\n'
>>> print(format(alignment))
1 CGGTTTTT 9
0 AGGTTT-- 6
0 AG-TTT-- 5
透過指定第 比對檔案格式 節中顯示的其中一種格式，format 將建立一個字串，顯示所要求格式的比對。
>>> format(alignment, "clustal")
'sequence_0 CGGTTTTT\nsequence_1 AGGTTT--\nsequence_2 AG-TTT--\n\n\n'
>>> print(format(alignment, "clustal"))
sequence_0 CGGTTTTT
sequence_1 AGGTTT--
sequence_2 AG-TTT--
>>> print(f"*** this is the alignment in Clustal format: ***\n{alignment:clustal}\n***")
*** this is the alignment in Clustal format: ***
sequence_0 CGGTTTTT
sequence_1 AGGTTT--
sequence_2 AG-TTT--
***
>>> format(alignment, "maf")
'a\ns sequence_0 1 8 + 9 CGGTTTTT\ns sequence_1 0 6 + 6 AGGTTT--\ns sequence_2 0 5 + 7 AG-TTT--\n\n'
>>> print(format(alignment, "maf"))
a
s sequence_0 1 8 + 9 CGGTTTTT
s sequence_1 0 6 + 6 AGGTTT--
s sequence_2 0 5 + 7 AG-TTT--
由於可選的關鍵字參數不能與 Python 內建的 format 函數或格式化字串一起使用，因此 Alignment 類別具有一個 format 方法，該方法帶有可選參數來自訂比對格式，如下面的子節所述。例如，我們可以使用特定數量的列，以 BED 格式列印比對（請參閱第 瀏覽器可延伸資料 (BED) 節）。
>>> print(pairwise_alignment)
target 1 CGGTTTTT 9
0 .|-|||-- 8
query 0 AG-TTT-- 5
>>> print(format(pairwise_alignment, "bed"))
target 1 7 query 0 + 1 7 0 2 2,3, 0,3,
>>> print(pairwise_alignment.format("bed"))
target 1 7 query 0 + 1 7 0 2 2,3, 0,3,
>>> print(pairwise_alignment.format("bed", bedN=3))
target 1 7
>>> print(pairwise_alignment.format("bed", bedN=6))
target 1 7 query 0 +
比對檔案格式
下表顯示可在 Bio.Align 中解析的比對格式。在 Bio.Align 函數中用於指定檔案格式的格式參數 fmt 不區分大小寫。如表中所示，Bio.Align 也可以寫入大多數這些檔案格式。
檔案格式 |
描述 |
文字 / 二進制 |
|
子節 |
|
A2M |
文字 |
是 |
|
|
瀏覽器可擴展資料 (BED) |
文字 |
是 |
|
|
bigBed |
二進制 |
是 |
|
|
bigMaf |
二進制 |
是 |
|
|
bigPsl |
二進制 |
是 |
|
|
UCSC chain 檔案 |
文字 |
是 |
|
|
ClustalW |
文字 |
是 |
|
|
EMBOSS |
文字 |
否 |
|
``exonerate`` |
Exonerate |
文字 |
是 |
|
|
已比對的 FASTA |
文字 |
是 |
|
|
HH-suite 輸出檔案 |
文字 |
否 |
|
|
多重比對格式 (MAF) |
文字 |
是 |
|
|
Mauve eXtended Multi-FastA (xmfa) 格式 |
文字 |
是 |
|
|
GCG 多序列格式 (MSF) |
文字 |
否 |
|
|
NEXUS |
文字 |
是 |
|
|
PHYLIP 輸出檔案 |
文字 |
是 |
|
|
模式空間佈局 (PSL) |
文字 |
是 |
|
|
序列比對/映射 (SAM) |
文字 |
是 |
|
``stockholm`` |
Stockholm |
文字 |
是 |
|
|
來自 BLAST 或 FASTA 的表格輸出 |
文字 |
否 |
對齊的 FASTA
對齊的 FASTA 格式的檔案正好儲存一個（成對或多個）序列比對，其中比對中的間隙以破折號 (-) 表示。使用 fmt="fasta" 以讀取或寫入對齊的 FASTA 格式的檔案。請注意，這與 William Pearson 的 FASTA 比對程式產生的輸出不同（解析此類輸出請參閱第 來自 BLAST 或 FASTA 的表格輸出 節）。
Biopython 測試套件中的檔案 probcons.fa 以對齊的 FASTA 格式儲存一個多重比對。此檔案的內容如下
>plas_horvu
D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-VD-VSKISQEEYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVTV
>plas_chlre
--VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-VN-ADAISRDDYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKIIV
>plas_anava
--VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKSLSHKQLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKITV
>plas_proho
VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-ES-APALSNTKLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTITV
>azup_achcy
VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-AE-A-------FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVEV
要讀取此檔案，請使用
>>> from Bio import Align
>>> alignment = Align.read("probcons.fa", "fasta")
>>> alignment
<Alignment object (5 rows x 101 columns) at ...>
我們可以列印比對以查看其預設表示形式
>>> print(alignment)
plas_horv 0 D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-VD-VSKISQE
plas_chlr 0 --VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-VN-ADAISRD
plas_anav 0 --VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKSLSHK
plas_proh 0 VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-ES-APALSNT
azup_achc 0 VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-AE-A------
plas_horv 57 EYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVTV 95
plas_chlr 56 DYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKIIV 94
plas_anav 58 QLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKITV 99
plas_proh 56 KLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTITV 94
azup_achc 51 -FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVEV 88
或者我們可以以對齊的 FASTA 格式列印它
>>> print(format(alignment, "fasta"))
>plas_horvu
D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-VD-VSKISQEEYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVTV
>plas_chlre
--VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-VN-ADAISRDDYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKIIV
>plas_anava
--VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKSLSHKQLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKITV
>plas_proho
VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-ES-APALSNTKLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTITV
>azup_achcy
VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-AE-A-------FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVEV
或任何其他可用的格式，例如 Clustal（請參閱第 ClustalW 節）
>>> print(format(alignment, "clustal"))
plas_horvu D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-
plas_chlre --VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-
plas_anava --VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKS
plas_proho VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-
azup_achcy VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-
plas_horvu VD-VSKISQEEYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVT
plas_chlre VN-ADAISRDDYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKII
plas_anava ADLAKSLSHKQLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKIT
plas_proho ES-APALSNTKLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTIT
azup_achcy AE-A-------FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVE
plas_horvu V
plas_chlre V
plas_anava V
plas_proho V
azup_achcy V
與比對關聯的序列是 SeqRecord 物件
>>> alignment.sequences
[SeqRecord(seq=Seq('DVLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSGVDVSKI...VTV'), id='plas_horvu', name='<unknown name>', description='', dbxrefs=[]), SeqRecord(seq=Seq('VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSGVNADAIS...IIV'), id='plas_chlre', name='<unknown name>', description='', dbxrefs=[]), SeqRecord(seq=Seq('VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKS...ITV'), id='plas_anava', name='<unknown name>', description='', dbxrefs=[]), SeqRecord(seq=Seq('VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDKVPAGESAPALS...ITV'), id='plas_proho', name='<unknown name>', description='', dbxrefs=[]), SeqRecord(seq=Seq('VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDKGHNVETIKGMIPDGAEAFKS...VEV'), id='azup_achcy', name='<unknown name>', description='', dbxrefs=[])]
請注意，這些序列不包含間隙（「-」字元），因為比對資訊儲存在 coordinates 屬性中。
>>> print(alignment.coordinates)
[[ 0 1 1 33 34 42 44 48 48 50 50 51 58 73 73 95]
[ 0 0 0 32 33 41 43 47 47 49 49 50 57 72 72 94]
[ 0 0 0 32 33 41 43 47 48 50 51 52 59 74 77 99]
[ 0 1 2 34 35 43 43 47 47 49 49 50 57 72 72 94]
[ 0 1 2 34 34 42 44 48 48 50 50 51 51 66 66 88]]
使用 Align.write 將此比對寫入檔案（在此，我們將使用 StringIO 物件而不是檔案）
>>> from io import StringIO
>>> stream = StringIO()
>>> Align.write(alignment, stream, "FASTA")
1
>>> print(stream.getvalue())
>plas_horvu
D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-VD-VSKISQEEYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVTV
>plas_chlre
--VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-VN-ADAISRDDYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKIIV
>plas_anava
--VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKSLSHKQLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKITV
>plas_proho
VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-ES-APALSNTKLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTITV
>azup_achcy
VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-AE-A-------FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVEV
請注意，Align.write 回傳寫入的比對數量（在此例中為 1）。
ClustalW
Clustal 是一組多重序列比對程式，既可以作為獨立程式也可以作為網路伺服器使用。檔案 opuntia.aln（可線上取得或在 Biopython 原始碼的 Doc/examples 子目錄中）是 Clustal 產生的輸出檔案。它的前幾行是
CLUSTAL 2.1 multiple sequence alignment
gi|6273285|gb|AF191659.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273284|gb|AF191658.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273287|gb|AF191661.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273286|gb|AF191660.1|AF191 TATACATAAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273290|gb|AF191664.1|AF191 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273289|gb|AF191663.1|AF191 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273291|gb|AF191665.1|AF191 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
******* **** *************************************
...
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("opuntia.aln", "clustal")
alignments 上的 metadata 屬性儲存檔案標頭中顯示的資訊
>>> alignments.metadata
{'Program': 'CLUSTAL', 'Version': '2.1'}
您可以在 alignments 上呼叫 next 以提取第一個（也是唯一的）比對
>>> alignment = next(alignments)
>>> print(alignment)
gi|627328 0 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627328 0 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627328 0 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627328 0 TATACATAAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627329 0 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627328 0 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627329 0 TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAAAAAAATGAAT
gi|627328 60 CTAAATGATATACGATTCCACTATGTAAGGTCTTTGAATCATATCATAAAAGACAATGTA
gi|627328 60 CTAAATGATATACGATTCCACTATGTAAGGTCTTTGAATCATATCATAAAAGACAATGTA
gi|627328 60 CTAAATGATATACGATTCCACTATGTAAGGTCTTTGAATCATATCATAAAAGACAATGTA
gi|627328 60 CTAAATGATATACGATTCCACTA...
如果您對元數據不感興趣，那麼使用 Align.read 函數會更方便，因為無論如何，每個 Clustal 檔案都只包含一個比對。
>>> from Bio import Align
>>> alignment = Align.read("opuntia.aln", "clustal")
Clustal 輸出檔案中每個比對區塊下方的共識行，如果序列在每個位置上是保守的，則包含一個星號。此資訊儲存在 alignment 的 column_annotations 屬性中
>>> alignment.column_annotations
{'clustal_consensus': '******* **** **********************************...
以 clustal 格式列印 alignment 將顯示序列比對，但不包括元數據
>>> print(format(alignment, "clustal"))
gi|6273285|gb|AF191659.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273284|gb|AF191658.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273287|gb|AF191661.1|AF191 TATACATT...
以 clustal 格式寫入 alignments 將包括元數據和序列比對
>>> from io import StringIO
>>> stream = StringIO()
>>> Align.write(alignments, stream, "clustal")
1
>>> print(stream.getvalue())
CLUSTAL 2.1 multiple sequence alignment
gi|6273285|gb|AF191659.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273284|gb|AF191658.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273287|gb|AF191661.1|AF191 TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAGAAAGAA
gi|6273286|gb|AF191660.1|AF191 TATACATAAAAGAAG...
如果您手動建立比對，並且想要在輸出中包含元數據資訊，請使用 Alignments（複數）物件（請參閱第 比對類別 節）。
Stockholm
這是 PFAM 使用的 Stockholm 檔案格式中蛋白質序列比對的範例
# STOCKHOLM 1.0
#=GF ID 7kD_DNA_binding
#=GF AC PF02294.20
#=GF DE 7kD DNA-binding domain
#=GF AU Mian N;0000-0003-4284-4749
#=GF AU Bateman A;0000-0002-6982-4660
#=GF SE Pfam-B_8148 (release 5.2)
#=GF GA 25.00 25.00;
#=GF TC 26.60 46.20;
#=GF NC 23.20 19.20;
#=GF BM hmmbuild HMM.ann SEED.ann
#=GF SM hmmsearch -Z 57096847 -E 1000 --cpu 4 HMM pfamseq
#=GF TP Domain
#=GF CL CL0049
#=GF RN [1]
#=GF RM 3130377
#=GF RT Microsequence analysis of DNA-binding proteins 7a, 7b, and 7e
#=GF RT from the archaebacterium Sulfolobus acidocaldarius.
#=GF RA Choli T, Wittmann-Liebold B, Reinhardt R;
#=GF RL J Biol Chem 1988;263:7087-7093.
#=GF DR INTERPRO; IPR003212;
#=GF DR SCOP; 1sso; fa;
#=GF DR SO; 0000417; polypeptide_domain;
#=GF CC This family contains members of the hyper-thermophilic
#=GF CC archaebacterium 7kD DNA-binding/endoribonuclease P2 family.
#=GF CC There are five 7kD DNA-binding proteins, 7a-7e, found as
#=GF CC monomers in the cell. Protein 7e shows the tightest DNA-binding
#=GF CC ability.
#=GF SQ 3
#=GS DN7_METS5/4-61 AC A4YEA2.1
#=GS DN7A_SACS2/3-61 AC P61991.2
#=GS DN7A_SACS2/3-61 DR PDB; 1SSO A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 1JIC A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 2CVR A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 1B4O A; 2-60;
#=GS DN7E_SULAC/3-60 AC P13125.2
DN7_METS5/4-61 KIKFKYKGQDLEVDISKVKKVWKVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLNMIGK
DN7A_SACS2/3-61 TVKFKYKGEEKQVDISKIKKVWRVGKMISFTYDEGGGKTGRGAVSEKDAPKELLQMLEK
#=GR DN7A_SACS2/3-61 SS EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT
DN7E_SULAC/3-60 KVRFKYKGEEKEVDTSKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELMDMLAR
#=GC SS_cons EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT
#=GC seq_cons KVKFKYKGEEKEVDISKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLsMLuK
//
這是 7kD_DNA_binding (PF02294.20) PFAM 條目的種子比對，從 InterPro 網站 (https://www.ebi.ac.uk/interpro/) 下載。此版本的 PFAM 條目也可在 Biopython 原始碼發行版中找到，作為子目錄 Tests/Stockholm/ 中的檔案 pfam2.seed.txt。我們可以如下載入此檔案
>>> from Bio import Align
>>> alignment = Align.read("pfam2.seed.txt", "stockholm")
>>> alignment
<Alignment object (3 rows x 59 columns) at ...>
我們可以列印比對的摘要
>>> print(alignment)
DN7_METS5 0 KIKFKYKGQDLEVDISKVKKVWKVGKMVSFTYDD-NGKTGRGAVSEKDAPKELLNMIGK
DN7A_SACS 0 TVKFKYKGEEKQVDISKIKKVWRVGKMISFTYDEGGGKTGRGAVSEKDAPKELLQMLEK
DN7E_SULA 0 KVRFKYKGEEKEVDTSKIKKVWRVGKMVSFTYDD-NGKTGRGAVSEKDAPKELMDMLAR
DN7_METS5 58
DN7A_SACS 59
DN7E_SULA 58
您也可以對比對物件呼叫 Python 內建的 format 函數，以特定檔案格式顯示它（有關詳細資訊，請參閱第 列印比對 節），例如以 Stockholm 格式重新產生檔案
>>> print(format(alignment, "stockholm"))
# STOCKHOLM 1.0
#=GF ID 7kD_DNA_binding
#=GF AC PF02294.20
#=GF DE 7kD DNA-binding domain
#=GF AU Mian N;0000-0003-4284-4749
#=GF AU Bateman A;0000-0002-6982-4660
#=GF SE Pfam-B_8148 (release 5.2)
#=GF GA 25.00 25.00;
#=GF TC 26.60 46.20;
#=GF NC 23.20 19.20;
#=GF BM hmmbuild HMM.ann SEED.ann
#=GF SM hmmsearch -Z 57096847 -E 1000 --cpu 4 HMM pfamseq
#=GF TP Domain
#=GF CL CL0049
#=GF RN [1]
#=GF RM 3130377
#=GF RT Microsequence analysis of DNA-binding proteins 7a, 7b, and 7e from
#=GF RT the archaebacterium Sulfolobus acidocaldarius.
#=GF RA Choli T, Wittmann-Liebold B, Reinhardt R;
#=GF RL J Biol Chem 1988;263:7087-7093.
#=GF DR INTERPRO; IPR003212;
#=GF DR SCOP; 1sso; fa;
#=GF DR SO; 0000417; polypeptide_domain;
#=GF CC This family contains members of the hyper-thermophilic
#=GF CC archaebacterium 7kD DNA-binding/endoribonuclease P2 family. There
#=GF CC are five 7kD DNA-binding proteins, 7a-7e, found as monomers in the
#=GF CC cell. Protein 7e shows the tightest DNA-binding ability.
#=GF SQ 3
#=GS DN7_METS5/4-61 AC A4YEA2.1
#=GS DN7A_SACS2/3-61 AC P61991.2
#=GS DN7A_SACS2/3-61 DR PDB; 1SSO A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 1JIC A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 2CVR A; 2-60;
#=GS DN7A_SACS2/3-61 DR PDB; 1B4O A; 2-60;
#=GS DN7E_SULAC/3-60 AC P13125.2
DN7_METS5/4-61 KIKFKYKGQDLEVDISKVKKVWKVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLNMIGK
DN7A_SACS2/3-61 TVKFKYKGEEKQVDISKIKKVWRVGKMISFTYDEGGGKTGRGAVSEKDAPKELLQMLEK
#=GR DN7A_SACS2/3-61 SS EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT
DN7E_SULAC/3-60 KVRFKYKGEEKEVDTSKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELMDMLAR
#=GC SS_cons EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT
#=GC seq_cons KVKFKYKGEEKEVDISKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLsMLuK
//
或者，以對齊的 FASTA 格式（請參閱第 對齊的 FASTA 節）
>>> print(format(alignment, "fasta"))
>DN7_METS5/4-61
KIKFKYKGQDLEVDISKVKKVWKVGKMVSFTYDD-NGKTGRGAVSEKDAPKELLNMIGK
>DN7A_SACS2/3-61
TVKFKYKGEEKQVDISKIKKVWRVGKMISFTYDEGGGKTGRGAVSEKDAPKELLQMLEK
>DN7E_SULAC/3-60
KVRFKYKGEEKEVDTSKIKKVWRVGKMVSFTYDD-NGKTGRGAVSEKDAPKELMDMLAR
或以 PHYLIP 格式（請參閱第 PHYLIP 輸出檔案 節）
>>> print(format(alignment, "phylip"))
3 59
DN7_METS5/KIKFKYKGQDLEVDISKVKKVWKVGKMVSFTYDD-NGKTGRGAVSEKDAPKELLNMIGK
DN7A_SACS2TVKFKYKGEEKQVDISKIKKVWRVGKMISFTYDEGGGKTGRGAVSEKDAPKELLQMLEK
DN7E_SULACKVRFKYKGEEKEVDTSKIKKVWRVGKMVSFTYDD-NGKTGRGAVSEKDAPKELMDMLAR
比對的一般資訊儲存在 Alignment 物件的 annotations 屬性下，例如
>>> alignment.annotations["identifier"]
'7kD_DNA_binding'
>>> alignment.annotations["clan"]
'CL0049'
>>> alignment.annotations["database references"]
[{'reference': 'INTERPRO; IPR003212;'}, {'reference': 'SCOP; 1sso; fa;'}, {'reference': 'SO; 0000417; polypeptide_domain;'}]
此比對中的個別序列以 SeqRecord 的形式儲存在 alignment.sequences 下，包括與每個序列記錄關聯的任何註解
>>> for record in alignment.sequences:
... print("%s %s %s" % (record.id, record.annotations["accession"], record.dbxrefs))
...
DN7_METS5/4-61 A4YEA2.1 []
DN7A_SACS2/3-61 P61991.2 ['PDB; 1SSO A; 2-60;', 'PDB; 1JIC A; 2-60;', 'PDB; 2CVR A; 2-60;', 'PDB; 1B4O A; 2-60;']
DN7E_SULAC/3-60 P13125.2 []
第二個序列 (DN7A_SACS2/3-61) 的二級結構儲存在 SeqRecord 的 letter_annotations 屬性中
>>> alignment.sequences[0].letter_annotations
{}
>>> alignment.sequences[1].letter_annotations
{'secondary structure': 'EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT'}
>>> alignment.sequences[2].letter_annotations
{}
共識序列和二級結構與整個序列比對相關聯，因此儲存在 Alignment 物件的 column_annotations 屬性中
>>> alignment.column_annotations
{'consensus secondary structure': 'EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT',
'consensus sequence': 'KVKFKYKGEEKEVDISKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLsMLuK'}
PHYLIP 輸出檔案
序列比對的 PHYLIP 格式源自 Joe Felsenstein 的 PHYLogeny Interference Package。PHYLIP 格式的檔案以兩個數字開頭，分別表示列印比對中的列數和欄數。序列比對本身可以是連續格式或交錯格式。前者的範例是 sequential.phy 檔案（在 Biopython 原始碼發行版的 Tests/Phylip/ 中提供）
3 384
CYS1_DICDI -----MKVIL LFVLAVFTVF VSS------- --------RG IPPEEQ---- --------SQ
FLEFQDKFNK KY-SHEEYLE RFEIFKSNLG KIEELNLIAI NHKADTKFGV NKFADLSSDE
FKNYYLNNKE AIFTDDLPVA DYLDDEFINS IPTAFDWRTR G-AVTPVKNQ GQCGSCWSFS
TTGNVEGQHF ISQNKLVSLS EQNLVDCDHE CMEYEGEEAC DEGCNGGLQP NAYNYIIKNG
GIQTESSYPY TAETGTQCNF NSANIGAKIS NFTMIP-KNE TVMAGYIVST GPLAIAADAV
E-WQFYIGGV F-DIPCN--P NSLDHGILIV GYSAKNTIFR KNMPYWIVKN SWGADWGEQG
YIYLRRGKNT CGVSNFVSTS II--
ALEU_HORVU MAHARVLLLA LAVLATAAVA VASSSSFADS NPIRPVTDRA ASTLESAVLG ALGRTRHALR
FARFAVRYGK SYESAAEVRR RFRIFSESLE EVRSTN---- RKGLPYRLGI NRFSDMSWEE
FQATRL-GAA QTCSATLAGN HLMRDA--AA LPETKDWRED G-IVSPVKNQ AHCGSCWTFS
TTGALEAAYT QATGKNISLS EQQLVDCAGG FNNF------ --GCNGGLPS QAFEYIKYNG
GIDTEESYPY KGVNGV-CHY KAENAAVQVL DSVNITLNAE DELKNAVGLV RPVSVAFQVI
DGFRQYKSGV YTSDHCGTTP DDVNHAVLAV GYGVENGV-- ---PYWLIKN SWGADWGDNG
YFKMEMGKNM CAIATCASYP VVAA
CATH_HUMAN ------MWAT LPLLCAGAWL LGV------- -PVCGAAELS VNSLEK---- --------FH
FKSWMSKHRK TY-STEEYHH RLQTFASNWR KINAHN---- NGNHTFKMAL NQFSDMSFAE
IKHKYLWSEP QNCSAT--KS NYLRGT--GP YPPSVDWRKK GNFVSPVKNQ GACGSCWTFS
TTGALESAIA IATGKMLSLA EQQLVDCAQD FNNY------ --GCQGGLPS QAFEYILYNK
GIMGEDTYPY QGKDGY-CKF QPGKAIGFVK DVANITIYDE EAMVEAVALY NPVSFAFEVT
QDFMMYRTGI YSSTSCHKTP DKVNHAVLAV GYGEKNGI-- ---PYWIVKN SWGPQWGMNG
YFLIERGKNM CGLAACASYP IPLV
在連續格式中，會先顯示一個序列的完整比對，然後再處理下一個序列。在交錯格式中，不同序列的比對彼此相鄰，例如在檔案 interlaced.phy 中（在 Biopython 原始碼發行版的 Tests/Phylip/ 中提供）
3 384
CYS1_DICDI -----MKVIL LFVLAVFTVF VSS------- --------RG IPPEEQ---- --------SQ
ALEU_HORVU MAHARVLLLA LAVLATAAVA VASSSSFADS NPIRPVTDRA ASTLESAVLG ALGRTRHALR
CATH_HUMAN ------MWAT LPLLCAGAWL LGV------- -PVCGAAELS VNSLEK---- --------FH
FLEFQDKFNK KY-SHEEYLE RFEIFKSNLG KIEELNLIAI NHKADTKFGV NKFADLSSDE
FARFAVRYGK SYESAAEVRR RFRIFSESLE EVRSTN---- RKGLPYRLGI NRFSDMSWEE
FKSWMSKHRK TY-STEEYHH RLQTFASNWR KINAHN---- NGNHTFKMAL NQFSDMSFAE
FKNYYLNNKE AIFTDDLPVA DYLDDEFINS IPTAFDWRTR G-AVTPVKNQ GQCGSCWSFS
FQATRL-GAA QTCSATLAGN HLMRDA--AA LPETKDWRED G-IVSPVKNQ AHCGSCWTFS
IKHKYLWSEP QNCSAT--KS NYLRGT--GP YPPSVDWRKK GNFVSPVKNQ GACGSCWTFS
TTGNVEGQHF ISQNKLVSLS EQNLVDCDHE CMEYEGEEAC DEGCNGGLQP NAYNYIIKNG
TTGALEAAYT QATGKNISLS EQQLVDCAGG FNNF------ --GCNGGLPS QAFEYIKYNG
TTGALESAIA IATGKMLSLA EQQLVDCAQD FNNY------ --GCQGGLPS QAFEYILYNK
GIQTESSYPY TAETGTQCNF NSANIGAKIS NFTMIP-KNE TVMAGYIVST GPLAIAADAV
GIDTEESYPY KGVNGV-CHY KAENAAVQVL DSVNITLNAE DELKNAVGLV RPVSVAFQVI
GIMGEDTYPY QGKDGY-CKF QPGKAIGFVK DVANITIYDE EAMVEAVALY NPVSFAFEVT
E-WQFYIGGV F-DIPCN--P NSLDHGILIV GYSAKNTIFR KNMPYWIVKN SWGADWGEQG
DGFRQYKSGV YTSDHCGTTP DDVNHAVLAV GYGVENGV-- ---PYWLIKN SWGADWGDNG
QDFMMYRTGI YSSTSCHKTP DKVNHAVLAV GYGEKNGI-- ---PYWIVKN SWGPQWGMNG
YIYLRRGKNT CGVSNFVSTS II--
YFKMEMGKNM CAIATCASYP VVAA
YFLIERGKNM CGLAACASYP IPLV
Bio.Align 中的解析器會從檔案內容偵測其為序列格式或交錯格式，然後進行適當的解析。
>>> from Bio import Align
>>> alignment = Align.read("sequential.phy", "phylip")
>>> alignment
<Alignment object (3 rows x 384 columns) at ...>
>>> alignment2 = Align.read("interlaced.phy", "phylip")
>>> alignment2
<Alignment object (3 rows x 384 columns) at ...>
>>> alignment == alignment2
True
這裡，如果兩個比對具有相同的序列內容和相同的比對座標，則認為它們相等。
>>> alignment.shape
(3, 384)
>>> print(alignment)
CYS1_DICD 0 -----MKVILLFVLAVFTVFVSS---------------RGIPPEEQ------------SQ
ALEU_HORV 0 MAHARVLLLALAVLATAAVAVASSSSFADSNPIRPVTDRAASTLESAVLGALGRTRHALR
CATH_HUMA 0 ------MWATLPLLCAGAWLLGV--------PVCGAAELSVNSLEK------------FH
CYS1_DICD 28 FLEFQDKFNKKY-SHEEYLERFEIFKSNLGKIEELNLIAINHKADTKFGVNKFADLSSDE
ALEU_HORV 60 FARFAVRYGKSYESAAEVRRRFRIFSESLEEVRSTN----RKGLPYRLGINRFSDMSWEE
CATH_HUMA 34 FKSWMSKHRKTY-STEEYHHRLQTFASNWRKINAHN----NGNHTFKMALNQFSDMSFAE
CYS1_DICD 87 FKNYYLNNKEAIFTDDLPVADYLDDEFINSIPTAFDWRTRG-AVTPVKNQGQCGSCWSFS
ALEU_HORV 116 FQATRL-GAAQTCSATLAGNHLMRDA--AALPETKDWREDG-IVSPVKNQAHCGSCWTFS
CATH_HUMA 89 IKHKYLWSEPQNCSAT--KSNYLRGT--GPYPPSVDWRKKGNFVSPVKNQGACGSCWTFS
CYS1_DICD 146 TTGNVEGQHFISQNKLVSLSEQNLVDCDHECMEYEGEEACDEGCNGGLQPNAYNYIIKNG
ALEU_HORV 172 TTGALEAAYTQATGKNISLSEQQLVDCAGGFNNF--------GCNGGLPSQAFEYIKYNG
CATH_HUMA 145 TTGALESAIAIATGKMLSLAEQQLVDCAQDFNNY--------GCQGGLPSQAFEYILYNK
CYS1_DICD 206 GIQTESSYPYTAETGTQCNFNSANIGAKISNFTMIP-KNETVMAGYIVSTGPLAIAADAV
ALEU_HORV 224 GIDTEESYPYKGVNGV-CHYKAENAAVQVLDSVNITLNAEDELKNAVGLVRPVSVAFQVI
CATH_HUMA 197 GIMGEDTYPYQGKDGY-CKFQPGKAIGFVKDVANITIYDEEAMVEAVALYNPVSFAFEVT
CYS1_DICD 265 E-WQFYIGGVF-DIPCN--PNSLDHGILIVGYSAKNTIFRKNMPYWIVKNSWGADWGEQG
ALEU_HORV 283 DGFRQYKSGVYTSDHCGTTPDDVNHAVLAVGYGVENGV-----PYWLIKNSWGADWGDNG
CATH_HUMA 256 QDFMMYRTGIYSSTSCHKTPDKVNHAVLAVGYGEKNGI-----PYWIVKNSWGPQWGMNG
CYS1_DICD 321 YIYLRRGKNTCGVSNFVSTSII-- 343
ALEU_HORV 338 YFKMEMGKNMCAIATCASYPVVAA 362
CATH_HUMA 311 YFLIERGKNMCGLAACASYPIPLV 335
當以 PHYLIP 格式輸出比對時，Bio.Align 會將每個比對的序列寫在單獨的一行上。
>>> print(format(alignment, "phylip"))
3 384
CYS1_DICDI-----MKVILLFVLAVFTVFVSS---------------RGIPPEEQ------------SQFLEFQDKFNKKY-SHEEYLERFEIFKSNLGKIEELNLIAINHKADTKFGVNKFADLSSDEFKNYYLNNKEAIFTDDLPVADYLDDEFINSIPTAFDWRTRG-AVTPVKNQGQCGSCWSFSTTGNVEGQHFISQNKLVSLSEQNLVDCDHECMEYEGEEACDEGCNGGLQPNAYNYIIKNGGIQTESSYPYTAETGTQCNFNSANIGAKISNFTMIP-KNETVMAGYIVSTGPLAIAADAVE-WQFYIGGVF-DIPCN--PNSLDHGILIVGYSAKNTIFRKNMPYWIVKNSWGADWGEQGYIYLRRGKNTCGVSNFVSTSII--
ALEU_HORVUMAHARVLLLALAVLATAAVAVASSSSFADSNPIRPVTDRAASTLESAVLGALGRTRHALRFARFAVRYGKSYESAAEVRRRFRIFSESLEEVRSTN----RKGLPYRLGINRFSDMSWEEFQATRL-GAAQTCSATLAGNHLMRDA--AALPETKDWREDG-IVSPVKNQAHCGSCWTFSTTGALEAAYTQATGKNISLSEQQLVDCAGGFNNF--------GCNGGLPSQAFEYIKYNGGIDTEESYPYKGVNGV-CHYKAENAAVQVLDSVNITLNAEDELKNAVGLVRPVSVAFQVIDGFRQYKSGVYTSDHCGTTPDDVNHAVLAVGYGVENGV-----PYWLIKNSWGADWGDNGYFKMEMGKNMCAIATCASYPVVAA
CATH_HUMAN------MWATLPLLCAGAWLLGV--------PVCGAAELSVNSLEK------------FHFKSWMSKHRKTY-STEEYHHRLQTFASNWRKINAHN----NGNHTFKMALNQFSDMSFAEIKHKYLWSEPQNCSAT--KSNYLRGT--GPYPPSVDWRKKGNFVSPVKNQGACGSCWTFSTTGALESAIAIATGKMLSLAEQQLVDCAQDFNNY--------GCQGGLPSQAFEYILYNKGIMGEDTYPYQGKDGY-CKFQPGKAIGFVKDVANITIYDEEAMVEAVALYNPVSFAFEVTQDFMMYRTGIYSSTSCHKTPDKVNHAVLAVGYGEKNGI-----PYWIVKNSWGPQWGMNGYFLIERGKNMCGLAACASYPIPLV
我們可以以 PHYLIP 格式寫入比對，解析結果，並確認它與原始比對物件相同。
>>> from io import StringIO
>>> stream = StringIO()
>>> Align.write(alignment, stream, "phylip")
1
>>> stream.seek(0)
0
>>> alignment3 = Align.read(stream, "phylip")
>>> alignment == alignment3
True
>>> [record.id for record in alignment.sequences]
['CYS1_DICDI', 'ALEU_HORVU', 'CATH_HUMAN']
>>> [record.id for record in alignment3.sequences]
['CYS1_DICDI', 'ALEU_HORVU', 'CATH_HUMAN']
EMBOSS
EMBOSS（European Molecular Biology Open Software Suite）是一套用於分子生物學和生物資訊學的開源軟體工具 [Rice2000]。它包括用於成對序列比對的軟體，如 needle 和 water。這是一個由 water 程式生成的 Smith-Waterman 局部成對序列比對輸出的範例（在 Biopython 發行版的 Tests/Emboss 目錄中以 water.txt 形式提供）。
########################################
# Program: water
# Rundate: Wed Jan 16 17:23:19 2002
# Report_file: stdout
########################################
#=======================================
#
# Aligned_sequences: 2
# 1: IXI_234
# 2: IXI_235
# Matrix: EBLOSUM62
# Gap_penalty: 10.0
# Extend_penalty: 0.5
#
# Length: 131
# Identity: 112/131 (85.5%)
# Similarity: 112/131 (85.5%)
# Gaps: 19/131 (14.5%)
# Score: 591.5
#
#
#=======================================
IXI_234 1 TSPASIRPPAGPSSRPAMVSSRRTRPSPPGPRRPTGRPCCSAAPRRPQAT 50
||||||||||||||| ||||||||||||||||||||||||||
IXI_235 1 TSPASIRPPAGPSSR---------RPSPPGPRRPTGRPCCSAAPRRPQAT 41
IXI_234 51 GGWKTCSGTCTTSTSTRHRGRSGWSARTTTAACLRASRKSMRAACSRSAG 100
|||||||||||||||||||||||| ||||||||||||||||
IXI_235 42 GGWKTCSGTCTTSTSTRHRGRSGW----------RASRKSMRAACSRSAG 81
IXI_234 101 SRPNRFAPTLMSSCITSTTGPPAWAGDRSHE 131
|||||||||||||||||||||||||||||||
IXI_235 82 SRPNRFAPTLMSSCITSTTGPPAWAGDRSHE 112
#---------------------------------------
#---------------------------------------
由於此輸出檔案僅包含一個比對，我們可以藉由 Align.read 直接提取它。這裡，我們將使用 Align.parse，以便查看此 water 運行的元數據。
>>> from Bio import Align
>>> alignments = Align.parse("water.txt", "emboss")
alignments 的 metadata 屬性會儲存檔案標頭中顯示的資訊，包括用於產生輸出的程式、程式運行的日期和時間、輸出檔案名稱，以及所使用的特定比對檔案格式（預設為 srspair）。
>>> alignments.metadata
{'Align_format': 'srspair', 'Program': 'water', 'Rundate': 'Wed Jan 16 17:23:19 2002', 'Report_file': 'stdout'}
要提取比對，我們使用
>>> alignment = next(alignments)
>>> alignment
<Alignment object (2 rows x 131 columns) at ...>
>>> alignment.shape
(2, 131)
>>> print(alignment)
IXI_234 0 TSPASIRPPAGPSSRPAMVSSRRTRPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTC
0 |||||||||||||||---------||||||||||||||||||||||||||||||||||||
IXI_235 0 TSPASIRPPAGPSSR---------RPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTC
IXI_234 60 TTSTSTRHRGRSGWSARTTTAACLRASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTG
60 ||||||||||||||----------||||||||||||||||||||||||||||||||||||
IXI_235 51 TTSTSTRHRGRSGW----------RASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTG
IXI_234 120 PPAWAGDRSHE 131
120 ||||||||||| 131
IXI_235 101 PPAWAGDRSHE 112
>>> print(alignment.coordinates)
[[ 0 15 24 74 84 131]
[ 0 15 15 65 65 112]]
我們可以使用索引來提取比對的特定部分。
>>> alignment[0]
'TSPASIRPPAGPSSRPAMVSSRRTRPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTCTTSTSTRHRGRSGWSARTTTAACLRASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTGPPAWAGDRSHE'
>>> alignment[1]
'TSPASIRPPAGPSSR---------RPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTCTTSTSTRHRGRSGW----------RASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTGPPAWAGDRSHE'
>>> alignment[1, 10:30]
'GPSSR---------RPSPPG'
alignment 的 annotations 屬性儲存與此比對相關的特定資訊。
>>> alignment.annotations
{'Matrix': 'EBLOSUM62', 'Gap_penalty': 10.0, 'Extend_penalty': 0.5, 'Identity': 112, 'Similarity': 112, 'Gaps': 19, 'Score': 591.5}
也可以藉由在 alignment 物件上呼叫 counts 方法來取得 gap、相同核苷酸和不匹配的數量。
>>> alignment.counts()
AlignmentCounts(gaps=19, identities=112, mismatches=0)
其中 AlignmentCounts 是 Python 標準函式庫中 collections 模組的 namedtuple。
兩個序列之間顯示的 consensus 行儲存在 column_annotations 屬性中。
>>> alignment.column_annotations
{'emboss_consensus': '||||||||||||||| |||||||||||||||||||||||||||||||||||||||||||||||||| |||||||||||||||||||||||||||||||||||||||||||||||'}
使用 format 函數（或 format 方法）以其他格式列印比對，例如以 PHYLIP 格式列印（請參閱PHYLIP 輸出檔案章節）。
>>> print(format(alignment, "phylip"))
2 131
IXI_234 TSPASIRPPAGPSSRPAMVSSRRTRPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTCTTSTSTRHRGRSGWSARTTTAACLRASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTGPPAWAGDRSHE
IXI_235 TSPASIRPPAGPSSR---------RPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTCTTSTSTRHRGRSGW----------RASRKSMRAACSRSAGSRPNRFAPTLMSSCITSTTGPPAWAGDRSHE
我們可以使用 alignment.sequences 來取得個別序列。但是，由於這是一個成對比對，我們也可以使用 alignment.target 和 alignment.query 來取得目標序列和查詢序列。
>>> alignment.target
SeqRecord(seq=Seq('TSPASIRPPAGPSSRPAMVSSRRTRPSPPGPRRPTGRPCCSAAPRRPQATGGWK...SHE'), id='IXI_234', name='<unknown name>', description='<unknown description>', dbxrefs=[])
>>> alignment.query
SeqRecord(seq=Seq('TSPASIRPPAGPSSRRPSPPGPRRPTGRPCCSAAPRRPQATGGWKTCSGTCTTS...SHE'), id='IXI_235', name='<unknown name>', description='<unknown description>', dbxrefs=[])
目前，Biopython 不支援以 EMBOSS 定義的輸出格式寫入序列比對。
GCG 多重序列格式 (MSF)
多重序列格式 (MSF) 是為儲存由 GCG（遺傳電腦組）程式集產生的多重序列比對而創建的。Biopython 發行版的 Tests/msf 目錄中的 W_prot.msf 檔案是 MSF 格式的序列比對檔案範例。該檔案顯示了 11 個蛋白質序列的比對。
!!AA_MULTIPLE_ALIGNMENT
MSF: 99 Type: P Oct 18, 2017 11:35 Check: 0 ..
Name: W*01:01:01:01 Len: 99 Check: 7236 Weight: 1.00
Name: W*01:01:01:02 Len: 99 Check: 7236 Weight: 1.00
Name: W*01:01:01:03 Len: 99 Check: 7236 Weight: 1.00
Name: W*01:01:01:04 Len: 99 Check: 7236 Weight: 1.00
Name: W*01:01:01:05 Len: 99 Check: 7236 Weight: 1.00
Name: W*01:01:01:06 Len: 99 Check: 7236 Weight: 1.00
Name: W*02:01 Len: 93 Check: 9483 Weight: 1.00
Name: W*03:01:01:01 Len: 93 Check: 9974 Weight: 1.00
Name: W*03:01:01:02 Len: 93 Check: 9974 Weight: 1.00
Name: W*04:01 Len: 93 Check: 9169 Weight: 1.00
Name: W*05:01 Len: 99 Check: 7331 Weight: 1.00
//
W*01:01:01:01 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:02 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:03 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:04 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:05 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:06 GLTPFNGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*02:01 GLTPSNGYTA ATWTRTAASS VGMNIPYDGA SYLVRNQELR SWTAADKAAQ
W*03:01:01:01 GLTPSSGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*03:01:01:02 GLTPSSGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*04:01 GLTPSNGYTA ATWTRTAASS VGMNIPYDGA SYLVRNQELR SWTAADKAAQ
W*05:01 GLTPSSGYTA ATWTRTAVSS VGMNIPYHGA SYLVRNQELR SWTAADKAAQ
W*01:01:01:01 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*01:01:01:02 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*01:01:01:03 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*01:01:01:04 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*01:01:01:05 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*01:01:01:06 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
W*02:01 MPWRRNMQSC SKPTCREGGR SGSAKSLRMG RRRCTAQNPK RLT
W*03:01:01:01 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK RLT
W*03:01:01:02 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK RLT
W*04:01 MPWRRNMQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK RLT
W*05:01 MPWRRNRQSC SKPTCREGGR SGSAKSLRMG RRGCSAQNPK DSHDPPPHL
要使用 Biopython 解析此檔案，請使用
>>> from Bio import Align
>>> alignment = Align.read("W_prot.msf", "msf")
解析器會略過所有直到（包含）以 “MSF:” 開頭的行。解析器會讀取後續行（直到 “//” 分隔符號），以驗證每個序列的長度。解析器會讀取比對區段（在 “//” 分隔符號之後），並將其儲存為 Alignment 物件。
>>> alignment
<Alignment object (11 rows x 99 columns) at ...>
>>> print(alignment)
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 0 GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*02:01 0 GLTPSNGYTAATWTRTAASSVGMNIPYDGASYLVRNQELRSWTAADKAAQMPWRRNMQSC
W*03:01:0 0 GLTPSSGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*03:01:0 0 GLTPSSGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*04:01 0 GLTPSNGYTAATWTRTAASSVGMNIPYDGASYLVRNQELRSWTAADKAAQMPWRRNMQSC
W*05:01 0 GLTPSSGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWRRNRQSC
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*01:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
W*02:01 60 SKPTCREGGRSGSAKSLRMGRRRCTAQNPKRLT------ 93
W*03:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKRLT------ 93
W*03:01:0 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKRLT------ 93
W*04:01 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKRLT------ 93
W*05:01 60 SKPTCREGGRSGSAKSLRMGRRGCSAQNPKDSHDPPPHL 99
序列及其名稱儲存在 alignment.sequences 屬性中。
>>> len(alignment.sequences)
11
>>> alignment.sequences[0].id
'W*01:01:01:01'
>>> alignment.sequences[0].seq
Seq('GLTPFNGYTAATWTRTAVSSVGMNIPYHGASYLVRNQELRSWTAADKAAQMPWR...PHL')
比對座標儲存在 alignment.coordinates 屬性中。
>>> print(alignment.coordinates)
[[ 0 93 99]
[ 0 93 99]
[ 0 93 99]
[ 0 93 99]
[ 0 93 99]
[ 0 93 99]
[ 0 93 93]
[ 0 93 93]
[ 0 93 93]
[ 0 93 93]
[ 0 93 99]]
目前，Biopython 不支援以 MSF 格式寫入序列比對。
Exonerate
Exonerate 是一個用於成對序列比對的通用程式 [Slater2005]。Exonerate 找到的序列比對可以以人類可讀的形式、以 “cigar”（Compact Idiosyncratic Gapped Alignment Report）格式或以 “vulgar”（Verbose Useful Labelled Gapped Alignment Report）格式輸出。使用者可以要求在輸出中包含其中一種或多種格式。Bio.Align 中的解析器只能解析 cigar 或 vulgar 格式的比對，而不會解析包含人類可讀格式比對的輸出。
Biopython 測試套件中的 exn_22_m_cdna2genome_vulgar.exn 檔案是一個以 vulgar 格式顯示比對的 Exonerate 輸出檔案範例。
Command line: [exonerate -m cdna2genome ../scer_cad1.fa /media/Waterloo/Downloads/genomes/scer_s288c/scer_s288c.fa --bestn 3 --showalignment no --showcigar no --showvulgar yes]
Hostname: [blackbriar]
vulgar: gi|296143771|ref|NM_001180731.1| 0 1230 + gi|330443520|ref|NC_001136.10| 1319275 1318045 - 6146 M 1 1 C 3 3 M 1226 1226
vulgar: gi|296143771|ref|NM_001180731.1| 1230 0 - gi|330443520|ref|NC_001136.10| 1318045 1319275 + 6146 M 129 129 C 3 3 M 1098 1098
vulgar: gi|296143771|ref|NM_001180731.1| 0 516 + gi|330443688|ref|NC_001145.3| 85010 667216 + 518 M 11 11 G 1 0 M 15 15 G 2 0 M 4 4 G 1 0 M 1 1 G 1 0 M 8 8 G 4 0 M 17 17 5 0 2 I 0 168904 3 0 2 M 4 4 G 0 1 M 8 8 G 2 0 M 3 3 G 1 0 M 33 33 G 0 2 M 7 7 G 0 1 M 102 102 5 0 2 I 0 96820 3 0 2 M 14 14 G 0 2 M 10 10 G 2 0 M 5 5 G 0 2 M 10 10 G 2 0 M 4 4 G 0 1 M 20 20 G 1 0 M 15 15 G 0 1 M 5 5 G 3 0 M 4 4 5 0 2 I 0 122114 3 0 2 M 20 20 G 0 5 M 6 6 5 0 2 I 0 193835 3 0 2 M 12 12 G 0 2 M 5 5 G 1 0 M 7 7 G 0 2 M 1 1 G 0 1 M 12 12 C 75 75 M 6 6 G 1 0 M 4 4 G 0 1 M 2 2 G 0 1 M 3 3 G 0 1 M 41 41
-- completed exonerate analysis
此檔案包含三個比對。要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("exn_22_m_cdna2genome_vulgar.exn", "exonerate")
字典 alignments.metadata 儲存關於這些比對的常規資訊，顯示在輸出檔案的頂部。
>>> alignments.metadata
{'Program': 'exonerate',
'Command line': 'exonerate -m cdna2genome ../scer_cad1.fa /media/Waterloo/Downloads/genomes/scer_s288c/scer_s288c.fa --bestn 3 --showalignment no --showcigar no --showvulgar yes',
'Hostname': 'blackbriar'}
現在我們可以迭代比對。第一個比對，比對分數為 6146.0，沒有 gap。
>>> alignment = next(alignments)
>>> alignment.score
6146.0
>>> print(alignment.coordinates)
[[1319275 1319274 1319271 1318045]
[ 0 1 4 1230]]
>>> print(alignment)
gi|330443 1319275 ????????????????????????????????????????????????????????????
0 ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
gi|296143 0 ????????????????????????????????????????????????????????????
...
gi|330443 1318075 ?????????????????????????????? 1318045
1200 |||||||||||||||||||||||||||||| 1230
gi|296143 1200 ?????????????????????????????? 1230
請注意，列印的比對中的目標（第一個序列）位於反向鏈，而查詢（第二個序列）位於正向鏈，目標座標遞減，而查詢座標遞增。使用 Python 內建的 format 函數以 exonerate 格式列印此比對會寫入 vulgar 行。
>>> print(format(alignment, "exonerate"))
vulgar: gi|296143771|ref|NM_001180731.1| 0 1230 + gi|330443520|ref|NC_001136.10| 1319275 1318045 - 6146 M 1 1 C 3 3 M 1226 1226
使用 format 方法允許我們請求 vulgar 行（預設）或 cigar 行。
>>> print(alignment.format("exonerate", "vulgar"))
vulgar: gi|296143771|ref|NM_001180731.1| 0 1230 + gi|330443520|ref|NC_001136.10| 1319275 1318045 - 6146 M 1 1 C 3 3 M 1226 1226
>>> print(alignment.format("exonerate", "cigar"))
cigar: gi|296143771|ref|NM_001180731.1| 0 1230 + gi|330443520|ref|NC_001136.10| 1319275 1318045 - 6146 M 1 M 3 M 1226
vulgar 行包含關於比對的資訊（在 M 1 1 C 3 3 M 1226 區段中），而 cigar 行 M 1 M 3 M 1226 中缺少這些資訊。vulgar 行指定比對以單個比對的核苷酸開始，後接三個形成密碼子 (C) 的比對核苷酸，後接 1226 個比對的核苷酸。在 cigar 行中，我們看到單個比對的核苷酸，後接三個比對的核苷酸，後接 1226 個比對的核苷酸；它未指定三個比對的核苷酸形成密碼子。來自 vulgar 行的此資訊儲存在 operations 屬性中。
>>> alignment.operations
bytearray(b'MCM')
有關其他操作碼的定義，請參閱 Exonerate 文件。
同樣，在呼叫 Bio.Align.write 寫入具有 vulgar 或 cigar 比對行的檔案時，可以使用 "vulgar" 或 "cigar" 引數。
我們可以以 BED 和 PSL 格式列印比對。
>>> print(format(alignment, "bed"))
gi|330443520|ref|NC_001136.10| 1318045 1319275 gi|296143771|ref|NM_001180731.1| 6146 - 1318045 1319275 0 3 1226,3,1, 0,1226,1229,
>>> print(format(alignment, "psl"))
1230 0 0 0 0 0 0 0 - gi|296143771|ref|NM_001180731.1| 1230 0 1230 gi|330443520|ref|NC_001136.10| 1319275 1318045 1319275 3 1226,3,1, 0,1226,1229, 1318045,1319271,1319274,
SAM 格式解析器定義其自己的（可選）operations 屬性（請參閱序列比對/對應 (SAM)章節），這與 Exonerate 格式解析器中定義的 operations 屬性不太一致。由於 operations 屬性是可選的，我們會在以 SAM 格式列印比對之前將其刪除。
>>> del alignment.operations
>>> print(format(alignment, "sam"))
gi|296143771|ref|NM_001180731.1| 16 gi|330443520|ref|NC_001136.10| 1318046 255 1226M3M1M * 0 0 * * AS:i:6146
第三個比對包含四個長 gap。
>>> alignment = next(alignments) # second alignment
>>> alignment = next(alignments) # third alignment
>>> print(alignment)
gi|330443 85010 ???????????-???????????????--????-?-????????----????????????
0 |||||||||||-|||||||||||||||--||||-|-||||||||----||||||||||||
gi|296143 0 ????????????????????????????????????????????????????????????
gi|330443 85061 ????????????????????????????????????????????????????????????
60 |||||-------------------------------------------------------
gi|296143 60 ?????-------------------------------------------------------
...
gi|330443 666990 ????????????????????????????????????????????????????????????
582000 --------------------------------------------------||||||||||
gi|296143 346 --------------------------------------------------??????????
gi|330443 667050 ?????????-??????????????????????????????????????????????????
582060 ||--|||||-|||||||--|-|||||||||||||||||||||||||||||||||||||||
gi|296143 356 ??--?????????????--?-???????????????????????????????????????
gi|330443 667109 ??????????????????????????????????????????????????????-?????
582120 ||||||||||||||||||||||||||||||||||||||||||||||||||||||-||||-
gi|296143 411 ???????????????????????????????????????????????????????????-
gi|330443 667168 ???????????????????????????????????????????????? 667216
582180 ||-|||-||||||||||||||||||||||||||||||||||||||||| 582228
gi|296143 470 ??-???-????????????????????????????????????????? 516
>>> print(format(alignment, "exonerate"))
vulgar: gi|296143771|ref|NM_001180731.1| 0 516 + gi|330443688|ref|NC_001145.3|
85010 667216 + 518 M 11 11 G 1 0 M 15 15 G 2 0 M 4 4 G 1 0 M 1 1 G 1 0 M 8 8
G 4 0 M 17 17 5 0 2 I 0 168904 3 0 2 M 4 4 G 0 1 M 8 8 G 2 0 M 3 3 G 1 0
M 33 33 G 0 2 M 7 7 G 0 1 M 102 102 5 0 2 I 0 96820 3 0 2 M 14 14 G 0 2 M 10 10
G 2 0 M 5 5 G 0 2 M 10 10 G 2 0 M 4 4 G 0 1 M 20 20 G 1 0 M 15 15 G 0 1 M 5 5
G 3 0 M 4 4 5 0 2 I 0 122114 3 0 2 M 20 20 G 0 5 M 6 6 5 0 2 I 0 193835 3 0 2
M 12 12 G 0 2 M 5 5 G 1 0 M 7 7 G 0 2 M 1 1 G 0 1 M 12 12 C 75 75 M 6 6 G 1 0
M 4 4 G 0 1 M 2 2 G 0 1 M 3 3 G 0 1 M 41 41
NEXUS
NEXUS 檔案格式 [Maddison1997] 由多個程式用於儲存親緣關係資訊。這是一個 NEXUS 格式檔案的範例（在 Biopython 發行版的 Tests/Nexus 子目錄中以 codonposset.nex 形式提供）。
#NEXUS
[MacClade 4.05 registered to Computational Biologist, University]
BEGIN DATA;
DIMENSIONS NTAX=2 NCHAR=22;
FORMAT DATATYPE=DNA MISSING=? GAP=- ;
MATRIX
[ 10 20 ]
[ . . ]
Aegotheles AAAAAGGCATTGTGGTGGGAAT [22]
Aerodramus ?????????TTGTGGTGGGAAT [13]
;
END;
BEGIN CODONS;
CODONPOSSET * CodonPositions =
N: 1-10,
1: 11-22\3,
2: 12-22\3,
3: 13-22\3;
CODESET * UNTITLED = Universal: all ;
END;
一般來說，NEXUS 格式的檔案可能更加複雜。Bio.Align 在很大程度上依賴於 Bio.Nexus 中的 NEXUS 解析器（請參閱使用 Bio.Phylo 進行親緣關係分析章節）從 NEXUS 檔案中提取 Alignment 物件。
要讀取此 NEXUS 檔案中的比對，請使用
>>> from Bio import Align
>>> alignment = Align.read("codonposset.nex", "nexus")
>>> print(alignment)
Aegothele 0 AAAAAGGCATTGTGGTGGGAAT 22
0 .........||||||||||||| 22
Aerodramu 0 ?????????TTGTGGTGGGAAT 22
>>> alignment.shape
(2, 22)
序列儲存在 sequences 屬性下。
>>> alignment.sequences[0].id
'Aegotheles'
>>> alignment.sequences[0].seq
Seq('AAAAAGGCATTGTGGTGGGAAT')
>>> alignment.sequences[0].annotations
{'molecule_type': 'DNA'}
>>> alignment.sequences[1].id
'Aerodramus'
>>> alignment.sequences[1].seq
Seq('?????????TTGTGGTGGGAAT')
>>> alignment.sequences[1].annotations
{'molecule_type': 'DNA'}
要以 NEXUS 格式列印此比對，請使用
>>> print(format(alignment, "nexus"))
#NEXUS
begin data;
dimensions ntax=2 nchar=22;
format datatype=dna missing=? gap=-;
matrix
Aegotheles AAAAAGGCATTGTGGTGGGAAT
Aerodramus ?????????TTGTGGTGGGAAT
;
end;
同樣，您可以使用 Align.write(alignment, "myfilename.nex", "nexus") 以 NEXUS 格式將比對寫入 myfilename.nex 檔案。
來自 BLAST 或 FASTA 的表格輸出
表格輸出的比對結果是由 FASTA 比對器 [Pearson1988] 搭配 -m 8CB 或 -m 8CC 參數，或是由 BLAST [Altschul1990] 搭配 -outfmt 7 參數執行產生。
Biopython 原始碼發佈版的 Tests/Fasta 子目錄中的 nucleotide_m8CC.txt 檔案，是 FASTA 使用 -m 8CC 參數所產生的輸出檔案範例。
# fasta36 -m 8CC seq/mgstm1.nt seq/gst.nlib
# FASTA 36.3.8h May, 2020
# Query: pGT875 - 657 nt
# Database: seq/gst.nlib
# Fields: query id, subject id, % identity, alignment length, mismatches, gap opens, q. start, q. end, s. start, s. end, evalue, bit score, aln_code
# 12 hits found
pGT875 pGT875 100.00 657 0 0 1 657 38 694 4.6e-191 655.6 657M
pGT875 RABGLTR 79.10 646 135 0 1 646 34 679 1.6e-116 408.0 646M
pGT875 BTGST 59.56 413 167 21 176 594 228 655 1.9e-07 45.7 149M1D7M1I17M3D60M5I6M1I13M2I13M4I30M2I6M2D112M
pGT875 RABGSTB 66.93 127 42 8 159 289 157 287 3.2e-07 45.0 15M2I17M2D11M1I58M1I11M1D7M1D8M
pGT875 OCDHPR 91.30 23 2 1 266 289 2303 2325 0.012 29.7 17M1D6M
...
# FASTA processed 1 queries
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("nucleotide_m8CC.txt", "tabular")
檔案標頭中顯示的資訊會儲存在 alignments 的 metadata 屬性中。
>>> alignments.metadata
{'Command line': 'fasta36 -m 8CC seq/mgstm1.nt seq/gst.nlib',
'Program': 'FASTA',
'Version': '36.3.8h May, 2020',
'Database': 'seq/gst.nlib'}
藉由迭代 alignments 來提取特定的比對結果。舉例來說，我們來看第四個比對結果。
>>> alignment = next(alignments)
>>> alignment = next(alignments)
>>> alignment = next(alignments)
>>> alignment = next(alignments)
>>> print(alignment)
RABGSTB 156 ??????????????????????????????????--????????????????????????
0 |||||||||||||||--|||||||||||||||||--|||||||||||-||||||||||||
pGT875 158 ???????????????--??????????????????????????????-????????????
RABGSTB 214 ??????????????????????????????????????????????????????????-?
60 ||||||||||||||||||||||||||||||||||||||||||||||-|||||||||||-|
pGT875 215 ??????????????????????????????????????????????-?????????????
RABGSTB 273 ??????-???????? 287
120 ||||||-|||||||| 135
pGT875 274 ??????????????? 289
>>> print(alignment.coordinates)
[[156 171 173 190 190 201 202 260 261 272 272 279 279 287]
[158 173 173 190 192 203 203 261 261 272 273 280 281 289]]
>>> alignment.aligned
array([[[156, 171],
[173, 190],
[190, 201],
[202, 260],
[261, 272],
[272, 279],
[279, 287]],
[[158, 173],
[173, 190],
[192, 203],
[203, 261],
[261, 272],
[273, 280],
[281, 289]]])
目標和查詢序列的序列資訊儲存在 target 和 query 屬性中（也儲存在 alignment.sequences 下）。
>>> alignment.target
SeqRecord(seq=Seq(None, length=287), id='RABGSTB', name='<unknown name>', description='<unknown description>', dbxrefs=[])
>>> alignment.query
SeqRecord(seq=Seq(None, length=657), id='pGT875', name='<unknown name>', description='<unknown description>', dbxrefs=[])
比對結果的資訊儲存在 alignment 的 annotations 屬性中。
>>> alignment.annotations
{'% identity': 66.93,
'mismatches': 42,
'gap opens': 8,
'evalue': 3.2e-07,
'bit score': 45.0}
BLAST 特別提供了許多選項，可以藉由包含或排除特定欄位來自訂表格輸出；詳細資訊請參閱 BLAST 文件。這些資訊會以字典形式儲存在 alignment.annotations、alignment.target.annotations 或 alignment.query.annotations 中，視情況而定。
HH-suite 輸出檔案
格式為 hhr 的比對檔案是由 HH-suite 中的 hhsearch 或 hhblits 產生 [Steinegger2019]。例如，請參閱 Biopython 測試套件中的 2uvo_hhblits.hhr 檔案。
Query 2UVO:A|PDBID|CHAIN|SEQUENCE
Match_columns 171
No_of_seqs 1560 out of 4005
Neff 8.3
Searched_HMMs 34
Date Fri Feb 15 16:34:13 2019
Command hhblits -i 2uvoAh.fasta -d /pdb70
No Hit Prob E-value P-value Score SS Cols Query HMM Template HMM
1 2uvo_A Agglutinin isolectin 1; 100.0 3.7E-34 4.8E-38 210.3 0.0 171 1-171 1-171 (171)
2 2wga ; lectin (agglutinin); 99.9 1.1E-30 1.4E-34 190.4 0.0 162 2-169 2-163 (164)
3 1ulk_A Lectin-C; chitin-bindin 99.8 5.2E-24 6.6E-28 148.2 0.0 120 1-124 2-121 (126)
...
31 4z8i_A BBTPGRP3, peptidoglycan 79.6 0.12 1.5E-05 36.1 0.0 37 1-37 9-54 (236)
32 1wga ; lectin (agglutinin); 40.4 2.6 0.00029 25.9 0.0 106 54-163 11-116 (164)
No 1
>2uvo_A Agglutinin isolectin 1; carbohydrate-binding protein, hevein domain, chitin-binding, GERM agglutinin, chitin-binding protein; HET: NDG NAG GOL; 1.40A {Triticum aestivum} PDB: 1wgc_A* 2cwg_A* 2x3t_A* 4aml_A* 7wga_A 9wga_A 2wgc_A 1wgt_A 1k7t_A* 1k7v_A* 1k7u_A 2x52_A* 1t0w_A*
Probab=99.95 E-value=3.7e-34 Score=210.31 Aligned_cols=171 Identities=100% Similarity=2.050 Sum_probs=166.9
Q 2UVO:A|PDBID|C 1 ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGATCTNNQCCSQYGYCGFGAEYCGAGCQG 80 (171)
Q Consensus 1 ~~cg~~~~~~~c~~~~CCs~~g~CG~~~~~c~~~c~~~~c~~~~~Cg~~~~~~~c~~~~CCs~~g~CG~~~~~c~~~c~~ 80 (171)
||||++.++..||++.|||+|+|||.+.+||+++||.+.|++..+|+++++.++|....|||.++||+.+.+||+.+||.
T Consensus 1 ~~cg~~~~~~~c~~~~CCS~~g~Cg~~~~~Cg~gC~~~~c~~~~~cg~~~~~~~c~~~~CCs~~g~Cg~~~~~c~~~c~~ 80 (171)
T 2uvo_A 1 ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGATCTNNQCCSQYGYCGFGAEYCGAGCQG 80 (171)
T ss_dssp CBCBGGGTTBBCGGGCEECTTSBEEBSHHHHSTTCCBSSCSSCCBCBGGGTTBCCSTTCEECTTSBEEBSHHHHSTTCCB
T ss_pred CCCCCCCCCcCCCCCCeeCCCCeECCCcccccCCccccccccccccCcccCCcccCCccccCCCceeCCCccccCCCccc
Confidence 79999999999999999999999999999999999999999999999999999999999999999999999999999999
Q 2UVO:A|PDBID|C 81 GPCRADIKCGSQAGGKLCPNNLCCSQWGFCGLGSEFCGGGCQSGACSTDKPCGKDAGGRVCTNNYCCSKWGSCGIGPGYC 160 (171)
Q Consensus 81 ~~~~~~~~Cg~~~~~~~c~~~~CCS~~G~CG~~~~~C~~~Cq~~~c~~~~~Cg~~~~~~~c~~~~CCS~~G~CG~~~~~C 160 (171)
+++++|+.|+...+++.||++.|||.|||||...+||+.+||+++|++|.+|++.+++++|..+.|||+++-||+...||
T Consensus 81 ~~~~~~~~cg~~~~~~~c~~~~CCs~~g~CG~~~~~C~~gCq~~~c~~~~~cg~~~~~~~c~~~~ccs~~g~Cg~~~~~C 160 (171)
T 2uvo_A 81 GPCRADIKCGSQAGGKLCPNNLCCSQWGFCGLGSEFCGGGCQSGACSTDKPCGKDAGGRVCTNNYCCSKWGSCGIGPGYC 160 (171)
T ss_dssp SSCSSCCBCBGGGTTBCCGGGCEECTTSBEEBSHHHHSTTCCBSSCSSCCCCBTTTTTBCCSTTCEECTTSCEEBSHHHH
T ss_pred ccccccccccccccCCCCCCCcccCCCCccCCCcccccCCCcCCccccccccccccccccCCCCCCcCCCCEecCchhhc
Confidence 99999999999988999999999999999999999999999999999999999999999999999999999999999999
Q 2UVO:A|PDBID|C 161 GAGCQSGGCDG 171 (171)
Q Consensus 161 ~~gCq~~~c~~ 171 (171)
+++||++.|||
T Consensus 161 ~~~cq~~~~~~ 171 (171)
T 2uvo_A 161 GAGCQSGGCDG 171 (171)
T ss_dssp STTCCBSSCC-
T ss_pred ccccccCCCCC
Confidence 99999999986
No 2
...
No 32
>1wga ; lectin (agglutinin); NMR {}
Probab=40.43 E-value=2.6 Score=25.90 Aligned_cols=106 Identities=20% Similarity=0.652 Sum_probs=54.7
Q 2UVO:A|PDBID|C 54 TCTNNQCCSQYGYCGFGAEYCGAGCQGGPCRADIKCGSQAGGKLCPNNLCCSQWGFCGLGSEFCGGGCQSGACSTDKPCG 133 (171)
Q Consensus 54 ~c~~~~CCs~~g~CG~~~~~c~~~c~~~~~~~~~~Cg~~~~~~~c~~~~CCS~~G~CG~~~~~C~~~Cq~~~c~~~~~Cg 133 (171)
.|....||.....|......|...|....|.....|... ...|....||.....|......|...|....+.....|.
T Consensus 11 ~c~~~~cc~~~~~c~~~~~~c~~~c~~~~c~~~~~c~~~--~~~c~~~~cc~~~~~c~~~~~~c~~~c~~~~c~~~~~c~ 88 (164)
T 1wga 11 XCXXXXCCXXXXXCXXXXXXCXXXCXXXXCXXXXXCXXX--XXXCXXXXCCXXXXXCXXXXXXCXXXCXXXXCXXXXXCX 88 (164)
T ss_pred ccccccccccccccccccccccccccccccccccccccc--ccccccccccccccccccccccccccccccccccccccc
Confidence 344556666666666666566555543333223333321 234666677777777777766666655544332223333
Q 2UVO:A|PDBID|C 134 KDAGGRVCTNNYCCSKWGSCGIGPGYCGAG 163 (171)
Q Consensus 134 ~~~~~~~c~~~~CCS~~G~CG~~~~~C~~g 163 (171)
.. ...|....||.....|......|...
T Consensus 89 ~~--~~~c~~~~cc~~~~~c~~~~~~c~~~ 116 (164)
T 1wga 89 XX--XXXCXXXXCCXXXXXCXXXXXXCXXX 116 (164)
T ss_pred cc--cccccccccccccccccccccccccc
Confidence 22 23344455555555555555544433
Done!
此檔案包含三個部分：
包含比對結果一般資訊的標頭；
包含每個取得的比對結果的一行摘要；
以及連續顯示的詳細比對結果。
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("2uvo_hhblits.hhr", "hhr")
大多數標頭資訊都儲存在 alignments 的 metadata 屬性中。
>>> alignments.metadata
{'Match_columns': 171,
'No_of_seqs': (1560, 4005),
'Neff': 8.3,
'Searched_HMMs': 34,
'Rundate': 'Fri Feb 15 16:34:13 2019',
'Command line': 'hhblits -i 2uvoAh.fasta -d /pdb70'}
除了查詢名稱，它是以屬性形式儲存。
>>> alignments.query_name
'2UVO:A|PDBID|CHAIN|SEQUENCE'
因為它會在每個比對結果中重複出現。
迭代比對結果。
>>> for alignment in alignments:
... print(alignment.target.id)
...
2uvo_A
2wga
1ulk_A
...
4z8i_A
1wga
讓我們更詳細地看看第一個比對結果。
>>> alignments = iter(alignments)
>>> alignment = next(alignments)
>>> alignment
<Alignment object (2 rows x 171 columns) at ...>
>>> print(alignment)
2uvo_A 0 ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGATCTNNQC
0 ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
2UVO:A|PD 0 ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGATCTNNQC
2uvo_A 60 CSQYGYCGFGAEYCGAGCQGGPCRADIKCGSQAGGKLCPNNLCCSQWGFCGLGSEFCGGG
60 ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
2UVO:A|PD 60 CSQYGYCGFGAEYCGAGCQGGPCRADIKCGSQAGGKLCPNNLCCSQWGFCGLGSEFCGGG
2uvo_A 120 CQSGACSTDKPCGKDAGGRVCTNNYCCSKWGSCGIGPGYCGAGCQSGGCDG 171
120 ||||||||||||||||||||||||||||||||||||||||||||||||||| 171
2UVO:A|PD 120 CQSGACSTDKPCGKDAGGRVCTNNYCCSKWGSCGIGPGYCGAGCQSGGCDG 171
目標和查詢序列儲存在 alignment.sequences 中。由於這些是成對比對，我們也可以透過 alignment.target 和 alignment.query 存取它們。
>>> alignment.target is alignment.sequences[0]
True
>>> alignment.query is alignment.sequences[1]
True
查詢的 ID 是從 alignments.query_name 設定（請注意，在 hhr 檔案的比對結果中印出的查詢 ID 已被縮寫）。
>>> alignment.query.id
'2UVO:A|PDBID|CHAIN|SEQUENCE'
目標的 ID 是從序列比對區塊取得（以 T 2uvo_A 開頭的行）。
>>> alignment.target.id
'2uvo_A'
目標和查詢的序列內容是從此比對結果中的可用資訊填入。
>>> alignment.target.seq
Seq('ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGAT...CDG')
>>> alignment.query.seq
Seq('ERCGEQGSNMECPNNLCCSQYGYCGMGGDYCGKGCQNGACWTSKRCGSQAGGAT...CDG')
如果比對結果沒有延伸到整個序列，則序列內容將會不完整（部分定義的序列；請參閱第 具有部分定義序列內容的序列 節）。
此比對區塊的第二行（以「>」開頭）顯示從中取得目標序列的隱馬可夫模型的名稱和描述。這些資訊會儲存在 alignment.target.annotations 字典中的 "hmm_name" 和 "hmm_description" 鍵下。
>>> alignment.target.annotations
{'hmm_name': '2uvo_A',
'hmm_description': 'Agglutinin isolectin 1; carbohydrate-binding protein, hevein domain, chitin-binding, GERM agglutinin, chitin-binding protein; HET: NDG NAG GOL; 1.40A {Triticum aestivum} PDB: 1wgc_A* 2cwg_A* 2x3t_A* 4aml_A* 7wga_A 9wga_A 2wgc_A 1wgt_A 1k7t_A* 1k7v_A* 1k7u_A 2x52_A* 1t0w_A*'}
alignment.target.letter_annotations 字典儲存了目標比對共識序列、PSIPRED 預測的二級結構，以及 DSSP 判定的目標二級結構。
>>> alignment.target.letter_annotations
{'Consensus': '~~cg~~~~~~~c~~~~CCS~~g~Cg~~~~~Cg~gC~~~~c~~~~~cg~~~~~~~c~~~~CCs~~g~Cg~~~~~c~~~c~~~~~~~~~~cg~~~~~~~c~~~~CCs~~g~CG~~~~~C~~gCq~~~c~~~~~cg~~~~~~~c~~~~ccs~~g~Cg~~~~~C~~~cq~~~~~~',
'ss_pred': 'CCCCCCCCCcCCCCCCeeCCCCeECCCcccccCCccccccccccccCcccCCcccCCccccCCCceeCCCccccCCCcccccccccccccccccCCCCCCCcccCCCCccCCCcccccCCCcCCccccccccccccccccCCCCCCcCCCCEecCchhhcccccccCCCCC',
'ss_dssp': 'CBCBGGGTTBBCGGGCEECTTSBEEBSHHHHSTTCCBSSCSSCCBCBGGGTTBCCSTTCEECTTSBEEBSHHHHSTTCCBSSCSSCCBCBGGGTTBCCGGGCEECTTSBEEBSHHHHSTTCCBSSCSSCCCCBTTTTTBCCSTTCEECTTSCEEBSHHHHSTTCCBSSCC '}
在此範例中，對於查詢序列，只有共識序列可用。
>>> alignment.query.letter_annotations
{'Consensus': '~~cg~~~~~~~c~~~~CCs~~g~CG~~~~~c~~~c~~~~c~~~~~Cg~~~~~~~c~~~~CCs~~g~CG~~~~~c~~~c~~~~~~~~~~Cg~~~~~~~c~~~~CCS~~G~CG~~~~~C~~~Cq~~~c~~~~~Cg~~~~~~~c~~~~CCS~~G~CG~~~~~C~~gCq~~~c~~'}
alignment.annotations 字典儲存了比對區塊第三行顯示的比對資訊。
>>> alignment.annotations
{'Probab': 99.95,
'E-value': 3.7e-34,
'Score': 210.31,
'Identities': 100.0,
'Similarity': 2.05,
'Sum_probs': 166.9}
成對比對的信賴值儲存在 alignment.column_annotations 字典中的 "Confidence" 鍵下。此字典也儲存了每個欄位的分數，顯示在每個比對區塊的查詢和目標部分之間。
>>> alignment.column_annotations
{'column score': '||||++.++..||++.|||+|+|||.+.+||+++||.+.|++..+|+++++.++|....|||.++||+.+.+||+.+||.+++++|+.|+...+++.||++.|||.|||||...+||+.+||+++|++|.+|++.+++++|..+.|||+++-||+...||+++||++.|||',
'Confidence': '799999999999999999999999999999999999999999999999999999999999999999999999999999999999999999998899999999999999999999999999999999999999999999999999999999999999999999999999986'}
A2M
A2M 檔案是由 SAM 序列比對和建模軟體系統中的 align2model 或 hmmscore 產生的比對檔案 [Krogh1994], [Hughey1996]。一個 A2M 檔案包含一個多序列比對。A2M 檔案格式類似於對齊的 FASTA（請參閱 對齊的 FASTA 節）。但是，為了區分插入和刪除，A2M 使用破折號和句點來表示間隙，並在對齊的序列中使用大寫和小寫字符。匹配以大寫字母表示，並且僅包含匹配或刪除的比對欄中的刪除以破折號表示。插入以小寫字母表示，與插入對齊的間隙顯示為句點。標頭行以「>」開頭，後接序列名稱，以及可選的描述。
Biopython 測試套件中的 probcons.a2m 檔案是 A2M 檔案的範例（關於相同比對結果的對齊 FASTA 格式，請參閱 對齊的 FASTA 節）。
>plas_horvu
D.VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG.VD.VSKISQEEYLTAPGETFSVTLTV...PGTYGFYCEPHAGAGMVGKVT
V
>plas_chlre
-.VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG.VN.ADAISRDDYLNAPGETYSVKLTA...AGEYGYYCEPHQGAGMVGKII
V
>plas_anava
-.VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKsADlAKSLSHKQLLMSPGQSTSTTFPAdapAGEYTFYCEPHRGAGMVGKIT
V
>plas_proho
VqIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG.ES.APALSNTKLRIAPGSFYSVTLGT...PGTYSFYCTPHRGAGMVGTIT
V
>azup_achcy
VhMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG.AE.A-------FKSKINENYKVTFTA...PGVYGVKCTPHYGMGMVGVVE
V
要解析此比對結果，請使用：
>>> from Bio import Align
>>> alignment = Align.read("probcons.a2m", "a2m")
>>> alignment
<Alignment object (5 rows x 101 columns) at ...>
>>> print(alignment)
plas_horv 0 D-VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG-VD-VSKISQE
plas_chlr 0 --VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG-VN-ADAISRD
plas_anav 0 --VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKSADLAKSLSHK
plas_proh 0 VQIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG-ES-APALSNT
azup_achc 0 VHMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG-AE-A------
plas_horv 57 EYLTAPGETFSVTLTV---PGTYGFYCEPHAGAGMVGKVTV 95
plas_chlr 56 DYLNAPGETYSVKLTA---AGEYGYYCEPHQGAGMVGKIIV 94
plas_anav 58 QLLMSPGQSTSTTFPADAPAGEYTFYCEPHRGAGMVGKITV 99
plas_proh 56 KLRIAPGSFYSVTLGT---PGTYSFYCTPHRGAGMVGTITV 94
azup_achc 51 -FKSKINENYKVTFTA---PGVYGVKCTPHYGMGMVGVVEV 88
解析器會分析 A2M 檔案中破折號、句點以及小寫和大寫字母的模式，以判斷欄位是匹配/不匹配/刪除（「D」）還是插入（「I」）。此資訊儲存在 alignment.column_annotations 字典的 match 鍵下。
>>> alignment.column_annotations
{'state': 'DIDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDIDDIDDDDDDDDDDDDDDDDDDDDDDDIIIDDDDDDDDDDDDDDDDDDDDDD'}
由於狀態資訊儲存在 alignment 中，我們可以以 A2M 格式列印比對結果。
>>> print(format(alignment, "a2m"))
>plas_horvu
D.VLLGANGGVLVFEPNDFSVKAGETITFKNNAGYPHNVVFDEDAVPSG.VD.VSKISQEEYLTAPGETFSVTLTV...PGTYGFYCEPHAGAGMVGKVTV
>plas_chlre
-.VKLGADSGALEFVPKTLTIKSGETVNFVNNAGFPHNIVFDEDAIPSG.VN.ADAISRDDYLNAPGETYSVKLTA...AGEYGYYCEPHQGAGMVGKIIV
>plas_anava
-.VKLGSDKGLLVFEPAKLTIKPGDTVEFLNNKVPPHNVVFDAALNPAKsADlAKSLSHKQLLMSPGQSTSTTFPAdapAGEYTFYCEPHRGAGMVGKITV
>plas_proho
VqIKMGTDKYAPLYEPKALSISAGDTVEFVMNKVGPHNVIFDK--VPAG.ES.APALSNTKLRIAPGSFYSVTLGT...PGTYSFYCTPHRGAGMVGTITV
>azup_achcy
VhMLNKGKDGAMVFEPASLKVAPGDTVTFIPTDK-GHNVETIKGMIPDG.AE.A-------FKSKINENYKVTFTA...PGVYGVKCTPHYGMGMVGVVEV
同樣地，可以使用 Align.write 將比對結果以 A2M 格式寫入輸出檔案（請參閱 寫入比對結果 節）。
Mauve 延伸多重 FASTA (xmfa) 格式
Mauve [Darling2004] 是一個用於建構多重基因體比對的軟體套件。這些比對結果會儲存在延伸多重 FASTA (xmfa) 格式中。根據呼叫 progressiveMauve （Mauve 中的比對程式）的方式，xmfa 格式會略有不同。
如果 progressiveMauve 是使用單一序列輸入檔案呼叫，例如：
progressiveMauve combined.fasta --output=combined.xmfa ...
其中 combined.fasta 包含基因體序列。
>equCab1
GAAAAGGAAAGTACGGCCCGGCCACTCCGGGTGTGTGCTAGGAGGGCTTA
>mm9
GAAGAGGAAAAGTAGATCCCTGGCGTCCGGAGCTGGGACGT
>canFam2
CAAGCCCTGCGCGCTCAGCCGGAGTGTCCCGGGCCCTGCTTTCCTTTTC
那麼輸出檔案 combined.xmfa 如下：
#FormatVersion Mauve1
#Sequence1File combined.fa
#Sequence1Entry 1
#Sequence1Format FastA
#Sequence2File combined.fa
#Sequence2Entry 2
#Sequence2Format FastA
#Sequence3File combined.fa
#Sequence3Entry 3
#Sequence3Format FastA
#BackboneFile combined.xmfa.bbcols
> 1:2-49 - combined.fa
AAGCCCTCCTAGCACACACCCGGAGTGG-CCGGGCCGTACTTTCCTTTT
> 2:0-0 + combined.fa
-------------------------------------------------
> 3:2-48 + combined.fa
AAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGCTTTCCTTTT
=
> 1:1-1 + combined.fa
G
=
> 1:50-50 + combined.fa
A
=
> 2:1-41 + combined.fa
GAAGAGGAAAAGTAGATCCCTGGCGTCCGGAGCTGGGACGT
=
> 3:1-1 + combined.fa
C
=
> 3:49-49 + combined.fa
C
=
數字 (1, 2, 3) 分別代表馬 (equCab1)、鼠 (mm9) 和狗 (canFam2) 的輸入基因體序列。此 xmfa 檔案由六個比對區塊組成，並以 = 字元分隔。使用 Align.parse 來提取這些比對結果。
>>> from Bio import Align
>>> alignments = Align.parse("combined.xmfa", "mauve")
檔案標頭資料儲存在 metadata 屬性中。
>>> alignments.metadata
{'FormatVersion': 'Mauve1',
'BackboneFile': 'combined.xmfa.bbcols',
'File': 'combined.fa'}
identifiers 屬性儲存三個序列的序列識別碼，在此案例中為三個數字。
>>> alignments.identifiers
['0', '1', '2']
這些識別碼會用於個別的比對結果中。
>>> for alignment in alignments:
... print([record.id for record in alignment.sequences])
... print(alignment)
... print("******")
...
['0', '1', '2']
0 49 AAGCCCTCCTAGCACACACCCGGAGTGG-CCGGGCCGTACTTTCCTTTT 1
1 0 ------------------------------------------------- 0
2 1 AAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGCTTTCCTTTT 48
******
['0']
0 0 G 1
******
['0']
0 49 A 50
******
['1']
1 0 GAAGAGGAAAAGTAGATCCCTGGCGTCCGGAGCTGGGACGT 41
******
['2']
2 0 C 1
******
['2']
2 48 C 49
******
請注意，只有第一個區塊是真正的比對結果；其他區塊只包含單一序列。藉由包含這些區塊，xmfa 檔案會包含 combined.fa 輸入檔案中提供的完整序列。
如果 progressiveMauve 是使用每個基因體個別的輸入檔案呼叫，例如：
progressiveMauve equCab1.fa canFam2.fa mm9.fa --output=separate.xmfa ...
其中每個 Fasta 檔案只包含一個物種的基因體序列，那麼輸出檔案 separate.xmfa 如下：
#FormatVersion Mauve1
#Sequence1File equCab1.fa
#Sequence1Format FastA
#Sequence2File canFam2.fa
#Sequence2Format FastA
#Sequence3File mm9.fa
#Sequence3Format FastA
#BackboneFile separate.xmfa.bbcols
> 1:1-50 - equCab1.fa
TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC
> 2:1-49 + canFam2.fa
CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC
> 3:1-19 - mm9.fa
---------------------------------GGATCTACTTTTCCTCTTC
=
> 3:20-41 + mm9.fa
CTGGCGTCCGGAGCTGGGACGT
=
馬的識別碼 equCab1、鼠的 mm9 和狗的 canFam2 現在會明確顯示在輸出檔案中。此 xmfa 檔案由兩個比對區塊組成，並以 = 字元分隔。使用 Align.parse 來提取這些比對結果。
>>> from Bio import Align
>>> alignments = Align.parse("separate.xmfa", "mauve")
檔案標頭資料現在不包含輸入檔案名稱。
>>> alignments.metadata
{'FormatVersion': 'Mauve1',
'BackboneFile': 'separate.xmfa.bbcols'}
identifiers 屬性儲存三個序列的序列識別碼。
>>> alignments.identifiers
['equCab1.fa', 'canFam2.fa', 'mm9.fa']
這些識別碼會用於個別的比對結果中。
>>> for alignment in alignments:
... print([record.id for record in alignment.sequences])
... print(alignment)
... print("******")
...
['equCab1.fa', 'canFam2.fa', 'mm9.fa']
equCab1.f 50 TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC 0
canFam2.f 0 CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC 49
mm9.fa 19 ---------------------------------GGATCTACTTTTCCTCTTC 0
******
['mm9.fa']
mm9.fa 19 CTGGCGTCCGGAGCTGGGACGT 41
******
要以 Mauve 格式輸出比對結果，請使用 Align.write。
>>> from io import StringIO
>>> stream = StringIO()
>>> alignments = Align.parse("separate.xmfa", "mauve")
>>> Align.write(alignments, stream, "mauve")
2
>>> print(stream.getvalue())
#FormatVersion Mauve1
#Sequence1File equCab1.fa
#Sequence1Format FastA
#Sequence2File canFam2.fa
#Sequence2Format FastA
#Sequence3File mm9.fa
#Sequence3Format FastA
#BackboneFile separate.xmfa.bbcols
> 1:1-50 - equCab1.fa
TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC
> 2:1-49 + canFam2.fa
CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC
> 3:1-19 - mm9.fa
---------------------------------GGATCTACTTTTCCTCTTC
=
> 3:20-41 + mm9.fa
CTGGCGTCCGGAGCTGGGACGT
=
在此，寫入器會使用儲存在 alignments.metadata 和 alignments.identifiers 中的資訊來建立此格式。如果您的 alignments 物件沒有這些屬性，您可以將它們作為關鍵字引數提供給 Align.write。
>>> stream = StringIO()
>>> alignments = Align.parse("separate.xmfa", "mauve")
>>> metadata = alignments.metadata
>>> identifiers = alignments.identifiers
>>> alignments = list(alignments) # this drops the attributes
>>> alignments.metadata
Traceback (most recent call last):
...
AttributeError: 'list' object has no attribute 'metadata'
>>> alignments.identifiers
Traceback (most recent call last):
...
AttributeError: 'list' object has no attribute 'identifiers'
>>> Align.write(alignments, stream, "mauve", metadata=metadata, identifiers=identifiers)
2
>>> print(stream.getvalue())
#FormatVersion Mauve1
#Sequence1File equCab1.fa
#Sequence1Format FastA
#Sequence2File canFam2.fa
#Sequence2Format FastA
#Sequence3File mm9.fa
#Sequence3Format FastA
#BackboneFile separate.xmfa.bbcols
> 1:1-50 - equCab1.fa
TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC
> 2:1-49 + canFam2.fa
CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC
> 3:1-19 - mm9.fa
---------------------------------GGATCTACTTTTCCTCTTC
=
> 3:20-41 + mm9.fa
CTGGCGTCCGGAGCTGGGACGT
=
Python 不允許您直接將這些屬性新增至 alignments 物件，如本範例所示，它已轉換為簡單的清單。但是，您可以建構 Alignments 物件並將屬性新增至其中（請參閱第 Alignments 類別 節）。
>>> alignments = Align.Alignments(alignments)
>>> alignments.metadata = metadata
>>> alignments.identifiers = identifiers
>>> stream = StringIO()
>>> Align.write(alignments, stream, "mauve", metadata=metadata, identifiers=identifiers)
2
>>> print(stream.getvalue())
#FormatVersion Mauve1
#Sequence1File equCab1.fa
#Sequence1Format FastA
#Sequence2File canFam2.fa
#Sequence2Format FastA
#Sequence3File mm9.fa
#Sequence3Format FastA
#BackboneFile separate.xmfa.bbcols
> 1:1-50 - equCab1.fa
TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC
> 2:1-49 + canFam2.fa
CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC
> 3:1-19 - mm9.fa
---------------------------------GGATCTACTTTTCCTCTTC
=
> 3:20-41 + mm9.fa
CTGGCGTCCGGAGCTGGGACGT
=
以 Mauve 格式列印單一比對結果時，請使用關鍵字引數來提供中繼資料和識別碼。
>>> alignment = alignments[0]
>>> print(alignment.format("mauve", metadata=metadata, identifiers=identifiers))
> 1:1-50 - equCab1.fa
TAAGCCCTCCTAGCACACACCCGGAGTGGCC-GGGCCGTAC-TTTCCTTTTC
> 2:1-49 + canFam2.fa
CAAGCCCTGC--GCGCTCAGCCGGAGTGTCCCGGGCCCTGC-TTTCCTTTTC
> 3:1-19 - mm9.fa
---------------------------------GGATCTACTTTTCCTCTTC
=
序列比對/地圖 (SAM)
序列比對/地圖 (SAM) 格式 [Li2009] 的檔案會儲存成對的序列比對結果，通常是新一代定序數據相對於參考基因體的比對結果。Biopython 測試套件中的 ex1.sam 檔案是 SAM 格式中最小檔案的範例。它的前幾行如下：
EAS56_57:6:190:289:82 69 chr1 100 0 * = 100 0 CTCAAGGTTGTTGCAAGGGGGTCTATGTGAACAAA <<<7<<<;<<<<<<<<8;;<7;4<;<;;;;;94<; MF:i:192
EAS56_57:6:190:289:82 137 chr1 100 73 35M = 100 0 AGGGGTGCAGAGCCGAGTCACGGGGTTGCCAGCAC <<<<<<;<<<<<<<<<<;<<;<<<<;8<6;9;;2; MF:i:64 Aq:i:0 NM:i:0 UQ:i:0 H0:i:1 H1:i:0
EAS51_64:3:190:727:308 99 chr1 103 99 35M = 263 195 GGTGCAGAGCCGAGTCACGGGGTTGCCAGCACAGG <<<<<<<<<<<<<<<<<<<<<<<<<<<::<<<844 MF:i:18 Aq:i:73 NM:i:0 UQ:i:0 H0:i:1 H1:i:0
...
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("ex1.sam", "sam")
>>> alignment = next(alignments)
第一行的 flag 值為 69。根據 SAM/BAM 檔案格式規範，flag 中包含位元旗標 4 的行表示未對應 (unmapped)。由於 69 對應此位置的位元設定為 True，因此此序列未對應，且未比對到基因組（儘管第一行顯示 chr1）。因此，此比對的目標（或 alignment.sequences 中的第一個項目）為 None。
>>> alignment.flag
69
>>> bin(69)
'0b1000101'
>>> bin(4)
'0b100'
>>> if alignment.flag & 4:
... print("unmapped")
... else:
... print("mapped")
...
unmapped
>>> alignment.sequences
[None, SeqRecord(seq=Seq('CTCAAGGTTGTTGCAAGGGGGTCTATGTGAACAAA'), id='EAS56_57:6:190:289:82', name='<unknown name>', description='', dbxrefs=[])]
>>> alignment.target is None
True
第二行代表比對到染色體 1 的結果。
>>> alignment = next(alignments)
>>> if alignment.flag & 4:
... print("unmapped")
... else:
... print("mapped")
...
mapped
>>> alignment.target
SeqRecord(seq=None, id='chr1', name='<unknown name>', description='', dbxrefs=[])
由於此 SAM 檔案未儲存每個比對的基因組序列資訊，因此我們無法列印比對結果。然而，我們可以以 SAM 格式或任何其他格式（例如 BED，請參閱 可延伸的瀏覽器資料 (BED) 章節），列印比對資訊，這些格式不需要目標序列資訊。
>>> format(alignment, "sam")
'EAS56_57:6:190:289:82\t137\tchr1\t100\t73\t35M\t=\t100\t0\tAGGGGTGCAGAGCCGAGTCACGGGGTTGCCAGCAC\t<<<<<<;<<<<<<<<<<;<<;<<<<;8<6;9;;2;\tMF:i:64\tAq:i:0\tNM:i:0\tUQ:i:0\tH0:i:1\tH1:i:0\n'
>>> format(alignment, "bed")
'chr1\t99\t134\tEAS56_57:6:190:289:82\t0\t+\t99\t134\t0\t1\t35,\t0,\n'
但是，我們無法以 PSL 格式列印比對結果（請參閱 模式空間佈局 (PSL) 章節），因為這需要知道目標序列 chr1 的大小。
>>> format(alignment, "psl")
Traceback (most recent call last):
...
TypeError: ...
如果您知道目標序列的大小，您可以手動設定它們。
>>> from Bio.Seq import Seq
>>> from Bio.SeqRecord import SeqRecord
>>> target = SeqRecord(Seq(None, length=1575), id="chr1")
>>> alignment.target = target
>>> format(alignment, "psl")
'35\t0\t0\t0\t0\t0\t0\t0\t+\tEAS56_57:6:190:289:82\t35\t0\t35\tchr1\t1575\t99\t134\t1\t35,\t0,\t99,\n'
Biopython 測試套件中的檔案 ex1_header.sam 包含相同的比對結果，但現在也包含標頭。其前幾行如下所示：
@HD\tVN:1.3\tSO:coordinate
@SQ\tSN:chr1\tLN:1575
@SQ\tSN:chr2\tLN:1584
EAS56_57:6:190:289:82 69 chr1 100 0 * = 100 0 CTCAAGGTTGTTGCAAGGGGGTCTATGTGAACAAA <<<7<<<;<<<<<<<<8;;<7;4<;<;;;;;94<; MF:i:192
...
標頭儲存關於比對結果的通用資訊，包括目標染色體的大小。目標資訊儲存在 alignments 物件的 targets 屬性中。
>>> from Bio import Align
>>> alignments = Align.parse("ex1_header.sam", "sam")
>>> len(alignments.targets)
2
>>> alignments.targets[0]
SeqRecord(seq=Seq(None, length=1575), id='chr1', name='<unknown name>', description='', dbxrefs=[])
>>> alignments.targets[1]
SeqRecord(seq=Seq(None, length=1584), id='chr2', name='<unknown name>', description='', dbxrefs=[])
標頭中提供的其他資訊儲存在 metadata 屬性中。
>>> alignments.metadata
{'HD': {'VN': '1.3', 'SO': 'coordinate'}}
有了目標資訊，我們現在也可以以 PSL 格式列印比對結果。
>>> alignment = next(alignments) # the unmapped sequence; skip it
>>> alignment = next(alignments)
>>> format(alignment, "psl")
'35\t0\t0\t0\t0\t0\t0\t0\t+\tEAS56_57:6:190:289:82\t35\t0\t35\tchr1\t1575\t99\t134\t1\t35,\t0,\t99,\n'
我們現在也可以以人類可讀的形式列印比對結果，但請注意，此檔案中無法取得目標序列內容。
>>> print(alignment)
chr1 99 ??????????????????????????????????? 134
0 ................................... 35
EAS56_57: 0 AGGGGTGCAGAGCCGAGTCACGGGGTTGCCAGCAC 35
Biopython 測試套件中 sam1.sam 檔案中的比對結果包含額外的 MD 標籤，該標籤顯示查詢序列與目標序列的差異。
@SQ SN:1 LN:239940
@PG ID:bwa PN:bwa VN:0.6.2-r126
HWI-1KL120:88:D0LRBACXX:1:1101:1780:2146 77 * 0 0 * * 0 0 GATGGGAAACCCATGGCCGAGTGGGAAGAAACCAGCTGAGGTCACATCACCAGAGGAGGGAGAGTGTGGCCCCTGACTCAGTCCATCAGCTTGTGGAGCTG @=?DDDDBFFFF7A;E?GGEGE8BB?FF?F>G@F=GIIDEIBCFF<FEFEC@EEEE2?8B8/=@((-;?@2<B9@##########################
...
HWI-1KL120:88:D0LRBACXX:1:1101:2852:2134 137 1 136186 25 101M = 136186 0 TCACGGTGGCCTGTTGAGGCAGGGGCTCACGCTGACCTCTCTCGGCGTGGGAGGGGCCGGTGTGAGGCAAGGGCTCACGCTGACCTCTCTCGGCGTGGGAG @C@FFFDFHGHHHJJJIJJJJIJJJGEDHHGGHGBGIIGIIAB@GEE=BDBBCCDD@D@B7@;@DDD?<A?DD728:>8()009>:>>C@>5??B###### XT:A:U NM:i:5 SM:i:25 AM:i:0 X0:i:1 X1:i:0 XM:i:5 XO:i:0 XG:i:0 MD:Z:25G14G2C34A12A9
解析器從 MD 標籤重建本地基因組序列，讓我們在列印比對結果時可以明確地看到目標序列。
>>> from Bio import Align
>>> alignments = Align.parse("sam1.sam", "sam")
>>> for alignment in alignments:
... if not alignment.flag & 4: # Skip the unmapped lines
... break
...
>>> alignment
<Alignment object (2 rows x 101 columns) at ...>
>>> print(alignment)
1 136185 TCACGGTGGCCTGTTGAGGCAGGGGGTCACGCTGACCTCTGTCCGCGTGGGAGGGGCCGG
0 |||||||||||||||||||||||||.||||||||||||||.||.||||||||||||||||
HWI-1KL12 0 TCACGGTGGCCTGTTGAGGCAGGGGCTCACGCTGACCTCTCTCGGCGTGGGAGGGGCCGG
1 136245 TGTGAGGCAAGGGCTCACACTGACCTCTCTCAGCGTGGGAG 136286
60 ||||||||||||||||||.||||||||||||.||||||||| 101
HWI-1KL12 60 TGTGAGGCAAGGGCTCACGCTGACCTCTCTCGGCGTGGGAG 101
SAM 檔案可能包含額外資訊，以區分簡單的序列插入和刪除與基因組的跳過區域（例如，內含子）、硬剪輯和軟剪輯以及填充序列區域。由於此資訊無法儲存在 Alignment 物件的 coordinates 屬性中，而是儲存在專用的 operations 屬性中。讓我們以這個 SAM 檔案中的第三個比對為例。
>>> from Bio import Align
>>> alignments = Align.parse("dna_rna.sam", "sam")
>>> alignment = next(alignments)
>>> alignment = next(alignments)
>>> alignment = next(alignments)
>>> print(format(alignment, "SAM"))
NR_111921.1 0 chr3 48663768 0 46M1827N82M3376N76M12H * 0 0 CACGAGAGGAGCGGAGGCGAGGGGTGAACGCGGAGCACTCCAATCGCTCCCAACTAGAGGTCCACCCAGGACCCAGAGACCTGGATTTGAGGCTGCTGGGCGGCAGATGGAGCGATCAGAAGACCAGGAGACGGGAGCTGGAGTGCAGTGGCTGTTCACAAGCGTGAAAGCAAAGATTAAAAAATTTGTTTTTATATTAAAAAA * AS:i:1000 NM:i:0
>>> print(alignment.coordinates)
[[48663767 48663813 48665640 48665722 48669098 48669174]
[ 0 46 46 128 128 204]]
>>> alignment.operations
bytearray(b'MNMNM')
>>> alignment.query.annotations["hard_clip_right"]
12
在這個比對中，cigar 字串 63M1062N75M468N43M 定義了 46 個對齊的核苷酸、一個 1827 個核苷酸的內含子、82 個對齊的核苷酸、一個 3376 個核苷酸的內含子、76 個對齊的核苷酸和 12 個硬剪輯的核苷酸。這些操作顯示在 operations 屬性中，硬剪輯除外，它儲存在 alignment.query.annotations["hard_clip_right"] 中（或如果適用，則儲存在 alignment.query.annotations["hard_clip_left"] 中）。
若要寫入從頭建立的比對的 SAM 檔案，請使用 Alignments（複數）物件（請參閱 Alignments 類別 章節）來儲存比對結果以及中繼資料和目標。
>>> from io import StringIO
>>> import numpy as np
>>> from Bio import Align
>>> from Bio.Seq import Seq
>>> from Bio.SeqRecord import SeqRecord
>>> alignments = Align.Alignments()
>>> seq1 = Seq(None, length=10000)
>>> target1 = SeqRecord(seq1, id="chr1")
>>> seq2 = Seq(None, length=15000)
>>> target2 = SeqRecord(seq2, id="chr2")
>>> alignments.targets = [target1, target2]
>>> alignments.metadata = {"HD": {"VN": "1.3", "SO": "coordinate"}}
>>> seqA = Seq(None, length=20)
>>> queryA = SeqRecord(seqA, id="readA")
>>> sequences = [target1, queryA]
>>> coordinates = np.array([[4300, 4320], [0, 20]])
>>> alignment = Align.Alignment(sequences, coordinates)
>>> alignments.append(alignment)
>>> seqB = Seq(None, length=25)
>>> queryB = SeqRecord(seqB, id="readB")
>>> sequences = [target1, queryB]
>>> coordinates = np.array([[5900, 5925], [25, 0]])
>>> alignment = Align.Alignment(sequences, coordinates)
>>> alignments.append(alignment)
>>> seqC = Seq(None, length=40)
>>> queryC = SeqRecord(seqC, id="readC")
>>> sequences = [target2, queryC]
>>> coordinates = np.array([[12300, 12318], [0, 18]])
>>> alignment = Align.Alignment(sequences, coordinates)
>>> alignments.append(alignment)
>>> stream = StringIO()
>>> Align.write(alignments, stream, "sam")
3
>>> print(stream.getvalue())
@HD VN:1.3 SO:coordinate
@SQ SN:chr1 LN:10000
@SQ SN:chr2 LN:15000
readA 0 chr1 4301 255 20M * 0 0 * *
readB 16 chr1 5901 255 25M * 0 0 * *
readC 0 chr2 12301 255 18M22S * 0 0 * *
可延伸的瀏覽器資料 (BED)
BED（可延伸的瀏覽器資料）檔案通常用於儲存基因轉錄本與基因組的比對結果。請參閱 UCSC 的說明，以取得 BED 格式的完整說明。
BED 檔案有三個必需欄位和九個可選欄位。子目錄 Tests/Blat 中的檔案 bed12.bed 是具有 12 個欄位的 BED 檔案範例。
chr22 1000 5000 mRNA1 960 + 1200 4900 255,0,0 2 567,488, 0,3512,
chr22 2000 6000 mRNA2 900 - 2300 5960 0,255,0 2 433,399, 0,3601,
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("bed12.bed", "bed")
>>> len(alignments)
2
>>> for alignment in alignments:
... print(alignment.coordinates)
...
[[1000 1567 4512 5000]
[ 0 567 567 1055]]
[[2000 2433 5601 6000]
[ 832 399 399 0]]
請注意，第一個序列（「mRNA1」）對應到正向鏈，而第二個序列（「mRNA2」）對應到反向鏈。
由於 BED 檔案不儲存每個染色體的長度，因此目標序列的長度會設定為其最大值。
>>> alignment.target
SeqRecord(seq=Seq(None, length=9223372036854775807), id='chr22', name='<unknown name>', description='', dbxrefs=[])
查詢序列的長度可以從其比對資訊推斷出來。
>>> alignment.query
SeqRecord(seq=Seq(None, length=832), id='mRNA2', name='<unknown name>', description='', dbxrefs=[])
比對分數（欄位 5）和欄位 7-9 中儲存的資訊（在 BED 格式規範中稱為 thickStart、thickEnd 和 itemRgb）會儲存為 alignment 物件的屬性。
>>> alignment.score
900.0
>>> alignment.thickStart
2300
>>> alignment.thickEnd
5960
>>> alignment.itemRgb
'0,255,0'
若要以 BED 格式列印比對結果，您可以使用 Python 的內建 format 函式。
>>> print(format(alignment, "bed"))
chr22 2000 6000 mRNA2 900 - 2300 5960 0,255,0 2 433,399, 0,3601,
或者，您可以使用 alignment 物件的 format 方法。這可讓您指定要寫入的欄位數，作為 bedN 關鍵字引數。
>>> print(alignment.format("bed"))
chr22 2000 6000 mRNA2 900 - 2300 5960 0,255,0 2 433,399, 0,3601,
>>> print(alignment.format("bed", 3))
chr22 2000 6000
>>> print(alignment.format("bed", 6))
chr22 2000 6000 mRNA2 900 -
相同的關鍵字引數可以與 Align.write 一起使用。
>>> Align.write(alignments, "mybed3file.bed", "bed", bedN=3)
2
>>> Align.write(alignments, "mybed6file.bed", "bed", bedN=6)
2
>>> Align.write(alignments, "mybed12file.bed", "bed")
2
bigBed
bigBed 檔案格式是 BED 檔案的索引二進位版本 可延伸的瀏覽器資料 (BED)。若要建立 bigBed 檔案，您可以從 UCSC 使用 bedToBigBed 程式 (`) <https://genome.ucsc.edu/goldenPath/help/bigBed.html>`__。或者，您可以透過呼叫 Bio.Align.write 函式並使用 fmt="bigbed" 來使用 Biopython 進行此操作。雖然這兩種方法應該產生相同的 bigBed 檔案，但使用 bedToBigBed 速度更快，而且可能更可靠，因為它是黃金標準。由於 bigBed 檔案具有內建索引，因此可讓您快速搜尋特定的基因組區域。
例如，讓我們解析 bigBed 檔案 dna_rna.bb，該檔案位於 Biopython 發行版本的 Tests/Blat 子目錄中。
>>> from Bio import Align
>>> alignments = Align.parse("dna_rna.bb", "bigbed")
>>> len(alignments)
4
>>> print(alignments.declaration)
table bed
"Browser Extensible Data"
(
string chrom; "Reference sequence chromosome or scaffold"
uint chromStart; "Start position in chromosome"
uint chromEnd; "End position in chromosome"
string name; "Name of item."
uint score; "Score (0-1000)"
char[1] strand; "+ or - for strand"
uint thickStart; "Start of where display should be thick (start codon)"
uint thickEnd; "End of where display should be thick (stop codon)"
uint reserved; "Used as itemRgb as of 2004-11-22"
int blockCount; "Number of blocks"
int[blockCount] blockSizes; "Comma separated list of block sizes"
int[blockCount] chromStarts; "Start positions relative to chromStart"
)
declaration 包含用於建立 bigBed 檔案的欄位規格，以 AutoSql 格式表示。目標序列（通常是對齊序列的染色體）儲存在 targets 屬性中。在 bigBed 格式中，只會儲存每個目標的識別碼和大小。在此範例中，只有一個染色體。
>>> alignments.targets
[SeqRecord(seq=Seq(None, length=198295559), id='chr3', name='<unknown name>', description='<unknown description>', dbxrefs=[])]
讓我們看看個別的比對。比對資訊的儲存方式與 BED 檔案相同（請參閱 可延伸的瀏覽器資料 (BED) 章節）。
>>> alignment = next(alignments)
>>> alignment.target.id
'chr3'
>>> alignment.query.id
'NR_046654.1'
>>> alignment.coordinates
array([[42530895, 42530958, 42532020, 42532095, 42532563, 42532606],
[ 181, 118, 118, 43, 43, 0]])
>>> alignment.thickStart
42530895
>>> alignment.thickEnd
42532606
>>> print(alignment)
chr3 42530895 ????????????????????????????????????????????????????????????
0 ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
NR_046654 181 ????????????????????????????????????????????????????????????
chr3 42530955 ????????????????????????????????????????????????????????????
60 |||---------------------------------------------------------
NR_046654 121 ???---------------------------------------------------------
...
chr3 42532515 ????????????????????????????????????????????????????????????
1620 ------------------------------------------------||||||||||||
NR_046654 43 ------------------------------------------------????????????
chr3 42532575 ??????????????????????????????? 42532606
1680 ||||||||||||||||||||||||||||||| 1711
NR_046654 31 ??????????????????????????????? 0
預設的 bigBed 格式不會儲存目標和查詢的序列內容。如果這些內容在其他地方可用（例如，Fasta 檔案），您可以設定 alignment.target.seq 和 alignment.query.seq，以在列印比對結果時顯示序列內容，或以需要序列內容的格式（例如 Clustal，請參閱 ClustalW 章節）寫入比對結果。Biopython 發行版本的 Tests 子目錄中的測試腳本 test_Align_bigbed.py 提供了一些關於如何執行此操作的範例。
現在讓我們看看如何搜尋序列區域。這些是以 BED 格式列印的 bigBed 檔案中儲存的序列（請參閱 可延伸的瀏覽器資料 (BED) 章節）。
>>> for alignment in alignments:
... print(format(alignment, "bed"))
...
chr3 42530895 42532606 NR_046654.1 1000 - 42530895 42532606 0 3 63,75,43, 0,1125,1668,
chr3 42530895 42532606 NR_046654.1_modified 978 - 42530895 42532606 0 5 27,36,17,56,43, 0,27,1125,1144,1668,
chr3 48663767 48669174 NR_111921.1 1000 + 48663767 48669174 0 3 46,82,76, 0,1873,5331,
chr3 48663767 48669174 NR_111921.1_modified 972 + 48663767 48669174 0 5 28,17,76,6,76, 0,29,1873,1949,5331,
使用 alignments 物件上的 search 方法，在 chr3 上找到位置 48000000 和 49000000 之間的區域。此方法會傳回一個迭代器。
>>> selected_alignments = alignments.search("chr3", 48000000, 49000000)
>>> for alignment in selected_alignments:
... print(alignment.query.id)
...
NR_111921.1
NR_111921.1_modified
染色體名稱可以是 None 以包含所有染色體，開始和結束位置可以是 None 以分別從位置 0 開始搜尋，或繼續搜尋到染色體的末尾。
以 bigBed 格式寫入比對結果與呼叫 Bio.Align.write 一樣簡單。
>>> Align.write(alignments, "output.bb", "bigbed")
您可以指定要包含在 bigBed 檔案中的 BED 欄位數。例如，若要寫入 BED6 檔案，請使用：
>>> Align.write(alignments, "output.bb", "bigbed", bedN=6)
與 bedToBigBed 相同，您可以在 bigBed 輸出中包含額外欄位。假設檔案 bedExample2.as（在 Biopython 發行版本的 Tests/Blat 子目錄中可用）以 AutoSql 格式儲存所包含 BED 欄位的宣告。我們可以讀取此宣告如下：
>>> from Bio.Align import bigbed
>>> with open("bedExample2.as") as stream:
... autosql_data = stream.read()
...
>>> declaration = bigbed.AutoSQLTable.from_string(autosql_data)
>>> type(declaration)
<class 'Bio.Align.bigbed.AutoSQLTable'>
>>> print(declaration)
table hg18KGchr7
"UCSC Genes for chr7 with color plus GeneSymbol and SwissProtID"
(
string chrom; "Reference sequence chromosome or scaffold"
uint chromStart; "Start position of feature on chromosome"
uint chromEnd; "End position of feature on chromosome"
string name; "Name of gene"
uint score; "Score"
char[1] strand; "+ or - for strand"
uint thickStart; "Coding region start"
uint thickEnd; "Coding region end"
uint reserved; "Green on + strand, Red on - strand"
string geneSymbol; "Gene Symbol"
string spID; "SWISS-PROT protein Accession number"
)
現在，我們可以透過呼叫以下程式碼，寫入包含 9 個 BED 欄位加上額外欄位 geneSymbol 和 spID 的 bigBed 檔案：
>>> Align.write(
... alignments,
... "output.bb",
... "bigbed",
... bedN=9,
... declaration=declaration,
... extraIndex=["name", "geneSymbol"],
... )
在此，我們也要求在大Bed檔案中，針對 name 和 geneSymbol 增加額外的索引。Align.write 預期在 alignment.annotations 字典中找到鍵 geneSymbol 和 spID。請參考 Biopython 發行版中 Tests 子目錄下的測試腳本 test_Align_bigbed.py，以取得更多以 bigBed 格式寫入比對檔案的範例。
可選參數包括 compress (預設值為 True)、blockSize (預設值為 256) 和 itemsPerSlot (預設值為 512)。請參閱 UCSC 的 bedToBigBed 程式的說明文件，以了解這些參數的描述。在建立 bigBed 檔案時，使用 compress=False 和 itemsPerSlot=1 可以加快搜尋 bigBed 檔案的速度。
模式空間佈局 (PSL)
PSL（模式空間佈局）檔案由 BLAST-Like Alignment Tool BLAT [Kent2002] 產生。與 BED 檔案（請參閱可擴展瀏覽器資料 (BED)章節）類似，PSL 檔案通常用於儲存轉錄本與基因組的比對。以下是一個簡短 BLAT 檔案的範例（在 Biopython 發行版的 Tests/Blat 子目錄中，以 dna_rna.psl 的名稱提供），其中包含標準的 5 行 PSL 標頭
psLayout version 3
match mis- rep. N's Q gap Q gap T gap T gap strand Q Q Q Q T T T T block blockSizes qStarts tStarts
match match count bases count bases name size start end name size start end count
---------------------------------------------------------------------------------------------------------------------------------------------------------------
165 0 39 0 0 0 2 5203 + NR_111921.1 216 0 204 chr3 198295559 48663767 48669174 3 46,82,76, 0,46,128, 48663767,48665640,48669098,
175 0 6 0 0 0 2 1530 - NR_046654.1 181 0 181 chr3 198295559 42530895 42532606 3 63,75,43, 0,63,138, 42530895,42532020,42532563,
162 2 39 0 1 2 3 5204 + NR_111921.1_modified 220 3 208 chr3 198295559 48663767 48669174 5 28,17,76,6,76, 3,31,48,126,132, 48663767,48663796,48665640,48665716,48669098,
172 1 6 0 1 3 3 1532 - NR_046654.1_modified 190 3 185 chr3 198295559 42530895 42532606 5 27,36,17,56,43, 5,35,71,88,144, 42530895,42530922,42532020,42532039,42532563,
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("dna_rna.psl", "psl")
>>> alignments.metadata
{'psLayout version': '3'}
迭代處理比對，以取得每行一個 Alignment 物件
>>> for alignment in alignments:
... print(alignment.target.id, alignment.query.id)
...
chr3 NR_046654.1
chr3 NR_046654.1_modified
chr3 NR_111921.1
chr3 NR_111921.1_modified
讓我們更詳細地查看最後一個比對。PSL 檔案中的前四欄顯示匹配數、不匹配數、與重複區域比對的核苷酸數，以及與 N（未知）字元比對的核苷酸數。這些值儲存在 Alignment 物件的屬性中
>>> alignment.matches
162
>>> alignment.misMatches
2
>>> alignment.repMatches
39
>>> alignment.nCount
0
由於目標和查詢的序列資料沒有明確儲存在 PSL 檔案中，因此 alignment.target 和 alignment.query 的序列內容未定義。然而，它們的序列長度是已知的
>>> alignment.target
SeqRecord(seq=Seq(None, length=198295559), id='chr3', ...)
>>> alignment.query
SeqRecord(seq=Seq(None, length=220), id='NR_111921.1_modified', ...)
我們可以以 BED 或 PSL 格式列印比對結果
>>> print(format(alignment, "bed"))
chr3 48663767 48669174 NR_111921.1_modified 0 + 48663767 48669174 0 5 28,17,76,6,76, 0,29,1873,1949,5331,
>>> print(format(alignment, "psl"))
162 2 39 0 1 2 3 5204 + NR_111921.1_modified 220 3 208 chr3 198295559 48663767 48669174 5 28,17,76,6,76, 3,31,48,126,132, 48663767,48663796,48665640,48665716,48669098,
在此，匹配數、不匹配數、重複區域匹配數和與未知核苷酸的匹配數取自 Alignment 物件的對應屬性。如果這些屬性不可用，例如比對結果不是來自 PSL 檔案，則會使用序列內容（如果可用）計算這些數字。目標序列中的重複區域會透過遮蔽序列為小寫或大寫字元來指示，這由 mask 關鍵字引數的下列值定義
False(預設值)：不分開計算與遮蔽序列的匹配；"lower"：計算並報告與小寫字元的匹配為與重複區域的匹配；"upper"：計算並報告與大寫字元的匹配為與重複區域的匹配；
用於未知核苷酸的字元由 wildcard 引數定義。為了與 BLAT 一致，萬用字元預設為 "N"。如果您不想單獨計算與任何未知核苷酸的匹配，請使用 wildcard=None。
>>> import numpy
>>> from Bio import Align
>>> query = "GGTGGGGG"
>>> target = "AAAAAAAggggGGNGAAAAA"
>>> coordinates = numpy.array([[0, 7, 15, 20], [0, 0, 8, 8]])
>>> alignment = Align.Alignment([target, query], coordinates)
>>> print(alignment)
target 0 AAAAAAAggggGGNGAAAAA 20
0 -------....||.|----- 20
query 0 -------GGTGGGGG----- 8
>>> line = alignment.format("psl")
>>> print(line)
6 1 0 1 0 0 0 0 + query 8 0 8 target 20 7 15 1 8, 0, 7,
>>> line = alignment.format("psl", mask="lower")
>>> print(line)
3 1 3 1 0 0 0 0 + query 8 0 8 target 20 7 15 1 8, 0, 7,
>>> line = alignment.format("psl", mask="lower", wildcard=None)
>>> print(line)
3 2 3 0 0 0 0 0 + query 8 0 8 target 20 7 15 1 8, 0, 7,
當使用 Bio.Align.write 將比對結果寫入 PSL 格式的輸出檔案時，可以使用相同的引數。此函數有一個額外的關鍵字 header（預設值為 True），用於指定是否應寫入 PSL 標頭。
除了 format 方法之外，您可以使用 Python 的內建 format 函數
>>> print(format(alignment, "psl"))
6 1 0 1 0 0 0 0 + query 8 0 8 target 20 7 15 1 8, 0, 7,
允許在 Python 中使用格式化的 (f-) 字串中使用 Alignment 物件
>>> line = f"The alignment in PSL format is '{alignment:psl}'."
>>> print(line)
The alignment in PSL format is '6 1 0 1 0 0 0 0 + query 8 0 8 target 20 7 15 1 8, 0, 7,
'
請注意，可選的關鍵字引數不能與 format 函數或格式化的字串一起使用。
bigPsl
bigPsl 檔案是一個 bigBed 檔案，具有 BED12+13 格式，由 12 個預定義的 BED 欄位和 UCSC 提供的 AutoSql 檔案 bigPsl.as 中定義的 13 個自訂欄位組成，產生 PSL 檔案的索引二進位版本（請參閱模式空間佈局 (PSL)章節）。若要建立 bigPsl 檔案，您可以使用 UCSC 的 pslToBigPsl 和 bedToBigBed 程式，也可以透過呼叫 Bio.Align.write 函數並使用 fmt="bigpsl" 來使用 Biopython。雖然這兩種方法應該會產生相同的 bigPsl 檔案，但 UCSC 工具的速度快得多，而且可能更可靠，因為它是黃金標準。由於 bigPsl 檔案是 bigBed 檔案，因此它們具有內建索引，允許您快速搜尋特定的基因組區域。
例如，讓我們剖析 Biopython 發行版的 Tests/Blat 子目錄中提供的 bigBed 檔案 dna_rna.psl.bb。此檔案是 bigPsl 檔案，相當於 bigBed 檔案 dna_rna.bb（請參閱bigBed 章節）和 PSL 檔案 dna_rna.psl（請參閱模式空間佈局 (PSL) 章節）。
>>> from Bio import Align
>>> alignments = Align.parse("dna_rna.psl.bb", "bigpsl")
>>> len(alignments)
4
>>> print(alignments.declaration)
table bigPsl
"bigPsl pairwise alignment"
(
string chrom; "Reference sequence chromosome or scaffold"
uint chromStart; "Start position in chromosome"
uint chromEnd; "End position in chromosome"
string name; "Name or ID of item, ideally both human readable and unique"
uint score; "Score (0-1000)"
char[1] strand; "+ or - indicates whether the query aligns to the + or - strand on the reference"
uint thickStart; "Start of where display should be thick (start codon)"
uint thickEnd; "End of where display should be thick (stop codon)"
uint reserved; "RGB value (use R,G,B string in input file)"
int blockCount; "Number of blocks"
int[blockCount] blockSizes; "Comma separated list of block sizes"
int[blockCount] chromStarts; "Start positions relative to chromStart"
uint oChromStart; "Start position in other chromosome"
uint oChromEnd; "End position in other chromosome"
char[1] oStrand; "+ or -, - means that psl was reversed into BED-compatible coordinates"
uint oChromSize; "Size of other chromosome."
int[blockCount] oChromStarts; "Start positions relative to oChromStart or from oChromStart+oChromSize depending on strand"
lstring oSequence; "Sequence on other chrom (or edit list, or empty)"
string oCDS; "CDS in NCBI format"
uint chromSize; "Size of target chromosome"
uint match; "Number of bases matched."
uint misMatch; "Number of bases that don't match"
uint repMatch; "Number of bases that match but are part of repeats"
uint nCount; "Number of 'N' bases"
uint seqType; "0=empty, 1=nucleotide, 2=amino_acid"
)
宣告包含 UCSC 的 bigPsl.as AutoSql 檔案定義的欄位規格。目標序列（通常是序列比對的染色體）儲存在 targets 屬性中。在 bigBed 格式中，僅儲存每個目標的識別碼和大小。在此範例中，只有一個染色體
>>> alignments.targets
[SeqRecord(seq=Seq(None, length=198295559), id='chr3', name='<unknown name>', description='<unknown description>', dbxrefs=[])]
迭代處理比對會為每一行提供一個 Alignment 物件
>>> for alignment in alignments:
... print(alignment.target.id, alignment.query.id)
...
chr3 NR_046654.1
chr3 NR_046654.1_modified
chr3 NR_111921.1
chr3 NR_111921.1_modified
讓我們看看各個比對。比對資訊的儲存方式與對應的 PSL 檔案相同（請參閱模式空間佈局 (PSL) 章節）
>>> alignment.coordinates
array([[48663767, 48663795, 48663796, 48663813, 48665640, 48665716,
48665716, 48665722, 48669098, 48669174],
[ 3, 31, 31, 48, 48, 124,
126, 132, 132, 208]])
>>> alignment.thickStart
48663767
>>> alignment.thickEnd
48669174
>>> alignment.matches
162
>>> alignment.misMatches
2
>>> alignment.repMatches
39
>>> alignment.nCount
0
我們可以以 BED 或 PSL 格式列印比對結果
>>> print(format(alignment, "bed"))
chr3 48663767 48669174 NR_111921.1_modified 1000 + 48663767 48669174 0 5 28,17,76,6,76, 0,29,1873,1949,5331,
>>> print(format(alignment, "psl"))
162 2 39 0 1 2 3 5204 + NR_111921.1_modified 220 3 208 chr3 198295559 48663767 48669174 5 28,17,76,6,76, 3,31,48,126,132, 48663767,48663796,48665640,48665716,48669098,
由於 bigPsl 檔案是 bigBed 檔案的特殊情況，您可以使用比對物件的 search 方法來尋找與特定基因組區域的比對。例如，我們可以尋找 chr3 上位置 48000000 和 49000000 之間的區域
>>> selected_alignments = alignments.search("chr3", 48000000, 49000000)
>>> for alignment in selected_alignments:
... print(alignment.query.id)
...
NR_111921.1
NR_111921.1_modified
染色體名稱可以是 None 以包含所有染色體，開始和結束位置可以是 None 以分別從位置 0 開始搜尋，或繼續搜尋到染色體的末尾。
若要使用 Biopython 寫入 bigPsl 檔案，請使用 Bio.Align.write(alignments, "myfilename.bb", fmt="bigpsl")，其中 myfilename.bb 是輸出 bigPsl 檔案的名稱。或者，您可以使用（二進位）串流進行輸出。其他選項包括
compress：如果為True（預設值），則使用 zlib 壓縮資料；如果為False，則不壓縮資料。extraIndex：要編製索引的額外欄位名稱字串清單。cds：如果為True，則尋找類型為 CDS 的查詢特徵，並以 NCBI 樣式將其寫入 PSL 檔案（預設值：False）。fa：如果為True，則在 PSL 檔案中包含查詢序列（預設值：False）。mask：指定是否遮蔽目標序列中的重複區域，並且應在repMatches欄位中報告，而不是在matches欄位中報告。可接受的值包括None：不遮蔽（預設值）；"lower"：以小寫字元遮蔽；"upper"：以大寫字元遮蔽。
wildcard：在nCount欄位中報告與目標或查詢序列中萬用字元（代表未知核苷酸）的比對，而不是在matches、misMatches或repMatches欄位中報告。預設值為"N"。
請參閱模式空間佈局 (PSL)章節，以瞭解如何取得匹配數、不匹配數、重複區域匹配數和與未知核苷酸的匹配數。
其他可選參數包括 blockSize (預設值為 256) 和 itemsPerSlot (預設值為 512)。請參閱 UCSC 的 bedToBigBed 程式的說明文件，以了解這些參數的描述。在建立 bigPsl 檔案時，使用 compress=False 和 itemsPerSlot=1 可以加快搜尋 bigPsl 檔案的速度。
多重比對格式 (MAF)
MAF（多重比對格式）檔案以人類可讀的格式儲存一系列多重序列比對。MAF 檔案通常用於儲存基因組彼此的比對。Tests/MAF 子目錄中的檔案 ucsc_test.maf 是簡單 MAF 檔案的範例
track name=euArc visibility=pack mafDot=off frames="multiz28wayFrames" speciesOrder="hg16 panTro1 baboon mm4 rn3" description="A sample alignment"
##maf version=1 scoring=tba.v8
# tba.v8 (((human chimp) baboon) (mouse rat))
# multiz.v7
# maf_project.v5 _tba_right.maf3 mouse _tba_C
# single_cov2.v4 single_cov2 /dev/stdin
a score=23262.0
s hg16.chr7 27578828 38 + 158545518 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
s panTro1.chr6 28741140 38 + 161576975 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
s baboon 116834 38 + 4622798 AAA-GGGAATGTTAACCAAATGA---GTTGTCTCTTATGGTG
s mm4.chr6 53215344 38 + 151104725 -AATGGGAATGTTAAGCAAACGA---ATTGTCTCTCAGTGTG
s rn3.chr4 81344243 40 + 187371129 -AA-GGGGATGCTAAGCCAATGAGTTGTTGTCTCTCAATGTG
a score=5062.0
s hg16.chr7 27699739 6 + 158545518 TAAAGA
s panTro1.chr6 28862317 6 + 161576975 TAAAGA
s baboon 241163 6 + 4622798 TAAAGA
s mm4.chr6 53303881 6 + 151104725 TAAAGA
s rn3.chr4 81444246 6 + 187371129 taagga
a score=6636.0
s hg16.chr7 27707221 13 + 158545518 gcagctgaaaaca
s panTro1.chr6 28869787 13 + 161576975 gcagctgaaaaca
s baboon 249182 13 + 4622798 gcagctgaaaaca
s mm4.chr6 53310102 13 + 151104725 ACAGCTGAAAATA
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("ucsc_test.maf", "maf")
檔案標頭中顯示的資訊（track 行和後續以「#」開頭的行）會儲存在 alignments 物件的 metadata 屬性中
>>> alignments.metadata
{'name': 'euArc',
'visibility': 'pack',
'mafDot': 'off',
'frames': 'multiz28wayFrames',
'speciesOrder': ['hg16', 'panTro1', 'baboon', 'mm4', 'rn3'],
'description': 'A sample alignment',
'MAF Version': '1',
'Scoring': 'tba.v8',
'Comments': ['tba.v8 (((human chimp) baboon) (mouse rat))',
'multiz.v7',
'maf_project.v5 _tba_right.maf3 mouse _tba_C',
'single_cov2.v4 single_cov2 /dev/stdin']}
透過迭代處理 alignments，我們可以為 MAF 檔案中的每個比對區塊取得一個 Alignment 物件
>>> alignment = next(alignments)
>>> alignment.score
23262.0
>>> {seq.id: len(seq) for seq in alignment.sequences}
{'hg16.chr7': 158545518,
'panTro1.chr6': 161576975,
'baboon': 4622798,
'mm4.chr6': 151104725,
'rn3.chr4': 187371129}
>>> print(alignment.coordinates)
[[27578828 27578829 27578831 27578831 27578850 27578850 27578866]
[28741140 28741141 28741143 28741143 28741162 28741162 28741178]
[ 116834 116835 116837 116837 116856 116856 116872]
[53215344 53215344 53215346 53215347 53215366 53215366 53215382]
[81344243 81344243 81344245 81344245 81344264 81344267 81344283]]
>>> print(alignment)
hg16.chr7 27578828 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG 27578866
panTro1.c 28741140 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG 28741178
baboon 116834 AAA-GGGAATGTTAACCAAATGA---GTTGTCTCTTATGGTG 116872
mm4.chr6 53215344 -AATGGGAATGTTAAGCAAACGA---ATTGTCTCTCAGTGTG 53215382
rn3.chr4 81344243 -AA-GGGGATGCTAAGCCAATGAGTTGTTGTCTCTCAATGTG 81344283
>>> print(format(alignment, "phylip"))
5 42
hg16.chr7 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
panTro1.chAAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
baboon AAA-GGGAATGTTAACCAAATGA---GTTGTCTCTTATGGTG
mm4.chr6 -AATGGGAATGTTAAGCAAACGA---ATTGTCTCTCAGTGTG
rn3.chr4 -AA-GGGGATGCTAAGCCAATGAGTTGTTGTCTCTCAATGTG
除了 “a” (對齊區塊) 和 “s” (序列) 行之外，MAF 檔案可能包含 “i” 行，其中包含此區塊之前和之後的基因組序列資訊；“e” 行，其中包含對齊中空白部分的資訊；以及 “q” 行，顯示每個對齊鹼基的品質。這是一個包含此類行的對齊區塊範例
a score=19159.000000
s mm9.chr10 3014644 45 + 129993255 CCTGTACC---CTTTGGTGAGAATTTTTGTTTCAGTGTTAAAAGTTTG
s hg18.chr6 15870786 46 - 170899992 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
i hg18.chr6 I 9085 C 0
s panTro2.chr6 16389355 46 - 173908612 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
q panTro2.chr6 99999999999999999999999-9999999999999999999-9999
i panTro2.chr6 I 9106 C 0
s calJac1.Contig6394 6182 46 + 133105 CCTATACCTTTCTTTCATGAGAA-TTTTGTTTGAATCCTAAAC-TTTT
i calJac1.Contig6394 N 0 C 0
s loxAfr1.scaffold_75566 1167 34 - 10574 ------------TTTGGTTAGAA-TTATGCTTTAATTCAAAAC-TTCC
q loxAfr1.scaffold_75566 ------------99999699899-9999999999999869998-9997
i loxAfr1.scaffold_75566 N 0 C 0
e tupBel1.scaffold_114895.1-498454 167376 4145 - 498454 I
e echTel1.scaffold_288249 87661 7564 + 100002 I
e otoGar1.scaffold_334.1-359464 181217 2931 - 359464 I
e ponAbe2.chr6 16161448 8044 - 174210431 I
這是檔案 ucsc_mm9_chr10.maf 中的第 10 個對齊區塊 (位於 Biopython 發行版的 Tests/MAF 子目錄中)
>>> from Bio import Align
>>> alignments = Align.parse("ucsc_mm9_chr10.maf", "maf")
>>> for i in range(10):
... alignment = next(alignments)
...
>>> alignment.score
19159.0
>>> print(alignment)
mm9.chr10 3014644 CCTGTACC---CTTTGGTGAGAATTTTTGTTTCAGTGTTAAAAGTTTG 3014689
hg18.chr6 155029206 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT 155029160
panTro2.c 157519257 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT 157519211
calJac1.C 6182 CCTATACCTTTCTTTCATGAGAA-TTTTGTTTGAATCCTAAAC-TTTT 6228
loxAfr1.s 9407 ------------TTTGGTTAGAA-TTATGCTTTAATTCAAAAC-TTCC 9373
“i” 行顯示目前對齊區塊中的序列與前一個和後續對齊區塊中的序列之間的關係。此資訊儲存在相應序列的 annotations 屬性中
>>> alignment.sequences[0].annotations
{}
>>> alignment.sequences[1].annotations
{'leftStatus': 'I', 'leftCount': 9085, 'rightStatus': 'C', 'rightCount': 0}
顯示在此區塊和前一個區塊之間插入了 9085 個鹼基 (”I”)，而此區塊與後續區塊是連續的 (”C”)。有關這些欄位和狀態字元的完整說明，請參閱 UCSC 文件。
“q” 行顯示序列品質，該品質儲存在相應序列的 annotations 屬性的 “quality” 字典鍵下
>>> alignment.sequences[2].annotations["quality"]
'9999999999999999999999999999999999999999999999'
>>> alignment.sequences[4].annotations["quality"]
'9999969989999999999999998699989997'
“e” 行顯示物種的資訊，這些物種在此對齊區塊之前和之後具有連續序列，但在這個對齊區塊中沒有對齊的核苷酸。此資訊儲存在 alignment.annotations 字典的 “empty” 鍵下
>>> alignment.annotations["empty"]
[(SeqRecord(seq=Seq(None, length=498454), id='tupBel1.scaffold_114895.1-498454', name='', description='', dbxrefs=[]), (331078, 326933), 'I'),
(SeqRecord(seq=Seq(None, length=100002), id='echTel1.scaffold_288249', name='', description='', dbxrefs=[]), (87661, 95225), 'I'),
(SeqRecord(seq=Seq(None, length=359464), id='otoGar1.scaffold_334.1-359464', name='', description='', dbxrefs=[]), (178247, 175316), 'I'),
(SeqRecord(seq=Seq(None, length=174210431), id='ponAbe2.chr6', name='', description='', dbxrefs=[]), (158048983, 158040939), 'I')]
例如，這顯示在 ponAbe2.chr6 基因組序列的相反鏈上，從位置 158040939 到 158048983 插入了不對齊的鹼基 (”I”)。同樣，有關 “e” 行的完整定義，請參閱 UCSC 文件。
若要以 MAF 格式列印對齊結果，您可以使用 Python 的內建 format 函式
>>> print(format(alignment, "MAF"))
a score=19159.000000
s mm9.chr10 3014644 45 + 129993255 CCTGTACC---CTTTGGTGAGAATTTTTGTTTCAGTGTTAAAAGTTTG
s hg18.chr6 15870786 46 - 170899992 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
i hg18.chr6 I 9085 C 0
s panTro2.chr6 16389355 46 - 173908612 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
q panTro2.chr6 99999999999999999999999-9999999999999999999-9999
i panTro2.chr6 I 9106 C 0
s calJac1.Contig6394 6182 46 + 133105 CCTATACCTTTCTTTCATGAGAA-TTTTGTTTGAATCCTAAAC-TTTT
i calJac1.Contig6394 N 0 C 0
s loxAfr1.scaffold_75566 1167 34 - 10574 ------------TTTGGTTAGAA-TTATGCTTTAATTCAAAAC-TTCC
q loxAfr1.scaffold_75566 ------------99999699899-9999999999999869998-9997
i loxAfr1.scaffold_75566 N 0 C 0
e tupBel1.scaffold_114895.1-498454 167376 4145 - 498454 I
e echTel1.scaffold_288249 87661 7564 + 100002 I
e otoGar1.scaffold_334.1-359464 181217 2931 - 359464 I
e ponAbe2.chr6 16161448 8044 - 174210431 I
若要寫入完整的 MAF 檔案，請使用 Bio.Align.write(alignments, "myfilename.maf", fmt="maf")，其中 myfilename.maf 是輸出 MAF 檔案的名稱。或者，您可以使用 (文字) 串流進行輸出。檔案標頭資訊將取自 alignments 物件的 metadata 屬性。如果您要從頭開始建立對齊結果，可以使用 Alignments (複數) 類別來建立類似列表的 alignments 物件 (請參閱章節 Alignments 類別)，並賦予它一個 metadata 屬性。
bigMaf
bigMaf 檔案是一個 bigBed 檔案，其格式為 BED3+1，由 3 個必要的 BED 欄位加上一個自訂欄位組成，該欄位將 MAF 對齊區塊儲存為字串，從而建立 MAF 檔案的索引二進位版本 (請參閱章節 多重對齊格式 (MAF))。相關聯的 AutoSql 檔案 bigMaf.as 由 UCSC 提供。若要建立 bigMaf 檔案，您可以使用 UCSC 的 mafToBigMaf 和 bedToBigBed 程式，或者您也可以使用 Biopython 並使用 fmt="bigmaf" 呼叫 Bio.Align.write 函式。雖然這兩種方法應產生相同的 bigMaf 檔案，但 UCSC 工具速度更快，可能更可靠，因為它是黃金標準。由於 bigMaf 檔案是 bigBed 檔案，因此它們具有內建索引，可讓您快速搜尋參考基因組的特定區域。
Biopython 發行版的 Tests/MAF 子目錄中的檔案 ucsc_test.bb 是 bigMaf 檔案的範例。此檔案等同於 MAF 檔案 ucsc_test.maf (請參閱章節 多重對齊格式 (MAF))。若要剖析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("ucsc_test.bb", "bigmaf")
>>> len(alignments)
3
>>> print(alignments.declaration)
table bedMaf
"Bed3 with MAF block"
(
string chrom; "Reference sequence chromosome or scaffold"
uint chromStart; "Start position in chromosome"
uint chromEnd; "End position in chromosome"
lstring mafBlock; "MAF block"
)
宣告包含 UCSC 的 bigMaf.as AutoSql 檔案定義的欄規格。
bigMaf 檔案不會儲存 MAF 檔案中的標頭資訊，但它確實定義了參考基因組。相應的 SeqRecord 儲存在 alignments 物件的 targets 屬性中
>>> alignments.reference
'hg16'
>>> alignments.targets
[SeqRecord(seq=Seq(None, length=158545518), id='hg16.chr7', ...)]
透過疊代 alignments，我們為 bigMaf 檔案中的每個對齊區塊取得一個 Alignment 物件
>>> alignment = next(alignments)
>>> alignment.score
23262.0
>>> {seq.id: len(seq) for seq in alignment.sequences}
{'hg16.chr7': 158545518,
'panTro1.chr6': 161576975,
'baboon': 4622798,
'mm4.chr6': 151104725,
'rn3.chr4': 187371129}
>>> print(alignment.coordinates)
[[27578828 27578829 27578831 27578831 27578850 27578850 27578866]
[28741140 28741141 28741143 28741143 28741162 28741162 28741178]
[ 116834 116835 116837 116837 116856 116856 116872]
[53215344 53215344 53215346 53215347 53215366 53215366 53215382]
[81344243 81344243 81344245 81344245 81344264 81344267 81344283]]
>>> print(alignment)
hg16.chr7 27578828 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG 27578866
panTro1.c 28741140 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG 28741178
baboon 116834 AAA-GGGAATGTTAACCAAATGA---GTTGTCTCTTATGGTG 116872
mm4.chr6 53215344 -AATGGGAATGTTAAGCAAACGA---ATTGTCTCTCAGTGTG 53215382
rn3.chr4 81344243 -AA-GGGGATGCTAAGCCAATGAGTTGTTGTCTCTCAATGTG 81344283
>>> print(format(alignment, "phylip"))
5 42
hg16.chr7 AAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
panTro1.chAAA-GGGAATGTTAACCAAATGA---ATTGTCTCTTACGGTG
baboon AAA-GGGAATGTTAACCAAATGA---GTTGTCTCTTATGGTG
mm4.chr6 -AATGGGAATGTTAAGCAAACGA---ATTGTCTCTCAGTGTG
rn3.chr4 -AA-GGGGATGCTAAGCCAATGAGTTGTTGTCTCTCAATGTG
“i”、“e” 和 “q” 行中的資訊以與相應 MAF 檔案相同的方式儲存 (請參閱章節 多重對齊格式 (MAF))
>>> from Bio import Align
>>> alignments = Align.parse("ucsc_mm9_chr10.bb", "bigmaf")
>>> for i in range(10):
... alignment = next(alignments)
...
>>> alignment.score
19159.0
>>> print(alignment)
mm9.chr10 3014644 CCTGTACC---CTTTGGTGAGAATTTTTGTTTCAGTGTTAAAAGTTTG 3014689
hg18.chr6 155029206 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT 155029160
panTro2.c 157519257 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT 157519211
calJac1.C 6182 CCTATACCTTTCTTTCATGAGAA-TTTTGTTTGAATCCTAAAC-TTTT 6228
loxAfr1.s 9407 ------------TTTGGTTAGAA-TTATGCTTTAATTCAAAAC-TTCC 9373
>>> print(format(alignment, "MAF"))
a score=19159.000000
s mm9.chr10 3014644 45 + 129993255 CCTGTACC---CTTTGGTGAGAATTTTTGTTTCAGTGTTAAAAGTTTG
s hg18.chr6 15870786 46 - 170899992 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
i hg18.chr6 I 9085 C 0
s panTro2.chr6 16389355 46 - 173908612 CCTATACCTTTCTTTTATGAGAA-TTTTGTTTTAATCCTAAAC-TTTT
q panTro2.chr6 99999999999999999999999-9999999999999999999-9999
i panTro2.chr6 I 9106 C 0
s calJac1.Contig6394 6182 46 + 133105 CCTATACCTTTCTTTCATGAGAA-TTTTGTTTGAATCCTAAAC-TTTT
i calJac1.Contig6394 N 0 C 0
s loxAfr1.scaffold_75566 1167 34 - 10574 ------------TTTGGTTAGAA-TTATGCTTTAATTCAAAAC-TTCC
q loxAfr1.scaffold_75566 ------------99999699899-9999999999999869998-9997
i loxAfr1.scaffold_75566 N 0 C 0
e tupBel1.scaffold_114895.1-498454 167376 4145 - 498454 I
e echTel1.scaffold_288249 87661 7564 + 100002 I
e otoGar1.scaffold_334.1-359464 181217 2931 - 359464 I
e ponAbe2.chr6 16161448 8044 - 174210431 I
>>> alignment.sequences[1].annotations
{'leftStatus': 'I', 'leftCount': 9085, 'rightStatus': 'C', 'rightCount': 0}
>>> alignment.sequences[2].annotations["quality"]
'9999999999999999999999999999999999999999999999'
>>> alignment.sequences[4].annotations["quality"]
'9999969989999999999999998699989997'
>>> alignment.annotations["empty"]
[(SeqRecord(seq=Seq(None, length=498454), id='tupBel1.scaffold_114895.1-498454', name='', description='', dbxrefs=[]), (331078, 326933), 'I'),
(SeqRecord(seq=Seq(None, length=100002), id='echTel1.scaffold_288249', name='', description='', dbxrefs=[]), (87661, 95225), 'I'),
(SeqRecord(seq=Seq(None, length=359464), id='otoGar1.scaffold_334.1-359464', name='', description='', dbxrefs=[]), (178247, 175316), 'I'),
(SeqRecord(seq=Seq(None, length=174210431), id='ponAbe2.chr6', name='', description='', dbxrefs=[]), (158048983, 158040939), 'I')]
若要寫入完整的 bigMaf 檔案，請使用 Bio.Align.write(alignments, "myfilename.bb", fmt="bigMaf")，其中 myfilename.bb 是輸出 bigMaf 檔案的名稱。或者，您可以使用 (二進位) 串流進行輸出。如果您要從頭開始建立對齊結果，可以使用 Alignments (複數) 類別來建立類似列表的 alignments 物件 (請參閱章節 Alignments 類別)，並賦予它一個 targets 屬性。後者必須是參考物種中染色體的 SeqRecord 物件列表，其順序與它們在對齊結果中出現的順序相同。或者，您可以在呼叫 Bio.Align.write 時使用 targets 關鍵字引數。每個 SeqRecord 的 id 必須採用 reference.chromosome 的形式，其中 reference 指的是參考物種。Bio.Align.write 具有附加的關鍵字引數 compress (預設為 True)，用於指定是否應使用 zlib 壓縮資料。其他可選引數為 blockSize (預設值為 256) 和 itemsPerSlot (預設值為 512)。有關這些引數的說明，請參閱 UCSC 的 bedToBigBed 程式的文件。
由於 bigMaf 檔案是 bigBed 檔案的特殊情況，您可以在 alignments 物件上使用 search 方法，以尋找參考物種特定區域的對齊結果。例如，我們可以尋找染色體 10 上位置 3018000 和 3019000 之間的 chr10 區域
>>> selected_alignments = alignments.search("mm9.chr10", 3018000, 3019000)
>>> for alignment in selected_alignments:
... start, end = alignment.coordinates[0, 0], alignment.coordinates[0, -1]
... print(start, end)
...
3017743 3018161
3018161 3018230
3018230 3018359
3018359 3018482
3018482 3018644
3018644 3018822
3018822 3018932
3018932 3019271
染色體名稱可以是 None 以包含所有染色體，開始和結束位置可以是 None 以分別從位置 0 開始搜尋或繼續搜尋到染色體末尾。請注意，我們只能搜尋參考物種的基因組位置。
在建立 bigMaf 檔案時，使用 compress=False 和 itemsPerSlot=1 可以加快搜尋 bigMaf 檔案的速度。
UCSC chain 檔案格式
Chain 檔案描述兩個核苷酸序列之間的成對對齊，允許兩個序列中存在間隙。Chain 檔案中僅儲存每個對齊子序列的長度和間隙長度；序列本身不儲存。Chain 檔案通常用於儲存兩個基因組組裝版本之間的對齊結果，允許將一個基因組組裝版本的對齊結果提升到另一個基因組組裝版本。這是一個 chain 檔案的範例 (在 Biopython 發行版的 Tests/Blat 子目錄中，檔案名為 psl_34_001.chain)
chain 16 chr4 191154276 + 61646095 61646111 hg18_dna 33 + 11 27 1
16
chain 33 chr1 249250621 + 10271783 10271816 hg18_dna 33 + 0 33 2
33
chain 17 chr2 243199373 + 53575980 53575997 hg18_dna 33 - 8 25 3
17
chain 35 chr9 141213431 + 85737865 85737906 hg19_dna 50 + 9 50 4
41
chain 41 chr8 146364022 + 95160479 95160520 hg19_dna 50 + 8 49 5
41
chain 30 chr22 51304566 + 42144400 42144436 hg19_dna 50 + 11 47 6
36
chain 41 chr2 243199373 + 183925984 183926028 hg19_dna 50 + 1 49 7
6 0 4
38
chain 31 chr19 59128983 + 35483340 35483510 hg19_dna 50 + 10 46 8
25 134 0
11
chain 39 chr18 78077248 + 23891310 23891349 hg19_dna 50 + 10 49 9
39
...
此檔案是透過在 PSL 檔案 psl_34_001.psl 上執行 UCSC 的 pslToChain 程式所產生。根據 chain 檔案格式規格，每個 chain 區塊之後應該要有一個空白行，但某些工具（包括 pslToChain）顯然沒有遵循此規則。
要解析此檔案，請使用
>>> from Bio import Align
>>> alignments = Align.parse("psl_34_001.chain", "chain")
迭代處理比對結果，為每個 chain 取得一個 Alignment 物件
>>> for alignment in alignments:
... print(alignment.target.id, alignment.query.id)
...
chr4 hg18_dna
chr1 hg18_dna
chr2 hg18_dna
chr9 hg19_dna
chr8 hg19_dna
chr22 hg19_dna
chr2 hg19_dna
...
chr1 hg19_dna
從頭開始迭代，直到抵達第七個比對結果
>>> alignments = iter(alignments)
>>> for i in range(7):
... alignment = next(alignments)
...
檢查比對分數和 chain ID（每個 chain 區塊的標頭行中的第一個和最後一個數字），以確認我們取得的是第七個比對結果
>>> alignment.score
41.0
>>> alignment.annotations["id"]
'7'
我們可以以 chain 檔案格式印出比對結果。比對坐標與 chain 區塊中的資訊一致，其中包含一段 6 個核苷酸的比對區段、一段 4 個核苷酸的間隙，以及一段 38 個核苷酸的比對區段
>>> print(format(alignment, "chain"))
chain 41 chr2 243199373 + 183925984 183926028 hg19_dna 50 + 1 49 7
6 0 4
38
>>> alignment.coordinates
array([[183925984, 183925990, 183925990, 183926028],
[ 1, 7, 11, 49]])
>>> print(alignment)
chr2 183925984 ??????----?????????????????????????????????????? 183926028
0 ||||||----|||||||||||||||||||||||||||||||||||||| 48
hg19_dna 1 ???????????????????????????????????????????????? 49
我們也可以用其他幾種比對檔案格式印出比對結果
>>> print(format(alignment, "BED"))
chr2 183925984 183926028 hg19_dna 41 + 183925984 183926028 0 2 6,38, 0,6,
>>> print(format(alignment, "PSL"))
44 0 0 0 1 4 0 0 + hg19_dna 50 1 49 chr2 243199373 183925984 183926028 2 6,38, 1,11, 183925984,183925990,
>>> print(format(alignment, "exonerate"))
vulgar: hg19_dna 1 49 + chr2 183925984 183926028 + 41 M 6 6 G 4 0 M 38 38
>>> print(alignment.format("exonerate", "cigar"))
cigar: hg19_dna 1 49 + chr2 183925984 183926028 + 41 M 6 I 4 M 38
>>> print(format(alignment, "sam"))
hg19_dna 0 chr2 183925985 255 1S6M4I38M1S * 0 0 * * AS:i:41 id:A:7