第20章高级

20.1 解析器的设计

过去很多Biopython解析器都是根据面向事件设计出来的，包括Scanner和Consumer。

Scanners是将输入的数据源进行逐行分析，只要识别出数据中的信息就会发送一个事件。例如，如果数据中包含生物名的信息，scanner只要读到某行包含名称信息时就会产生一个“生物名”的事件。

Consumers是用来接收Scanners所发出事件的对象。接着上面的例子，当consumer收到了“生物名”事件，在当前应用程序中无论以何种方式都会运行。

这是一个非常灵活的构架，如果你想要将一个文件解析成多种其他格式的话，这会是很有优势的。例如，“Bio.GenBank”模块可以运用这种方式构建“SeqRecord”或者其他独特的文件格式记录对象。

最近，很多添加了“Bio.SeqIO”和“Bio.AlignIO”的解析器在使用一种更为简单的方法，但是只能产生单一形式的文件格式（分别是“SeqRecord”和“MultipleSeqAlignment”）。在某些情况，“Bio.SeqIO”解析器实际上包含了另一种Biopython解析器。例如，“Bio.SwissProt”解析器产生了特定的SwissProt格式对象，又转换成了“SeqRecord”格式对象。

20.2 替换矩阵

20.2.1 SubsMat

这个模块提供了一个类和一些固定的方法来产生替换矩阵，类似于BLOSUM或者PAM矩阵，但是是基于用户提供的数据。此外，你还可以从已建立的替换矩阵集合MatrixInfo.py中选择一个矩阵。 “SeqMat”类来自于一个dictionary库:

class SeqMat(dict)

这个dictionary的格式是 {(i1,j1):n1, (i1,j2):n2,...,(ik,jk):nk} i和j是字母，而n是一个值。

属性
1. self.alphabet: Bio.Alphabet中定义的一个类
2. self.ab_list: 排列好的字母列表。主要是内部原因而需要。
方法
1. ```
__init__(self,data=None,alphabet=None, mat_name='', build_later=0):
```
  1. data: 可以是一个dictionary，也可以是另一个SeqMat实例。
  2. alphabet: 一个Bio.Alphabet的实例。如果没有提供，可以从数据构建一个alphabet。
  3. mat_name: 矩阵名，例如 "BLOSUM62" 或者 "PAM250"
  4. build_later: 默认值为假。如果是真，用户应该只提供alphabet和空的dictionary。
    
    如果想要之后再构建矩阵，这样会跳过alphabet大小和矩阵大小的检查。
2. ```
entropy(self,obs_freq_mat)
```
  1. obs_freq_mat: 一个观测频率矩阵。基于“obs_freq_mat”的频率返回矩阵的熵值。矩阵实例应该是LO或者SUBS。
3. ```
sum(self)
```
  计算矩阵的alphabet中每个字母值的总和，返回值是以dictionary的形式 {i1: s1, i2: s2,...,in:sn}, 其中:
  - i: 一个字母;
  - s: 半矩阵中某个字母值的总和;
  - n: alphabet中字母的个数。
4. ```
print_mat(self,f,format="%4d",bottomformat="%4s",alphabet=None)
```
  将矩阵打印到文件句柄f。“format”是矩阵值的格式；“bottomformat”是底部行的格式，包括矩阵字母。下面是一个3字母矩阵的例子：
```
A 23
B 12 34
C 7  22  27
  A   B   C
```
  ``alphabet``可选的自变量是alphabet中所有字符的一个字符串。
  
  如果用户提供了数据，则轴上字母的顺序是根据字符串中的顺序，而不是字母表的顺序。
用法

安排下面这个部分是因为大多数的读者希望能够知道如何产生一个对数机率矩阵。当然，也有人研究产生的中间过渡矩阵。但是大部分的人只是想要一个对数机率矩阵，仅此而已。
1. 产生一个可接受的代替矩阵
  
  首先，你应该从数据中产生出一个可接受代替矩阵（ARM）。
  
  ARM中的数值是根据你的数据中替换的个数决定的。数据可以是一对或者多对的序列比对结果。例如，丙氨酸被半胱氨酸替换了10次，而半胱氨酸被丙氨酸替换了12次，其相对应的ARM为：
```
('A','C'): 10, ('C','A'): 12
```
  由于顺序并不重要，用户也可以只用一个输入:
```
('A','C'): 22
```
  一个SeqMat实例的初始化可以用全矩阵（第一种计数方法：10,12），也可以用半矩阵（后一种方法，22）。
  
  一个蛋白字母全矩阵的大小应该是20x20 = 400。而一个这样的半矩阵大小是20x20/2 + 20/2 = 210。这是因为相同字母的输入并没有改变（矩阵的对角线）。如果一个大小为N的alphabet：
  1. 全矩阵大小: N*N
  2. 半矩阵大小: N(N+1)/2
  如果已经通过了全矩阵，SeqMat的构造器会自动产生半矩阵。
  
  如果通过了半矩阵，则键的字母将按照字母表顺序排列('A','C')，而不是('C','A')。
  
  讲到这里，如果你想知道的仅仅只是怎样产生一个对数机率矩阵的话，请直接看用法示例那个章节。
  
  对于想要更加深入地知道核苷酸/氨基酸频率数据的读者，接下来要讲的是内部功能的细节。
