蛋白质数据库(PDB)生物信息学数据库是世界上最大的实验确定的蛋白质结构储存库, nucleic acids, 以及复杂的装配. 所有数据都是通过x射线等实验方法收集的, spectroscopy, crystallography, NMR, etc.

本文解释了如何从PDB中提取、过滤和清理数据. 这反过来又支持本文中解释的分析类型 蛋白质二硫键在不同生命领域的发生:蛋白质数据库中蛋白质的比较, published in 蛋白质工程，设计与选择，第27卷，第3期，2014年3月1日，页. 65–72.

PDB有许多具有不同分辨率、方法、突变等的重复结构. 用相同或相似的蛋白质做实验，在任何组分析中都可能产生偏差, 所以我们需要从任何一组重复序列中选择正确的结构. 为此，我们需要使用一组非冗余(NR)蛋白质.

为了实现规范化，我建议下载 化合物词典 用于将原子名称导入到使用3NF或使用星型模式和维度建模的数据库中. (I’ve also used DSSP 帮助消除有问题的结构. 我不会在本文中讨论这一点，但请注意，我没有使用任何其他DSSP特性.)

本研究中使用的数据包含含有至少一个来自不同物种的二硫键的单单位蛋白质. 执行分析, 所有二硫键首先分为连续键和非连续键, by domain (archaea, prokaryote, viral, eukaryote, or other), and also by length.

Output

准备好输入 R， SPSS或其他工具，分析师将需要将数据放在 database 表具有这种结构:

Materials and Methods

为了实现这一目标，第一步是 gather data from rcsb.org. 该网站包含各种格式的可下载PDB实验结构.

尽管数据以多种格式存储, in this example, 只使用格式化的固定空格分隔文本格式(PDB). PDB文本格式的另一种替代方案是它的XML版本, PDBML, 但它有时包含不正确的原子位置命名项, 哪些会给数据分析带来问题. 也可以使用较旧的mmCIF和其他格式, 但本文不会对它们进行解释.

The PDB Format

PDB格式是一种碎片化的固定宽度文本格式，可以很容易地被SQL查询解析, Java plugins, or Perl modules, for example. 文件容器中的每种数据类型都表示为以适当的标记开头的一行，我们将在以下小节中详细介绍每种标记类型. 行长度小于等于80个字符, 标签需要六个或更少的字符加上一个或多个空格，总共需要八个字节. 也有标记和数据之间没有空格的情况，通常为 CONECT tags.

TITLE

The TITLE 标签将一行标记为实验标题的一部分, 包含分子名称和其他相关数据，如插入, mutation, 或者是特定氨基酸的缺失.

12345678901234567890123456789012345678901234567890123456789012345678901234567890
标题为一种含有二硫化物的双二硫衍生物            
标题2 13-33被两个α -氨基丁酸取代，nmr, 30              
TITLE    3 STRUCTURES                                                           

在有多条线到a的情况下 TITLE record, 然后标题必须连接起来, 按连续号排序, which is placed, right-aligned, 在第9和第10字节上(取决于这些行的数量).

ATOM

The data stored in ATOM 线是实验中每个原子的坐标数据. 有时一个实验有插入、突变、替代位置或多个模型. 这导致同一个原子重复多次. 选择正确的原子将在后面解释.

12345678901234567890123456789012345678901234567890123456789012345678901234567890
原子390和原子26 -1.120  -2.842   4.624  1.00  0.00           N  
原子391可以达到26 -0.334  -2.509   3.469  1.00  0.00           C  
Atom 392 c gly 26 0.682  -1.548   3.972  1.00  0.00           C  
原子393，原子26，原子26.420  -0.786   4.898  1.00  0.00           O  
原子394的高度为26 -0.832  -2.438   5.489  1.00  0.00           H  
原子395的直径为26英寸.163  -3.399   3.111  1.00  0.00           H  
原子396的直径为26 -0.955  -2.006   2.739  1.00  0.00           H  

