序列比对

序列比对是两个或多个已经相互比对的序列的集合-通常插入空位并添加引导或尾部空位-使得所有序列串的长度相同。

对齐可以延伸到每个序列的全长，或者可以限制到每个序列的一个小节。在Biopython中，所有序列比对均由 Alignment 对象，在部分中描述路线对象 . Alignment 对象可以通过解析比对软件（例如Classal或BLAT）的输出来获得（在部分中描述阅读和写作对齐 .或者通过使用Biopython的成对序列对齐器，它可以将两个序列相互对齐（参见第章双序列比对 ).

参见Chapter 多个序列比对对象了解老年人的描述 MultipleSeqAlignment 类和中的解析器 Bio.AlignIO 分析序列比对软件的输出，生成 MultipleSeqAlignment 对象

路线对象

的 Alignment 类定义在 Bio.Align .通常你会得到一个 Alignment 通过解析对齐程序的输出来创建对象（部分阅读和写作对齐）或通过运行Biopython的成对对齐器（第章双序列比对 ).然而，为了本节的好处，我们将创建一个 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]]

初始 G 核苷酸 seqA 和最终 CC 核苷酸 seqB 未包含在对齐中，因此此处缺失。但这很容易解决：

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

which identical to the Alignment object created above in section 根据序列和坐标创建对齐对象.

默认情况下 coordinates 论点 Alignment 初始值设定项 None ，它假设对齐中没有间隙。无缺口比对中的所有序列必须具有相同的长度。如果 coordinates 论点是 None ，那么初始化器将填写 coordinates 属性 Alignment 对象为您：

>>> 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 ;
- bytes ;
- bytearray ;
- 具有数据类型的NumPy数组 numpy.int32 ;
- 具有连续格式缓冲区的任何其他对象 "c" , "B" , "i" ，或者 "I" ;
- 中定义的对象列表或数组 alphabet 属性 PairwiseAligner 创建路线的对象（请参阅部分广义成对排列 ).
对于成对比对（意味着两个序列的比对），属性 target 和 query 是别名 sequences[0] 和 sequences[1] ，分别。
coordinates ：NumPy的整元数组，存储定义序列如何相互对齐的序列索引;
score ：对齐分数，由解析器在对齐文件中找到，或由 PairwiseAligner (see部分基本用法 );
annotations ：存储与对齐相关的大多数其他注释的字典;
column_annotations ：存储注释的字典，这些注释沿着比对延伸并具有与比对相同的长度，例如共有序列（参见第 ClustalW 例如）。

一个 Alignment 由解析器创建的对象 Bio.Align 可能具有其他属性，具体取决于读取对齐的对齐文件格式。

对路线进行切片和索引

形式的切片 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 == alignment2 计算结果为 True ）如果每个序列在 alignment1.sequences 和 alignment2.sequences 彼此相等，并且 alignment1.coordinates 和 alignment2.coordinates 都有相同的坐标如果这两个条件中的任何一个不满足， 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 对象，从中可以获得计数作为属性。

>>> print(pairwise_alignment)
target            1 CGGTTTTT 9
                  0 .|-|||-- 8
query             0 AG-TTT-- 5

>>> counts = pairwise_alignment.counts()
>>> counts.aligned
5
>>> counts.identities
4
>>> counts.mismatches
1
>>> counts.gaps
3
>>> counts.insertions
0
>>> counts.deletions
3
>>> counts.internal_deletions
1
>>> counts.right_deletions
2

使用 wildcard 参数，用于指定在计算身份、阳性和不匹配时应忽略的字母（例如 wildcard="?" 或 wildcard="N" 是常见的选择）。

对于两个以上序列的比对，计算并总和比对中所有序列对的同一性、错配和缺口的数量。

>>> print(alignment)
                  1 CGGTTTTT 9
                  0 AG-TTT-- 5
                  0 AGGTTT-- 6

>>> counts = alignment.counts()
>>> counts.aligned
16
>>> counts.identities
14
>>> counts.mismatches
2
>>> counts.insertions
1
>>> counts.deletions
5

这里，插入被定义为在比对中将较后序列插入到较早序列中的序列。与上面的成对比对相反，插入和删除之间的区别对于多序列比对可能没有意义，并且您可能会对缺口的数量（=插入+删除）更感兴趣：

>>> counts.gaps
6
>>> counts.left_gaps
0
>>> counts.right_gaps
4
>>> counts.internal_gaps
2

为了加快计算速度，您可以使用 ignore_sequences=True 跳过计算匹配和不匹配的数量（这仍然会计算对齐序列的数量）：

>>> counts = alignment.counts(ignore_sequences=True)
>>> counts.aligned
16
>>> print(counts.identities)
None
>>> print(counts.mismatches)
None
>>> counts.insertions
1
>>> counts.deletions
5