2. 产生观察频率矩阵(OFM)
  
  用法:
```
OFM = SubsMat._build_obs_freq_mat(ARM)
```
  观察频率矩阵是由可接受代替矩阵产生的，只是将替换的个数换成了替换频率。
3. 产生一个期望频率矩阵(EFM)
  
  用法:
```
EFM = SubsMat._build_exp_freq_mat(OFM,exp_freq_table)
```
  1. exp_freq_table: 应该是一个FreqTable的实例。
    
    更加详细的信息可以查看章节20.2.2 简单地说，期望频率表表示字母表中每个元素显示的频率。
    
    这个表相当于一个dictionary，字母是键，字母对应的频率是值。总值是1。
  期望频率表可以（理论上说应该）从观察频率矩阵得到。
  
  所以大多数情况你可以用下面的代码产生``exp_freq_table``:
```
>>> exp_freq_table = SubsMat._exp_freq_table_from_obs_freq(OFM)
>>> EFM = SubsMat._build_exp_freq_mat(OFM,exp_freq_table)
```
  如果你想的话，也可以使用自己提供的``exp_freq_table``
4. 产生一个替换频率矩阵(SFM)
  
  用法:
```
SFM = SubsMat._build_subs_mat(OFM,EFM)
```
  使用观察频率矩阵(OFM)和期望频率矩阵(EFM)。
  
  得到相应值的除法结果。
5. 产生一个对数机率矩阵(LOM)
  
  用法:
```
LOM=SubsMat._build_log_odds_mat(SFM[,logbase=10,factor=10.0,round_digit=1])
```
  1. 使用一个替换频率矩阵(SFM)。
  2. logbase: 用来产生对数机率值的对数的底。
  3. factor: 因子是对数机率值的乘数。
    
    每个数通过log(LOM[key])*factor产生，如果需要还可以通过``round_digit``四舍五入到小数点后。
用法示例

正如大部分的人都想用最简单的方法产生对数机率矩阵，SubsMat提供了一个函数可以完成所有的事情：
```
make_log_odds_matrix(acc_rep_mat,exp_freq_table=None,logbase=10,
                      factor=10.0,round_digit=0):
```
1. acc_rep_mat: 用户提供可接受代替矩阵
2. exp_freq_table: 期望频率表。如果提供了就使用期望频率表，否则就从``acc_rep_mat``产生。
3. logbase: 对数机率矩阵的对数的底。默认底为10。
4. round_digit: 四舍五入的小数点位数。默认为0。

20.2.2 FreqTable

FreqTable.FreqTable(UserDict.UserDict)

属性:
1. alphabet: 一个Bio.Alphabet的实例。
2. data: 频率dictionary
3. count: 计数dictionary(如果有计数的话)。
功能:
1. read_count(f): 从字符串f读入一个计数文件。然后将其转换成频率。
2. read_freq(f): 从字符串f读入一个频率数据文件。当然，我们不用计数，我们感兴趣的是字母频率。
用法示例: 文件中有残基的个数，用空格分格，形式如下（以3个字母为例）：
```
A   35
B   65
C   100
```
用``FreqTable.read_count(file_handle)``方法读入。

一个等价的频率文件:
```
A  0.175
B  0.325
C  0.5
```
相反地，残基频率或者计数也可以作为dictionary输入。一个计数dictionary的例子（3个字母）：
```
{'A': 35, 'B': 65, 'C': 100}
```
这也意味着'C'的频率是0.5，'B'的频率是0.325，'A'的频率是0.175。A、B、C的总和为200。

一个相同数据的频率dictionary如下：
```
{'A': 0.175, 'B': 0.325, 'C': 0.5}
```
总和为1。

当传入一个dictionary数据作为参数，应该指出这是一个计数还是频率的dictionary。因此FreqTable类的构建需要两个参数： dictionary本身和FreqTable.COUNT或者FreqTable.FREQ，分别代表计数或者频率。

读入期望的计数，readCount会产生频率。下面的任意一个都可以用来产生频率表（ftab）：
```
>>> from SubsMat import *
>>> ftab = FreqTable.FreqTable(my_frequency_dictionary,FreqTable.FREQ)
>>> ftab = FreqTable.FreqTable(my_count_dictionary,FreqTable.COUNT)
>>> ftab = FreqTable.read_count(open('myCountFile'))
>>> ftab = FreqTable.read_frequency(open('myFrequencyFile'))
```

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

chr20.rst

chr20.rst

第20章高级

20.1 解析器的设计

20.2 替换矩阵

20.2.1 SubsMat

20.2.2 FreqTable

Files

chr20.rst

Latest commit

History

chr20.rst

File metadata and controls

第20章 高级

20.1 解析器的设计

20.2 替换矩阵

20.2.1 SubsMat

20.2.2 FreqTable

第20章高级