上面的例子取自实验 1BAH. 第一列标记记录的类型，第二列是原子的序列号. 结构中的每个原子都有自己的序号.

序列号旁边是原子位置标签，占用四个字节. 从那个原子位置开始, 提取这种元素的化学符号是可能的, 哪一个不总是存在于记录数据的单独列中.

原子名之后是三个字母的剩余代码. 在蛋白质的情况下，残基对应于氨基酸. 接下来，用一个字母对链进行编码. By chain 我们指的是单链氨基酸, with or without gaps, although sometimes ligands can be assigned to a chain; this case is detectable through very large gaps in an amino acid sequence, 下一列是哪个. 有时可以扫描包含突变的相同结构, 在这种情况下，插入代码在序列列后面的额外列中可用. 插入代码包含一个字母，以帮助区分哪些残留物受到影响.

接下来的三列是每个原子的空间坐标，单位是埃(Å)。. 这些坐标旁边是占用列, 也就是说原子在那个位置的概率是多少, 通常的等级是0到1.

倒数第二列是温度因子, 它携带着晶体无序性的信息, measured in Ų. 大于60Ų的值表示无序，小于30Ų的值表示自信. 它并不总是出现在PDB文件中，因为它取决于实验方法.

下一列——符号和电荷——通常是缺失的. 如上所述，化学符号可以从原子位置列中得到. 当电荷存在时，符号的后缀是一个整数，后面跟着 + or -, e.g. N1+.

TER

这标志着链条的末端. 即使没有这条线，也很容易分辨出链条的终点. Thus, often the TER line is missing.

MODEL and ENDMDL

A MODEL 线条标志着结构模型的起始位置，它包含模型的序列号.

在该模型中的所有原子行之后，它以 ENDMDL line.

SSBOND

这些线含有半胱氨酸氨基酸之间的二硫键. 二硫键可以存在于其他残留物类型中, 但本文只分析氨基酸, 所以只有半胱氨酸被包括在内. 下面的示例来自代码实验 132L:

12345678901234567890123456789012345678901234567890123456789012345678901234567890
Ssbond 1 cys a 6 cys a 127 1555 1555 2.01
2 . a . 1555 a . 1555 a . 1555.05
3 . a . 64 . a . 80 . 1555.02
4天a 76天a 94 1555 1555.02

在这个例子中，有四个二硫键在文件中被标记，它们的序列号在第二列中. 所有这些键都使用半胱氨酸(第3列和第6列)，它们都是链状的 A (columns 4 and 7). 每条链后面都有一个残基序列号，表示该键在肽链中的位置. 插入码在每个残基序列旁边, 但在这个例子中，它们不存在因为没有氨基酸插入到那个区域. 最后一列之前的两列保留给 symmetry operations，最后一列是硫原子之间的距离，用Å来测量.

让我们花点时间来了解一下这些数据的背景.

下面的图片是用rcsb拍摄的.org NGL查看器，展示实验结构 132L. 特别是，它们显示了一个没有配体的蛋白质. 第一张图使用棍子表示, 用CPK着色，硫磺和它们的化学键呈黄色. v形硫连接代表蛋氨酸连接, 而z形连接是半胱氨酸之间的二硫键.

在下一张图中，一种简化的蛋白质可视化方法叫做 骨干可视化 是用氨基酸来表示颜色的，其中半胱氨酸是黄色的. 它代表相同的蛋白质, 不包括其侧链, 在这个例子中，只有一部分的肽基, the protein backbone. 它由三个原子组成:n端，c端和c端. 这张图没有显示二硫键, 但为了显示蛋白质的空间排列，它被简化了:

管道是通过将肽键原子与c原子连接而形成的. 半胱氨酸的颜色与CPK着色法中硫的颜色相同. 当半胱氨酸足够接近时, 它们的硫形成二硫键, 这加强了结构. 否则蛋白质就会过度结合, 而且它的结构在高温下不太稳定.

CONECT