对于蛋白质比对，除了身份和错配的数量之外，您还可以通过提供替代矩阵来计算阳性匹配的数量（参见第章取代矩阵 ):

>>> from Bio.Align import substitution_matrices
>>> substitution_matrix = substitution_matrices.load("BLOSUM62")
>>> protein_alignment = Alignment(["GQCGSCWSFS", "GACGSCWTFS"])
>>> print(protein_alignment)
target            0 GQCGSCWSFS 10
                  0 |.|||||.|| 10
query             0 GACGSCWTFS 10

>>> counts = protein_alignment.counts(substitution_matrix)
>>> counts.aligned
10
>>> counts.identities
8
>>> counts.mismatches
2
>>> counts.positives
9

其中阳性数量比身份数量高一个，因为S与T的不匹配具有正分数（而Q与A的不匹配分数不为正）：

>>> substitution_matrix["S", "T"]
1.0
>>> substitution_matrix["Q", "A"]
-1.0

一个 AlignmentCounts 对象具有以下属性：

property	description
`aligned`	对齐中对齐字母的数量。如果部分或所有序列未定义，也会计算此量。如果所有序列都是已知的，那么 `aligned` = `identities` + `mismatches` .如果某些序列是未定义的，则 `aligned` > `identities` + `mismatches` .
`identities`	对齐中相同字母的数量
`mismatches`	对齐中不匹配的字母的数量
`positives`	具有正值的对齐字母的数量
`left_insertions`	对齐左侧的插入数量
`left_deletions`	对齐左侧的删除数量
`right_insertions`	对齐右侧的插入数量
`right_deletions`	对齐右侧的删除数量
`internal_insertions`	对齐内部的插入数量
`internal_deletions`	对齐内部的删除数量
`insertions`	插入总数，等于 `left_insertions` + `right_insertions` + `internal_insertions`
`deletions`	删除总数，等于 `left_deletions` + `right_deletions` + `internal_deletions`
`left_gaps`	路线左侧的间隙数量，等于 `left_insertions` + `left_deletions`
`right_gaps`	对齐右侧的间隙数量，等于 `right_insertions` + `right_deletions`
`internal_gaps`	对齐内部的间隙数量，等于 `internal_insertions` + `internal_deletions`
`gaps`	对齐中的间隙总数，等于 `left_gaps` + `right_gaps` + `internal_gaps`

信件频率

的 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替换的次数。例如， C 与 A 在后面的序列中是

>>> 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 对齐中的' s是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' ).您可以使用创建Unicode（Python字符串）字符数组 dtype='U' :

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

(the印刷 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 函数应用于完整序列（包括初始序列 A 和最终 GG 核苷酸 seqB ），而不仅仅是对齐的部分。

的 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 对象，可以添加在一起。参阅Chapter 顺序注释对象了解如何处理注释的详细信息。此示例是一个特殊情况，因为两个原始路线共享相同的名称，这意味着当添加行时，它们也会获得相同的名称。

任何常见的注释都会保留，但不同的注释会丢失。这与中使用的行为相同 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

以及转录物和序列之间的成对比对（例如，通过RN-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

使用 map 方法对 alignment1 ， 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 .

地图方法还可用于提升不同基因组组装之间的对齐。在这种情况下，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]]

表明黑猩猩基因组组装panTro 5上chr1的范围122250000：122909835与黑猩猩基因组组装panTro 6上chr1的范围111776384：112019978对齐。见章节链文件格式有关链文件格式的更多信息。

>>> 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]]

这表明表达序列标签DC 525629的核苷酸范围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 对应于目标 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

现在我们可以通过调用来获得DC 525629针对黑猩猩基因组组装panTro 6的坐标 chain.map ， transcript 作为论点：

>>> 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]]

这表明表达的序列标签DC 525629的第32：344个碱基范围与黑猩猩基因组组装panTro 6的chr 1的第111982717：112009302个碱基范围对齐。请注意，DC 525629在黑猩猩基因组组装panTro 5上的基因组跨度为122907314 - 122835789 = 71525个碱基对，而在panTro 6上的基因组跨度为112009302 - 111982717 = 26585个碱基对。

绘制多序列比对

考虑黑猩猩、人类、猕猴、绒猴、小鼠和大鼠基因组序列的多重比对：

>>> 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作为链文件提供，通常用于UCSC liftOver 工具.的 .chain 中的文件 Tests/Align Biopython来源分布中的条目是从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 (see部分阅读对齐）或由成对对齐器返回的迭代器（参见第双序列比对 ).调用 iter 在迭代器上将始终返回 Alignments 对象本身。相比之下，呼叫 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 ...>

相反， alignment_iterator1 和 alignment_iterator2 通过致电获得 iter 上 Alignments 对象彼此相同：

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

调用 iter 上 Alignments 对象将迭代器重置为它的第一个项，这样你就可以再次循环它。您还可以使用 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 (see部分阅读对齐）或通过成对对齐器（请参阅部分双序列比对 ).

元数据可以作为属性存储在 Alignments 对象，而普通的 list 不接受属性：

>>> alignment_list.score = 100
Traceback (most recent call last):
 ...
AttributeError: 'list' object has no attribute 'score'...
>>> alignments.score = 100
>>> alignments.score
100

阅读和写作对齐

序列比对软件（例如Classal）的输出可以解析为 Alignment 物体由 Bio.Align.read 和 Bio.Align.parse 功能协调发展的它们的用途类似于 read 和 parse 功能 Bio.SeqIO (see部分解析或读取序列）： read 函数用于读取包含单个对齐的输出文件并返回 Alignment 对象，而 parse 函数返回迭代器，以迭代存储在包含一个或多个对齐的输出文件中的对齐。部分对齐文件格式描述了可以在其中解析的对齐格式 Bio.Align . Bio.Align 还提供了一种 write 可以以大多数格式编写对齐的函数。

阅读对齐

使用 Bio.Align.parse 解析序列比对文件。例如文件 ucsc_mm9_chr10.maf 包含48个以APM（多重比对格式）格式进行的多重序列比对（请参阅部分多重对齐格式（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 可能包含存储文件中找到的元数据的属性，例如用于创建对齐的软件的版本号。为每种文件格式存储的特定属性在部分中描述对齐文件格式 .对于APM文件，我们可以获取所使用的文件格式版本和评分方案：

>>> 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 很快就会回来。

如果对齐的数量不是太大，并且内存中可以容纳，那么可以将对齐迭代器转换为对齐列表。为此，您可以调用 list 上 alignments :

>>> 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课程），允许您将元数据字典存储为 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 对对齐方式进行函数以获取以请求的格式显示对齐方式的字符串，或使用 Alignment 格式化（f-）字符串中的对象。如果在没有论据的情况下调用， 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.Alignn中解析的对齐格式。format参数 fmt 用于 Bio.Align 指定文件格式不区分大小写的函数。这些文件格式中的大多数也可以由 Bio.Align ，如表中所示。

文件格式 `fmt`	描述	文本/二进制	支持 `write`	款
`a2m`	A2m	文本	是的	1.7.11
`bed`	浏览器可扩展数据（BED）	文本	是的	1.7.14
`bigbed`	bigBed	二进制	是的	1.7.15
`bigmaf`	bigMaf	二进制	是的	1.7.19
`bigpsl`	bigPsl	二进制	是的	1.7.17
`chain`	UCSC链文件	文本	是的	1.7.20
`clustal`	ClustalW	文本	是的	1.7.2
`emboss`	EMBOSS	文本	没有	1.7.5
`exonerate`	免除	文本	是的	1 .7.7
`fasta`	对齐FASTA	文本	是的	1.7.1
`hhr`	HH套件输出文件	文本	没有	1.7.10
`maf`	多重对齐格式（MAF）	文本	是的	1.7.18
`mauve`	Mauve eXtented Multi-FastA（xmfa）格式	文本	是的	1.7.12
`msf`	GCG多序列格式（MSF）	文本	没有	1.7.6
`nexus`	NEXUS	文本	是的	1.7.8
`phylip`	PHYLIP输出文件	文本	是的	1.7.4
`psl`	图案空间布局（PSL）	文本	是的	1.7.16
`sam`	序列比对/图谱（Sam）	文本	是的	1.7.13
`stockholm`	斯德哥尔摩	文本	是的	1 .7.3
`tabular`	来自AMPS或FASTA的表格输出	文本	没有	1.7.9

对齐FASTA

对齐FASTA格式的文件恰好存储一个（成对或多个）序列对齐，其中对齐中的间隙由虚线表示 (- ).使用 fmt="fasta" 以对齐的FASTA格式读取或写入文件。请注意，这与William Pearson的FASTA对齐程序生成的输出不同（解析此类输出在部分中描述来自AMPS或FASTA的表格输出相反）。

文件 probcons.fa 在Biopython的测试套件中，以对齐的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

或任何其他可用的格式，例如Classal（请参阅部分 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

Classal是一组多序列比对程序，既可以作为独立程序也可以作为网络服务器使用。文件 opuntia.aln （在线或在 Doc/examples Biopython源代码的副本）是由Classal生成的输出文件。它的前几行是

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

的 metadata 属性对 alignments 存储文件头中显示的信息：

>>> alignments.metadata
{'Program': 'CLUSTAL', 'Version': '2.1'}

你可以叫 next 上 alignments 要拉出第一个（也是唯一一个）对齐：

>>> 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 函数，因为无论如何每个Classal文件只包含一个对齐：

>>> from Bio import Align
>>> alignment = Align.read("opuntia.aln", "clustal")

如果序列在每个位置都保守，则Classal输出文件中每个比对块下方的共识线包含星号。此信息存储在 column_annotations 属性 alignment :

>>> alignment.column_annotations
{'clustal_consensus': '******* **** **********************************...

打印 alignment 在 clustal 格式将显示序列比对，但不包括元数据：

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

写 alignments 在 clustal 格式将包括元数据和序列比对：

>>> 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 （复数）对象（请参阅部分 Alignments课程）如果您正在手工创建对齐，并且希望在输出中包括元数据信息。

斯德哥尔摩

这是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
//

这是从InterPro网站（https：//www.ebi.ac.uk/interpro/）下载的7kD_DNA_binding（PF 02294.20）PFAM条目的种子对齐。PFAM条目的此版本也可以在Biopython源发行版中作为文件提供 pfam2.seed.txt 子目录中 Tests/Stockholm/ .我们可以如下方式加载此文件：

>>> 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 对对齐对象进行函数，以特定文件格式显示它（请参阅部分印刷对齐详情），例如以斯德哥尔摩格式重新生成文件：

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

路线的一般信息存储在 annotations 属性 Alignment 例如，对象

>>> 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;'}]

该比对中的各个序列存储在 alignment.sequences 作为 SeqRecord s，包括与每个序列记录相关的任何注释：

>>> 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 ）存储在 letter_annotations 属性 SeqRecord :

>>> alignment.sequences[0].letter_annotations
{}
>>> alignment.sequences[1].letter_annotations
{'secondary structure': 'EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT'}
>>> alignment.sequences[2].letter_annotations
{}

共有序列和二级结构与整个序列比对相关，因此存储在 column_annotations 属性 Alignment 对象：

>>> alignment.column_annotations
{'consensus secondary structure': 'EEEEESSSSEEEEETTTEEEEEESSSSEEEEEE-SSSSEEEEEEETTTS-CHHHHHHTT',
 'consensus sequence': 'KVKFKYKGEEKEVDISKIKKVWRVGKMVSFTYDD.NGKTGRGAVSEKDAPKELLsMLuK'}

PHYLIP输出文件

用于序列比对的PHYLIP格式源自Joe Felsenstein的PHYLogogeny Interference Pack。PHYLIP格式的文件以两个数字开头，代表印刷对齐中的行和列数。序列比对本身可以是顺序格式或交错格式。前者的一个例子是 sequential.phy 文件（提供在 Tests/Phylip/ 在Biopython源分发版中）：

 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 （提供于 Tests/Phylip/ 在Biopython源分发版中）：

 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

CLARSS（欧洲分子生物学开放软件套件）是一套用于分子生物学和生物信息学的开源软件工具 [Rice2000]. 它包括以下软件 needle 和 water 用于成对序列比对。这是由 water Smith-Waterman局部成对序列比对程序（可作为 water.txt 在 Tests/Emboss Biopython发行版目录）：

########################################
# 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")

的 metadata 属性 alignments 存储文件头中显示的信息，包括用于生成输出的程序、运行程序的日期和时间、输出文件名以及使用的特定对齐文件格式（假定为 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'

的 annotations 属性 alignment 专门存储与此对齐相关的信息：

>>> alignment.annotations
{'Matrix': 'EBLOSUM62', 'Gap_penalty': 10.0, 'Extend_penalty': 0.5, 'Identity': 112, 'Similarity': 112, 'Gaps': 19, 'Score': 591.5}

还可以通过调用 counts 方法对 alignment 对象：

>>> counts = alignment.counts()

其中 counts 变量是一个 AlignmentCounts 对象收集有关对齐库中间隙、匹配和不匹配数量的信息（请参阅计算身份、不匹配和差距 )):