这些记录用于标记原子之间的连接. 有时这些标记根本不存在，而其他时候则填充了所有数据. 在分析二硫键的情况下，这部分是有用的，但不是必需的. That’s because, in this project, 通过测量距离来添加无标记的键, 这就是开销, 还需要检查.

SOURCE

这些记录包含了提取分子的来源生物的信息. 它们包含子记录，便于在分类法中定位, 并且具有与标题记录相同的多行结构.

SOURCE    MOL_ID: 1;                                  
来源2生物科学:冈比亚按蚊;    
来源3生物常见:非洲疟蚊; 
来源4 organm_taxid: 7165;                      
源5基因:gst1-6;                              
来源6表达系统:大肠杆菌;       
来源7 expression_system_taxid: 562;              
来源8 expression_system_strain: bl21 (de3) plyss;  
来源9 expression_system_vector_type:质粒;    
来源10 expression_system_plasmid: pxaggst1-6      

NR Format

这是一个非冗余(NR)链PDB集列表. 它的快照可以在 ftp.ncbi.nih.gov/mmdb/nrtable/. 其目的是为了避免由于蛋白质相似性造成的不必要的偏差. 通过比较所有PDB结构，NR有三个具有不同身份p值水平的集合. 结果被添加到文本文件中，稍后将对此进行解释. 本项目并不需要所有的列，因此将只解释重要的列.

如前所述，前两列包含唯一的PDB实验代码和链标识符 ATOM records above. Columns 6, 9, 和C包含关于p值代表性的信息, 用哪一种方法计算序列的相似度 BLAST. If that value is zero, then it is not accepted to be part of a set; if the value is 1, then it is. 上述列表示p值截止值为10e-7的集合的可接受性, 10e-40, and 10e-80, respectively. 只有p值截止值为10e-7的集合才会用于分析.

最后一列包含有关结构可接受性的信息，其中 a is acceptable and n is not.

#---------------------------------------------------------------------------------------------------------------------------
# 1 2 3 4 5 6 7 8 9 a b c d d d d d I j k I m n p q
#---------------------------------------------------------------------------------------------------------------------------
3f8v a 69715 1 1 1 1 1 1 1 1 1 1 9427 1 1 0.00   0.00   0.00   0.00  1.08   1   6   5   164 X a
3dke x 68132 1 2 0 1 1 2 0 39139 1 1 0 0.00   0.00   0.00   0.00  1.25   1  11   7   164 X a
3hh3 a 77317 1 3 0 1 3 0 1 3 0 90 1 0 0.00   0.00   0.00   0.00  1.25   1   5   4   164 X a
3hh5 a 77319 1 4 0 1 4 0 1 4 0 90 2 0 0.00   0.00   0.00   0.00  1.25   1   4   4   164 X a

数据库构建和数据解析

现在我们已经了解了我们要处理的问题和我们需要做的事情，让我们开始吧.

Downloading Data

此分析的所有数据都可以在以下三个地址找到:

• ftp.wwpdb.org/pub/pdb/data/structures/divided/pdb/
• ftp.ncbi.nih.gov/mmdb/nrtable/
• ftp.ncbi.nih.gov / pub /分类/ taxdmp.zip

前两个链接包含一个存档列表. 应该使用每个文件的最新存档，以避免由于缺乏解决方案或质量而产生的问题. 第三个链接直接包含最新的分类法存档.

Parsing Data

通常PDB文件的解析是由Java、Perl或Python中的插件或模块完成的. 在这项研究中, 我编写了一个自定义Perl应用程序，没有使用预先编写的pdb解析模块. 其原因是在解析大量数据时, in my experience, 使用实验数据最常见的问题是数据中的错误. 有时候坐标会有误差, distances, line lengths, 注释出现在不该出现的地方, etc.

处理这个问题最有效的方法是最初将所有内容以原始文本的形式存储在数据库中. 编写通用解析器是为了处理完全符合规范的理想数据. But in practice, data is not ideal, 这将在过滤一节中解释，您将在其中找到Perl导入脚本.