>>> counts.identities
112
>>> counts.mismatches
0
>>> counts.gaps
19

两个序列之间显示的共识线存储在 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不支持以CLARSS定义的输出格式编写序列对齐。

GCG多序列格式（MSF）

创建多序列格式（MSF）是为了存储GCG（遗传学计算机组）程序集生成的多个序列比对。文件 W_prot.msf 在 Tests/msf Biopython发行版的目录是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格式编写序列对齐。

免除

Exconerate是一个用于成对序列比对的通用程序 [Slater2005]. Exonerate找到的序列比对可以以人类可读的形式、“雪茄”（紧凑的独特间隙比对报告）格式或“粗俗”（Verbose Useful Labeled Gapped比对报告）格式输出。用户可以请求在输出中包括一种或多种这些格式。中的解析器 Bio.Align 只能分析雪茄或粗俗格式的对齐，并且不会分析包括人类可读格式的对齐的输出。

文件 exn_22_m_cdna2genome_vulgar.exn Biopython测试套件中有一个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，没有间隙：

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

请注意，打印比对中的目标（第一个序列）在反向链上，而查询（第二个序列）在正向链上，目标坐标减少，查询坐标增加。打印此对齐方式 exonerate 使用Python内置格式 format 函数写了一句粗俗的行：

>>> 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 方法允许我们请求粗俗行（默认）或雪茄行：

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

粗俗行包含有关路线的信息（在部分中 M 1 1 C 3 3 M 1226 ）雪茄系列中缺少的 M 1 M 3 M 1226 .粗俗的行指定对齐从单个对齐的核酸开始，然后是形成密码子的三个对齐的核酸 (C ），然后是1226个对齐的核酸。在雪茄系列中，我们看到一个对齐的核苷，然后是三个对齐的核苷，然后是1226个对齐的核苷;它没有指定三个对齐的核苷构成密码子。来自粗俗行的此信息存储在 operations 属性：

>>> alignment.operations
bytearray(b'MCM')

有关其他操作代码的定义，请参阅免责文档。

同样， "vulgar" 或 "cigar" 调用时可以使用参数 Bio.Align.write 用粗俗或雪茄对齐线编写文件。

我们可以以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）），这与 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

第三对线包含四个长间隙：

>>> 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格式的文件示例（可用为 codonposset.nex 在 Tests/Nexus Biopython发行版中的格式）：

#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 严重依赖NEXUS解析器 Bio.Nexus (see章系统发生学与Bio.Phylo ）提取 Alignment NEXUS文件中的对象。

要读取此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 .

来自AMPS或FASTA的表格输出

表格输出中的对齐输出由FASTA对齐器生成 [Pearson1988] run with the -m 8CB 或 -m 8CC 论点，或通过BLAST [Altschul1990] run with the -outfmt 7 论点

文件 nucleotide_m8CC.txt 在 Tests/Fasta Biopython源发行版的收件箱是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")

文件头中显示的信息存储在 metadata 属性 alignments :

>>> 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=[])

对齐信息存储在 annotations 属性 alignment :

>>> alignment.annotations
{'% identity': 66.93,
 'mismatches': 42,
 'gap opens': 8,
 'evalue': 3.2e-07,
 'bit score': 45.0}

BST特别提供了许多通过包括或排除特定列来自定义表格输出的选项;有关详细信息，请参阅BST文档。此信息存储在字典中 alignment.annotations , alignment.target.annotations ，或者 alignment.query.annotations ，视情况而定。

HH套件输出文件

中的对齐文件 hhr 格式由以下人员生成 hhsearch 或 hhblits 在卫生间 [Steinegger2019]. 作为示例，请参阅该文件 2uvo_hhblits.hhr 在Biopython的测试套件中：

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

大部分标头信息存储在 metadata 属性 alignments :

>>> 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 (note查询ID在 hhr 文件缩写）：

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

序列内容将不完整（部分定义的序列;请参阅部分具有部分定义序列内容的序列）如果对齐没有延伸到整个序列上。

该对齐块的第二行以“''开始，显示了获取目标序列的隐马尔科夫模型的名称和描述。这些都存储在钥匙下 "hmm_name" 和 "hmm_description" 在 alignment.target.annotations 词典：

>>> 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预测的二级结构和DS SP确定的目标二级结构：

>>> 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}

成对对齐的置信度值存储在 "Confidence" 钥匙 alignment.column_annotations 字典此字典还存储每列的分数，显示在查询和每个对齐块的目标部分之间：

>>> alignment.column_annotations
{'column score': '||||++.++..||++.|||+|+|||.+.+||+++||.+.|++..+|+++++.++|....|||.++||+.+.+||+.+||.+++++|+.|+...+++.||++.|||.|||||...+||+.+||+++|++|.+|++.+++++|..+.|||+++-||+...||+++||++.|||',
 'Confidence': '799999999999999999999999999999999999999999999999999999999999999999999999999999999999999999998899999999999999999999999999999999999999999999999999999999999999999999999999986'}

A2m

A2 M文件是由以下人员创建的对齐文件 align2model 或 hmmscore 在Sam序列比对和建模软件系统中 [Krogh1994], [Hughey1996]. A2 M文件包含一个多重对齐。A2 M文件格式类似于对齐的FASTA（请参阅部分对齐FASTA ).然而，为了区分插入和缺失，A2 M使用破折号和句号来代表空白，以及对齐序列中的大写和大写字符。在仅包含匹配或删除的对齐列中，匹配用大写字母表示，删除用虚线表示。插入用大写字母表示，与插入对齐的间隙显示为句号。标题行以“''''开头，后面是序列的名称，还可以选择描述。

文件 probcons.a2m Biopython的测试套件中是A2 M文件的示例（请参阅部分对齐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

解析器分析A2 M文件中的破折号、句点以及大小写字母的模式，以确定列是匹配/不匹配/删除（“' D '）还是插入（“' I ''''）。此信息存储在 match 密钥 alignment.column_annotations 词典：