Database Construction

在构建数据库时，请注意该数据库是为处理数据而构建的. 之后的分析将在SPSS或R中完成. 出于我们的目的，建议至少使用版本8的PostgreSQL.4.

表结构直接从下载的文件中复制，只做了一些小的修改. In this case, 记录的数量太少了，不值得花时间进行规范化. 如前所述，此数据库仅是一次性使用的，而这些表则不是 为在网站上提供服务而构建的它们只是用来处理数据的. Once that is finished, they can be dropped, 或备份为补充数据, 也许是因为其他研究人员重复了这个过程.

In this case, 最终的结果将是一个表格，然后可以导出到一个文件中，以便与SPSS或R等统计工具一起使用.

Tables

Data extraction from ATOM 必须连接到记录 HEADER or TITLE records. 数据层次结构如下图所示.

由于这张图是第三范式(3NF)中数据库的简化表示。, 对于我们的目的来说，它包含了太多的开销. 目的:计算原子间的距离，用于二硫键检测, 我们需要进行连接. In this case, 我们将有一个表连接到它自己两次, 也分别连接到二级和一级结构上两次, 这是一个非常缓慢的过程. 因为不是每个分析都需要二级结构信息, 如果您需要重用数据或分析更大量的二硫键，则提出另一种模式:

二硫键不像其他共价键那么常见, 因此不需要仓库模型, 虽然它可以被使用. 下面的星型模式和维度建模将花费太多时间来开发, 并且会使查询更加复杂:

在需要处理所有键的情况下，我推荐使用星型模式.

(否则就不需要了，因为二硫键不像其他键那么常见. 在这项工作中, 二硫键的数目接近30个,000, 在3NF中哪个可能足够快, 但我决定通过一个非规范化的表来处理它, 所以这里没有画出来.)

所有共价键的预期总数至少是三级结构中原子数的两倍, 在这种情况下，3NF会非常慢, 因此，需要使用星型模式形式进行反规范化. In that schema, 有些表有两个外键检查, 这是因为两个原子之间形成了键, 所以每个原子都需要有自己的 primary_structure_id, atom_name_id and residue_id.

有两种方法来填充 d_atom_name 维度表:来自数据, 从另一个来源, 我之前提到的化学成分词典. 其格式类似于PDB格式:Only RESIDUE and CONECT lines are useful. This is because RESIDUE的第一列包含剩余的三个字母代码，并且 CONECT 包含原子及其连接的名称，它们也是原子名称. So from this file, 我们可以解析所有原子名称并将它们包含在数据库中, 尽管我建议您允许数据库包含未列出的原子名称.

RESIDUE   PRO     17
连接3 ca CD h   
连接c4c4ha  
连接c和c 
连接c   
连接cb4 ca cghb2 hb3 
连接cg4 cb CD hg2 hg3 
连接cd4和cghd2 hd3 
连接oxt 2c HXT 
连接h和n   
连接1个ca  
连接hb21 - 1cb  
连接hb31 1cb  
连接hg2 - 1cg  
连接hg3 - 1cg  
连接hd2 1 CD  
连接hd31光盘  
连接hxt1oxt 
END   
Het pro 17
南脯氨酸
公式为c5h9n1o2

在这个项目中，编码速度比执行速度和存储消耗更重要. 我决定不去正常化——毕竟, 我们的目标是生成一个包含介绍中提到的列的表.

在这一部分中，只会解释最重要的表格.

The main tables are:

• proteins:实验名称和实验代码表.
• ps:主结构表，它将包含 sequence, chain_id, and code.
• ts:从原始数据中提取并转换成的三级/四级结构表 ATOM record format. 这将用作staging表，并且可以在提取后删除. Ligands are excluded.
• sources获得实验数据的生物体列表.
• tax_names, taxonomy_path, taxonomy_paths:来自NCBI分类法数据库的林奈分类法名称, 用于从列表中列出的生物体中获取分类路径 sources.
• nr:从NR集合中提取的NCBI非冗余蛋白列表.
• pdb_ssbond:给定PDB文件中的二硫键列表.