>>> alignment.column_annotations
{'state': 'DIDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDDIDDIDDDDDDDDDDDDDDDDDDDDDDDIIIDDDDDDDDDDDDDDDDDDDDDD'}

由于状态信息存储在 alignment ，我们可以以A2 M格式打印对齐：

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

同样，可以使用A2 M格式将对齐写入输出文件 Align.write (see部分书写对齐 ).

Mauve eXtented Multi-FastA（xmfa）格式

淡紫色 [Darling2004] 是一个用于构建多基因组比对的软件包。这些比对以eXtented Multi-FastA（xmfa）格式存储。具体取决于如何 progressiveMauve (the 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 for dog现在显式显示在输出文件中。此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] 存储成对序列比对，通常是针对参考基因组的下一代测序数据。文件 ex1.sam Biopython的测试套件中的一个是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文件格式规范，标志包含逐位标志4的行将被取消映射。由于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'

文件 ex1_header.sam Biopython的测试套件包含相同的路线，但现在还包括一个标题。它的前几行如下：

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

标头存储有关比对的一般信息，包括目标染色体的大小。目标信息存储在 targets 属性 alignments 对象：

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

文件中的对齐 sam1.sam Biopython测试套件中包含额外的 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文件可能包括额外信息，以区分基因组跳过区域（例如插入子）、硬剪切和软剪切以及填充序列区域的简单序列插入和缺失。由于此信息无法存储在 coordinates 的属性 Alignment 对象，并存储在专用的 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

在这种排列中，雪茄串 63M1062N75M468N43M 定义了46个对齐的核苷、一个由1827个核苷酸组成的Intron、82个对齐的核苷、一个由3376个核苷酸组成的Intron、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（浏览器可扩展数据）文件通常用于存储基因转录本与基因组的比对。看到 description from UCSC 了解床位格式的完整解释。

BED文件有三个必填字段和九个可选字段。文件 bed12.bed 在 Tests/Blat 是具有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]]

请注意，第一个序列（“mRNA 1”）被定位到正向链，而第二个序列（“mRNA 2”）被定位到反向链。

由于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中的信息（称为 thickStart , thickEnd ，而且 itemRgb 在BED格式规范中）作为属性存储在 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,

也可以使用 format 方法 alignment object.这允许您指定要写为 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文件，您可以使用 bedToBigBed program from UCSC (`) <https://genome.ucsc.edu/goldenPath/help/bigBed.html>`__ .或者您可以通过调用 Bio.Align.write 函数 fmt="bigbed" .虽然这两种方法应该会产生相同的bigBed文件，但使用 bedToBigBed 速度要快得多，而且可能更可靠，因为它是黄金标准。由于bigBed文件带有内置索引，因此可以快速搜索特定的基因组区域。

例如，让我们解析bigBed文件 dna_rna.bb ，可在 Tests/Blat Biopython发行版中的收件箱：

>>> 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 在打印对齐时显示序列内容，或以需要序列内容的格式编写对齐（例如Classal，请参阅部分 ClustalW ).测试脚本 test_Align_bigbed.py 在 Tests Biopython发行版中的收件箱提供了一些有关如何做到这一点的示例。

现在让我们看看如何搜索序列区域。这些是存储在bigBed文件中的序列，以BED格式打印（请参阅部分浏览器可扩展数据（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,

使用 search 方法对 alignments 目标在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字段的数量。例如，要编写BED 6文件，请使用

>>> Align.write(alignments, "output.bb", "bigbed", bedN=6)

一样的 bedToBigBed ，您可以在bigBed输出中包括其他列。假设文件 bedExample2.as （可在 Tests/Blat Biopython发行版的收件箱）以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字段以及其他字段的bigBed文件 geneSymbol 和 spID 通过调用

>>> Align.write(
...     alignments,
...     "output.bb",
...     "bigbed",
...     bedN=9,
...     declaration=declaration,
...     extraIndex=["name", "geneSymbol"],
... )

在这里，我们还要求纳入有关 name 和 geneSymbol 在bigBed文件中。 Align.write 希望找到钥匙 geneSymbol 和 spID 在 alignment.annotations 字典请参阅测试脚本 test_Align_bigbed.py 在 Tests 请在Biopython发行版中查找以bigBed格式编写对齐文件的更多示例。

可选参数 compress （默认值为 True ), blockSize （默认值为256），以及 itemsPerSlot （默认值为512）。请参阅UCSC的文档 bedToBigBed 描述这些论点的程序。搜索 bigBed 文件可以使用更快 compress=False 和 itemsPerSlot=1 创建bigBed文件时。

图案空间布局（PSL）

PSL（图案空间布局）文件是由类RST对齐工具BLAT生成的 [Kent2002]. 喜欢BED文件（请参阅部分浏览器可扩展数据（BED）），PSL文件通常用于存储转录本与基因组的比对。这是短BLAT文件的示例（提供为 dna_rna.psl 在 Tests/Blat Biopython发行版的收件箱），标准PSL标题由5行组成：

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 object.如果这些属性不可用，例如，如果对齐不是来自PSL文件，则使用序列内容（如果可用）计算这些数字。目标序列中的重复区域通过将序列掩蔽为大写或大写字符来指示，如 mask 关键字参数：

False （默认）：不单独计算与屏蔽序列的匹配项;
"lower" ：计数并报告与大小写字符的匹配，作为与重复区域的匹配;
"upper" ：计数并报告与大写字符的匹配项作为与重复区域的匹配项;

用于未知核苷的字符由 wildcard 论点为了与BLAT保持一致，RST字符是 "N" 在默认情况下使用 wildcard=None 如果您不想单独计算与任何未知核苷的匹配项。

>>> import numpy as np
>>> from Bio import Align
>>> query = "GGTGGGGG"
>>> target = "AAAAAAAggggGGNGAAAAA"
>>> coordinates = np.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,

The same arguments can be used when writing alignments to an output file in PSL format using Bio.Align.write. This function has an additional keyword header (True by default) specifying if the PSL header should be written.

除了有 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,

允许 Alignment Python中用于格式化（f-）字符串的对象：

>>> 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格式由AutoSql文件中定义的12个预定义BED字段和13个自定义字段组成 bigPsl.as 由UCSC提供，创建PSL文件的索引二进制版本（请参阅部分图案空间布局（PSL） ).要创建bigPsl文件，您可以使用 pslToBigPsl 和 bedToBigBed 来自UCSC的节目。或者您可以通过调用 Bio.Align.write 函数 fmt="bigpsl" .虽然这两种方法应该会产生相同的bigPsl文件，但UCSC工具速度要快得多，而且可能更可靠，因为它是黄金标准。由于bigPsl文件是bigBed文件，因此它们带有内置索引，允许您快速搜索特定的基因组区域。

例如，让我们解析bigBed文件 dna_rna.psl.bb ，可在 Tests/Blat 在Biopython发行版中删除。该文件是bigPsl相当于bigBed文件 dna_rna.bb (see部分 bigBed ）和PSL文件 dna_rna.psl (see部分图案空间布局（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"
)

该声明包含由 bigPsl.as 来自UCSC的AutoSQL文件。目标序列（通常是序列比对的染色体）存储在 targets 属性在bigBed格式中，仅存储每个目标的标识符和大小。在这个例子中，只有一条染色体：

>>> alignments.targets
[SeqRecord(seq=Seq(None, length=198295559), id='chr3', name='<unknown name>', description='<unknown description>', dbxrefs=[])]

迭代路线会为每条线提供一个路线对象：

>>> 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 ：报告与目标或查询序列中的RST字符（代表未知核苷）的比对 nCount 字段而不是在 matches , misMatches ，或者 repMatches 领域的默认值为 "N" .

见章节图案空间布局（PSL）解释如何获得匹配、错配、重复区匹配以及与未知核苷匹配的数量。

进一步的可选论点是 blockSize （默认值为256），以及 itemsPerSlot （默认值为512）。请参阅UCSC的文档 bedToBigBed 描述这些论点的程序。搜索 bigPsl 文件可以使用更快 compress=False 和 itemsPerSlot=1 创建bigPsl文件时。

多重对齐格式（MAF）

APM（多重比对格式）文件以人类可读的格式存储一系列多重序列比对。APM文件通常用于存储基因组彼此的对齐。文件 ucsc_test.maf 在 Tests/MAF Biopython发行版的收件箱是简单的APM文件的一个示例：

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

文件头中显示的信息（以“#'#'开始的轨道行和后续行）存储在 metadata 属性 alignments 对象：

>>> 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 我们得到一个 Alignment APM文件中每个路线块的对象：

>>> 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 ''（序列）行外，APM文件还可能包含包含有关该块之前和之后基因组序列的信息的“' i '”行，包含有关比对空部分的信息，以及显示每个比对碱基质量的“''行。这是包括此类线的对齐块的示例：

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

这是文件中的第10个对齐块 ucsc_mm9_chr10.maf （可在 Tests/MAF Biopython发行版的副本）：

>>> 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 documentation 获取这些字段和状态字符的完整描述。

“' q '”行显示序列质量，该质量存储在的“' quality '”字典项下 annotations 相应序列的属性：

>>> alignment.sequences[2].annotations["quality"]
'9999999999999999999999999999999999999999999999'
>>> alignment.sequences[4].annotations["quality"]
'9999969989999999999999998699989997'

“e”线显示了有关在该比对区块之前和之后具有连续序列的物种的信息，但在该比对区块中没有比对核苷。这存储在 alignment.annotations 词典：

>>> 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')]

例如，这表明在相对链上从位置158040939到158048983插入了非对齐碱基（“' I '”） ponAbe2.chr6 genomic sequence. Again, see the UCSC documentation 了解“e”行的完整定义。

要以APM格式打印对齐，您可以使用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

要编写完整的APM文件，请使用 Bio.Align.write(alignments, "myfilename.maf", fmt="maf") ，在哪里 myfilename.maf 是输出MAF文件的名称。或者，您可以使用（文本）流进行输出。文件头信息将从 metadata 属性 alignments object.如果您要从头开始创建路线，则可以使用 Alignments （复数）类创建类似列表 alignments 对象（请参阅部分 Alignments课程）并给它一个 metadata 属性

bigMaf

bigMaf文件是BED 3 +1格式的bigBed文件，包括3个必需的BED字段和一个自定义字段，该字段将APM对齐块存储为字符串，从而创建APM文件的索引二进制版本（请参阅部分多重对齐格式（MAF） ).关联的AutoSql文件 bigMaf.as 由UCSC提供。要创建bigMaf文件，您可以使用 mafToBigMaf 和 bedToBigBed 来自UCSC的节目。或者您可以通过调用Bio.Align. writer函数来使用Biopython fmt="bigmaf" .虽然这两种方法应该会产生相同的bigMaf文件，但UCSC工具速度要快得多，而且可能更可靠，因为它是黄金标准。由于bigMaf文件是bigBed文件，因此它们带有内置索引，允许您快速搜索参考基因组的特定区域。