数据的过滤和处理

数据是从RCSB PDB存储库的快照中检索的.

每个文件被导入到单个表中 raw_pdb 在我们的PostgreSQL数据库中使用Perl脚本. 该脚本使用每个块10,000个插入的事务.

The structure of raw_pdb is this:

导入脚本看起来像这样:

#!/usr/bin/perl -w

use FindBin '$Bin';
use DBI;

$dbName = 'bioinf'; 
$dbLogin = 'ezop'; 
$dbPass = 'XYZ';

$conn = DBI->connect("DBI:Pg:database=$dbName;host=localhost", "$dbLogin", "$dbPass", {'RaiseError' => 1, 'AutoCommit' => 0});

die "./pdb_lines_unos.pl " if not defined($ARGV[0]);

$recordCount = 0;
$table = $ARGV[0];
$fName = $ARGV[1];
打开F， "zcat $fName|";
while () {
    chomp;
    $linija = $_;
    $recordCount += 1;
    $sql = "insert into $table (code, line_num, line_cont) values "?, ?, ?)";

    $conn->do($sql, undef, $fName, $recordCount, $linija);
    $conn->commit() if ($recordCount%10000 == 0);
}
close F;
$conn->commit();

1;


在导入行之后，将使用下面定义的函数对它们进行解析.

From raw_pdb 数据，我们生成表格 ts, ps, proteins, sources, sources_organela, and ss_bond 通过解析相应的记录.

The ps 该表有三个重要列: chain, length, and sequence. 利用c - α原子对每条链和按残基序列排序的每个残基生成蛋白质序列, 只取第一个插入和第一个替代位置. chain is taken from the TS.chain column, and length 仅仅是预先计算的长度 sequence string. 因为本文只分析单链和链内连接, 在这里，我们的分析跳过了多链蛋白.

Within SSBOND 记录，所有的二硫键都存储在 pdb_ssbond 表，它继承自 pdb_ssbond_extended table.  pdb_ssbond looks like this:

Column
Type
Nullable
Default
Description
id
integer
not null
nextval(“pdb_ssbond_id_seq”::regclass)
 
code
character(4)
 
 
four-letter code
ser_num
integer
 
 
ssbond的序列号
residue1
character(3)
 
 
first residue in bond
chain_id1
character(1)
 
 
first chain in bond
res_seq1
integer
 
 
第一余数的序号
i_code1
character(1)
 
 
第一残基的插入码
residue2
character(3)
 
 
键合中的第二残基
chain_id2
character(1)
 
 
second chain in bond
res_seq2
integer
 
 
第二余数的序号
i_code2
character(1)
 
 
第二残基的插入码
sym1
character(6)
 
 
第一对称算子
sym2
character(6)
 
 
第二对称算子
dist
numeric(5,2)
 
 
原子间距离
is_reactive
boolean
not null
false
反应性标记(待定)
is_consecutive
boolean
not null
true
连续债券保证金(待定)
rep7
boolean
not null
false
set-7结构标记(待设置)
rep40
boolean
not null
false
set-40结构标记(待设置)
rep80
boolean
not null
false
set-80结构标记(待设置)
is_from_pdb
boolean
 
true
从PDB中获取作为源(待设置)

"pdb_ssbond_pkey"主键，btree (id)
"ndxcode1" btree (code, chain_id1, res_seq1)
"ndxcode2" btree (code, chain_id2, res_seq2)

选择count(1)作为CNT
    ， std_dev (dist)为std_dev  
    ， avg(dist)为avg_val
    ， -stddev(dist) + avg(dist) as left1
    ， stddev(dist) + avg(dist)为right1
    ， -2 * stddev(dist) + avg(dist) as left2
    ， 2 * stddev(dist) + avg(dist) = right2
    ， -3 * stddev(dist) + avg(dist) as left3
    ， 3 * stddev(dist) + avg(dist) = right3
    ， -4 * stddev(dist) + avg(dist) as left4
    ， 4 * stddev(dist) + avg(dist) = right4
    ,         min(dist)
    ,         max(dist)