文件 ucsc_test.bb 在 Tests/MAF Biopython发行版的收件箱是bigMaf文件的一个示例。该文件相当于APM文件 ucsc_test.maf (see部分多重对齐格式（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文件不存储在APM文件中找到的标头信息，但它确实定义了参考基因组。相应的 SeqRecord 存储在 targets 属性 alignments 对象：

>>> alignments.reference
'hg16'
>>> alignments.targets
[SeqRecord(seq=Seq(None, length=158545518), id='hg16.chr7', ...)]

通过迭代 alignments 我们得到一个 Alignment bigMaf文件中每个对齐块的对象：

>>> 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”行中的信息的存储方式与相应的APM文件相同（请参阅部分多重对齐格式（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')]

To write a complete bigMaf file, use Bio.Align.write(alignments, "myfilename.bb", fmt="bigMaf"), where myfilename.bb is the name of the output bigMaf file. Alternatively, you can use a (binary) stream for output. If you are creating the alignments from scratch, you can use the Alignments (plural) class to create a list-like alignments object (see Section Alignments课程) and give it a targets attribute. The latter must be a list of SeqRecord objects for the chromosomes for the reference species in the order in which they appear in the alignments. Alternatively, you can use the targets keyword argument when calling Bio.Align.write. The id of each SeqRecord must be of the form reference.chromosome, where reference refers to the reference species. Bio.Align.write has the additional keyword argument compress (True by default) specifying whether the data should be compressed using zlib. Further optional arguments are blockSize (default value is 256), and itemsPerSlot (default value is 512). See the documentation of UCSC's bedToBigBed program for a description of these arguments.

由于bigMaf文件是bigBed文件的特殊情况，因此您可以使用 search 方法对 alignments 目标是找到与参考物种特定区域的对齐。例如，我们可以寻找chr10上10号染色体3018000和3019000位置之间的区域：

>>> 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文件时。

链文件格式

链文件描述了两个核苷酸序列之间的成对比对，允许两个序列中存在缺口。只有每个对齐序列的长度和间隙长度存储在链文件中;序列本身不存储。链文件通常用于存储两个基因组组装版本之间的比对，允许将一个基因组组装版本的比对转移到另一个基因组组装。这是链文件的示例（可用为 psl_34_001.chain 在 Tests/Blat Biopython发行版的副本）：

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

此文件是通过运行UCSC的 pslToChain PSL文件上的程序 psl_34_001.psl .根据链文件格式规范，每个链块后面应该有一行，但有些工具（包括 pslToChain ）显然不遵守这个规则。

要解析此文件，请使用

>>> from Bio import Align
>>> alignments = Align.parse("psl_34_001.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)
...

检查对齐分数和链ID（分别是每个链块标题行中的第一个和最后一个数字）以确认我们获得了第七个对齐：

>>> alignment.score
41.0
>>> alignment.annotations["id"]
'7'

我们可以以链文件格式打印对齐方式。对齐坐标与链块中的信息一致，对齐部分为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