from pdb_ssbond_tmp
where dist > 0
    and dist < 20;

 创建或替换函数pdb_ssbond_clean2 
    clean_icodes boolean,
    clean_altloc boolean,
    mark_reactive布尔,
    mark_consecutive布尔)
RETURNS void AS $$
declare _res integer;
BEGIN  
    从pdb_ssbond b中删除 
    where exists (
        select a.id
        from pdb_ssbond a
        where a.code = b.code
        and a.id > b.id
        and (
            (   a.chain_id1 = b.chain_id1 and a.res_seq1 = b.res_seq1
            and a.chain_id2 = b.chain_id2 and a.res_seq2 = b.res_seq2)
            or
            (   a.chain_id1 = b.chain_id2 and a.res_seq1 = b.res_seq2
            and a.chain_id2 = b.chain_id1 and a.res_seq2 = b.res_seq1)
        )
    ) ; 

    delete
    from pdb_ssbond
    chain_id1 = chain_id2
    res_seq1 = res_seq2
    和i_code1 = i_code2;
        
 
         
    update pdb_ssbond 
    设置is_continuous = true
    ， is_reactive = false;
        
    从pdb_ssbond中删除 
    不等于1的地方.6175344456264 and 2.48801947642267;
     
          
    
    if clean_icodes then  
        从pdb_ssbond中删除 
        where i_code1 not in ("， ''， 'A') 
        or i_code2 not in ("， ''， 'A');
    end if;
    
    if clean_altloc then  
        从pdb_ssbond a中删除
        where exists (
            select 1
            from pdb_ssbond  b
            where   b.code = a.code
                and b.chain_id1 = a.chain_id1
                and b.res_seq1 = a.res_seq1
                and b.i_code1 = a.i_code1
                and b.ser_num 

Select * from sources_organela where code = '1rav';

Select * from taxonomy_path order by length(path) limit 20 offset 2000

创建或替换函数len2class
    RETURNS integer AS
$BODY$
BEGIN
return 
    case 
        when len <= 100                then 1
        when len >  100 and len <= 200 then 2
        when len >  200 and len <= 300 then 3
        when len >  300 and len <= 400 then 4
                                       else 5
    end;
END;
$BODY$
    语言plpgsql VOLATILE
    COST 100;

 SELECT p.code,
    p.title,
    p.ss_bonds,
    p.ssbonds_overlap,
    p.intra_count,
    p.sci_name_src,
    p.len,
    p.tax_path,
    p.pfam_families,
    p.src_class,
    ( SELECT s.id
           FROM src_classes s
          WHERE s.src_class::text = p.src_class::text) AS src_class_id，
    p.len_class7,
    ( SELECT s.val
           FROM sources_organela s
          WHERE s.code = p.code::bpchar AND s.tag::text = 'EXPRESSION_SYSTEM_CELLULAR_LOCATION'::text)作为EXPRESSION_SYSTEM_CELLULAR_LOCATION，
    ( SELECT s.val
           FROM sources_organela s
          WHERE s.code = p.code::bpchar AND s.标签::text = 'CELLULAR_LOCATION'::text)
    ps.sequence,
    ps.uniprot_code,
    ps.accession_code,
    ps.cc,
    ps.seq_uniprot,
    ps.chain_id
   FROM proteins p
     JOIN nr n ON n.code::text = p.code::text AND n.rep7 = true
     JOIN ps ps ON ps.code::text = n.code::text AND ps.chain_id = n.chain_id AND ps.het = false
  WHERE p.is_excluded = false.chain_cnt = 1 AND p.is_set7 = true AND p.reactive_cnt7 = 0
  ORDER BY p.code;