ジスルフィド結合研究のためのバイオインフォマティクスデータベースの開発

Protein Data Bank（PDB）バイオインフォマティクスデータベースは、実験的に決定されたタンパク質、核酸、および複雑なアセンブリの構造の世界最大のリポジトリです。 すべてのデータは、X線、分光法、結晶学、NMRなどの実験方法を使用して収集されます。

この記事では、PDBからデータを抽出、フィルタリング、およびクリーンアップする方法について説明します。 これにより、「さまざまな生命領域におけるタンパク質ジスルフィド結合の発生：タンパク質工学、設計および選択、第27巻、第3号に掲載された、タンパク質データバンクのタンパク質の比較」の記事で説明されているタイプの分析が可能になります。 2014年3月1日、65〜72ページ。

PDBには、解像度、方法、変異などが異なる多くの繰り返し構造があります。同じまたは類似のタンパク質で実験を行うと、グループ分析に偏りが生じる可能性があるため、重複セットの中から正しい構造を選択する必要があります。 。 そのためには、非冗長（NR）タンパク質セットを使用する必要があります。

正規化の目的で、3NFを使用するデータベース、またはスタースキーマと次元モデリングを使用するデータベースに原子名をインポートするための化合物辞書をダウンロードすることをお勧めします。 （問題のある構造を排除するためにDSSPも使用しました。この記事ではそれについては説明しませんが、他のDSSP機能は使用しなかったことに注意してください。）

この研究で使用されたデータには、異なる種から得られた少なくとも1つのジスルフィド結合を含む単一ユニットのタンパク質が含まれています。 分析を実行するために、すべてのジスルフィド結合は、最初に、ドメイン（古細菌、原核生物、ウイルス、真核生物、またはその他）によって、また長さによって、連続または非連続として分類されます。

出力

R、SPSS、またはその他のツールに入力できるようにするには、アナリストは次の構造のデータベーステーブルにデータを配置する必要があります。

材料および方法

この目標を達成するための最初のステップは、rcsb.orgからデータを収集することです。 そのサイトには、さまざまな形式の実験のダウンロード可能なPDB構造が含まれています。

データは複数の形式で保存されますが、この例では、フォーマットされた固定スペース区切りのテキスト形式（PDB）のみが使用されます。 PDBテキスト形式の代替はそのXMLバージョンであるPDBMLですが、データ分析で問題を引き起こす可能性のある不正な形式の原子位置命名エントリが含まれている場合があります。 古いmmCIFやその他の形式も利用できる場合がありますが、この記事では説明しません。

PDBフォーマット

PDB形式は、断片化された固定幅のテキスト形式であり、SQLクエリ、Javaプラグイン、Perlモジュールなどで簡単に解析できます。 ファイルコンテナ内の各データ型は、適切なタグで始まる行として表されます。次のサブセクションで各タグタイプについて説明します。 行の長さは80文字以下です。ここで、タグは6文字以下で、1つ以上のスペースを合わせて8バイトを使用します。 通常はCONECTタグの場合、タグとデータの間にスペースがない場合もあります。

TITLE

TITLEタグは、分子名と特定のアミノ酸の挿入、変異、削除などの他の関連データを含む、実験のタイトル（の一部）として行をマークします。

 12345678901234567890123456789012345678901234567890123456789012345678901234567890 TITLE A TWO DISULFIDE DERIVATIVE OF CHARYBDOTOXIN WITH DISULFIDE TITLE 2 13-33 REPLACED BY TWO ALPHA-AMINOBUTYRIC ACIDS, NMR, 30 TITLE 3 STRUCTURES

TITLEレコードに複数の行がある場合は、タイトルを連結する必要があります。連続番号の順に、バイト9と10に右揃えで配置されます（これらの行の数によって異なります）。

ATOM

ATOMラインに保存されているデータは、実験における各原子の座標データです。 実験には、挿入、突然変異、代替の場所、または複数のモデルが含まれる場合があります。 これにより、同じ原子が複数回繰り返されます。 適切な原子の選択については、後で説明します。

 12345678901234567890123456789012345678901234567890123456789012345678901234567890 ATOM 390 N GLY A 26 -1.120 -2.842 4.624 1.00 0.00 N ATOM 391 CA GLY A 26 -0.334 -2.509 3.469 1.00 0.00 C ATOM 392 C GLY A 26 0.682 -1.548 3.972 1.00 0.00 C ATOM 393 O GLY A 26 0.420 -0.786 4.898 1.00 0.00 O ATOM 394 H GLY A 26 -0.832 -2.438 5.489 1.00 0.00 H ATOM 395 HA2 GLY A 26 0.163 -3.399 3.111 1.00 0.00 H ATOM 396 HA3 GLY A 26 -0.955 -2.006 2.739 1.00 0.00 H

上記の例は、実験1BAHから抜粋したものです。 最初の列はレコードのタイプを示し、2番目の列はアトムのシリアル番号です。 構造内のすべての原子には、独自のシリアル番号があります。

シリアル番号の横には、4バイトかかる原子位置ラベルがあります。 その原子位置から、元素の化学記号を抽出することができます。これは、レコードデータの別の列に常に存在するとは限りません。

原子名の後に3文字の剰余コードがあります。 タンパク質の場合、その残基はアミノ酸に対応します。 次に、チェーンは1文字でコード化されます。 鎖とは、ギャップの有無にかかわらず、アミノ酸の単一鎖を意味しますが、リガンドを鎖に割り当てることができる場合もあります。 このケースは、次の列にあるアミノ酸配列の非常に大きなギャップから検出できます。 変異を含めて同じ構造をスキャンできる場合もあります。その場合、挿入コードはシーケンス列の後の追加の列で使用できます。 挿入コードには、影響を受ける残基を区別するのに役立つ文字が含まれています。

次の3つの列は、オングストローム（Å）で測定された各原子の空間座標です。 これらの座標の隣には、占有列があります。これは、原子がその場所にある確率を、通常の0から1のスケールで示します。

最後から2番目の列は温度係数で、Ųで測定された結晶の乱れに関する情報を伝達します。 60Ųより大きい値は無秩序を意味し、30Ųより小さい値は自信を意味します。 実験方法に依存するため、PDBファイルに常に存在するとは限りません。

次の列（記号と料金）は通常欠落しています。 化学記号は、前述のように、原子位置の列から収集できます。 電荷が存在する場合、整数の後に+または-が続く記号として接尾辞が付けられます（例： N1+ ）。

TER

これはチェーンの終わりを示します。 この線がなくても、チェーンがどこで終わっているかを簡単に見分けることができます。 したがって、多くの場合、 TER行が欠落しています。

ENDMDLとMODEL

MODEL線は、構造のモデルが開始する場所を示し、モデルのシリアル番号が含まれています。

そのモデルのすべてのアトミックラインの後、 ENDMDLラインで終わります。

SSBOND

これらの系統には、システインアミノ酸間にジスルフィド結合が含まれています。 ジスルフィド結合は他の種類の残基に存在する可能性がありますが、この記事ではアミノ酸のみが分析されるため、システインのみが含まれます。 次の例は、コード132Lを使用した実験からのものです。

 12345678901234567890123456789012345678901234567890123456789012345678901234567890 SSBOND 1 CYS A 6 CYS A 127 1555 1555 2.01 SSBOND 2 CYS A 30 CYS A 115 1555 1555 2.05 SSBOND 3 CYS A 64 CYS A 80 1555 1555 2.02 SSBOND 4 CYS A 76 CYS A 94 1555 1555 2.02

この例では、ファイルに4つのジスルフィド結合がタグ付けされており、2番目の列にシーケンス番号が示されています。 これらの結合はすべてシステインを使用し（列3および6）、それらはすべて鎖Aに存在します（列4および7）。 各鎖の後に、ペプチド鎖における結合の位置を示す残基配列番号があります。 挿入コードは各残基配列の隣にありますが、この例では、その領域にアミノ酸が挿入されていないため、挿入コードは存在しません。 最後の列の前の2つの列は対称操作用に予約されており、最後の列は硫黄原子間の距離であり、Åで測定されます。

このデータにいくつかのコンテキストを与えるために少し時間を取ってみましょう。

rcsb.org NGLビューアを使用して撮影した以下の写真は、実験132Lの構造を示しています。 特に、それらはリガンドのないタンパク質を示しています。 最初の画像はスティック表現を使用しており、CPKカラーは硫黄とその結合を黄色で示しています。 V字型の硫黄結合はメチオニン結合を表し、Z字型の結合はシステイン間のジスルフィド結合です。

次の図では、バックボーンの視覚化と呼ばれるタンパク質の視覚化の簡略化された方法がアミノ酸で色分けされて示されています。ここで、システインは黄色です。 これは同じタンパク質を表しており、その側鎖は除外され、そのペプチドグループの一部のみが含まれています。この場合はタンパク質のバックボーンです。 N末端、Cアルファ、C末端の3つの原子で構成されています。 この写真はジスルフィド結合を示していませんが、タンパク質の空間配置を示すために簡略化されています。

パイプは、ペプチド結合原子をC-アルファ原子と結合することによって作成されます。 システインの色は、CPK配色法の硫黄の色と同じです。 システンが十分に近づくと、それらの硫黄がジスルフィド結合を形成し、それが構造を強化します。 そうしないと、このタンパク質は結合しすぎて、その構造は高温で不安定になります。

CONECT

これらのレコードは、原子間の接続にタグを付けるために使用されます。 これらのタグがまったく存在しない場合もあれば、すべてのデータが入力される場合もあります。ジスルフィド結合を分析する場合、この部分は役立つ場合がありますが、必須ではありません。 これは、このプロジェクトでは、距離を測定することでタグなしの結合が追加されるため、これはオーバーヘッドになり、チェックする必要があるためです。

SOURCE

これらのレコードには、分子が抽出されたソース生物に関する情報が含まれています。 これらには、分類法で簡単に配置できるサブレコードが含まれており、タイトルレコードで見たのと同じ複数行構造になっています。

 SOURCE MOL_ID: 1; SOURCE 2 ORGANISM_SCIENTIFIC: ANOPHELES GAMBIAE; SOURCE 3 ORGANISM_COMMON: AFRICAN MALARIA MOSQUITO; SOURCE 4 ORGANISM_TAXID: 7165; SOURCE 5 GENE: GST1-6; SOURCE 6 EXPRESSION_SYSTEM: ESCHERICHIA COLI; SOURCE 7 EXPRESSION_SYSTEM_TAXID: 562; SOURCE 8 EXPRESSION_SYSTEM_STRAIN: BL21(DE3)PLYSS; SOURCE 9 EXPRESSION_SYSTEM_VECTOR_TYPE: PLASMID; SOURCE 10 EXPRESSION_SYSTEM_PLASMID: PXAGGST1-6

NRフォーマット

これは、非冗長（NR）チェーンPDBセットのリストです。 そのスナップショットはftp.ncbi.nih.gov/mmdb/nrtable/にあります。 その目的は、タンパク質の類似性によって引き起こされる不必要なバイアスを回避することです。 NRには、すべてのPDB構造の比較によって作成された異なるIDp値レベルの3つのセットがあります。 結果は、後で説明するテキストファイルに追加されます。 このプロジェクトではすべての列が必要なわけではないため、重要な列のみを説明します。

最初の2つの列には、上記のATOMレコードで説明したように、一意のPDB実験コードとチェーン識別子が含まれています。 列6、9、およびCには、BLASTによって計算されたシーケンスの類似性のレベルであるp値の代表性に関する情報が含まれています。 その値がゼロの場合、セットの一部として受け入れられません。 値が1の場合、それはです。 上記の列は、それぞれ10e-7、10e-40、および10e-80のp値カットオフを持つセットの受け入れを表しています。 分析には、p値のカットオフが10e-7のセットのみが使用されます。

最後の列には、構造の受け入れ可能性に関する情報が含まれています。ここでaは受け入れ可能であり、 nは受け入れられません。

 #--------------------------------------------------------------------------------------------------------------------------- # 1 2 3 4 5 6 7 8 9 ABCDEFGHIJKLMNOPQ #--------------------------------------------------------------------------------------------------------------------------- 3F8V A 69715 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 2 0 1 2 0 39139 1 1 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

データベースの構築とデータの解析

何を扱っているのか、何をする必要があるのか​​がわかったので、始めましょう。

データのダウンロード

この分析のすべてのデータは、次の3つのアドレスにあります。

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

最初の2つのリンクには、アーカイブのリストが含まれています。 それぞれの最新のアーカイブを使用して、解像度や品質の欠如に起因する問題を回避する必要があります。 3番目のリンクには、最新の分類法アーカイブが直接含まれています。

データの解析

通常、PDBファイルの解析は、Java、Perl、またはPythonのプラグインまたはモジュールによって行われます。 この調査の場合、事前に作成されたPDB解析モジュールを使用せずにカスタムPerlアプリケーションを作成しました。 その理由は、大量のデータを解析する場合、私の経験では、実験データを使用する際の最も一般的な問題はデータのエラーです。 座標、距離、行の長さ、あるべきではない場所のコメントなどにエラーがある場合があります。

これに対処する最も効果的な方法は、最初にすべてを生のテキストとしてデータベースに保存することです。 一般的なパーサーは、仕様に完全に準拠する理想的なデータを処理するように作成されています。 しかし実際には、データは理想的ではありません。これについては、Perlインポートスクリプトが掲載されているフィルタリングのセクションで説明します。

データベース構築

データベースを構築するとき、このデータベースはデータを処理するために構築されていることに注意してください。 後で分析はSPSSまたはRで行われます。ここでの目的のために、少なくともバージョン8.4のPostgreSQLを使用することをお勧めします。

テーブル構造は、ダウンロードしたファイルから直接コピーされますが、わずかな変更が加えられています。 この場合、レコードの数が少なすぎるため、正規化に時間を費やす価値はありません。 前述のように、このデータベースは使い捨てです。これらのテーブルはWebサイトで提供されるように構築されておらず、データを処理するためだけにあります。 それが終了したら、おそらく他の研究者によるプロセスを繰り返すために、それらをドロップするか、補足データとしてバックアップすることができます。

この場合、最終結果は1つのテーブルになり、SPSSやRなどの統計ツールで使用するためにファイルにエクスポートできます。

テーブル

ATOMレコードからのデータ抽出は、 HEADERまたはTITLEレコードに接続する必要があります。 データ階層は次の図で説明されています。

この図は、第3正規形（3NF）のデータベースを簡略化して表したものであるため、私たちの目的では、オーバーヘッドが多すぎます。 理由：ジスルフィド結合を検出するための原子間の距離を計算するには、結合を行う必要があります。 この場合、テーブルはそれ自体に2回結合され、2次構造と1次構造にもそれぞれ2回結合されます。これは、非常に遅いプロセスです。 すべての分析に二次構造情報が必要なわけではないため、データを再利用したり、大量のジスルフィド結合を分析したりする必要がある場合に備えて、別のスキーマを提案します。

ジスルフィド結合は他の共有結合ほど頻繁ではないため、ウェアハウスモデルは使用できますが、必要ありません。 以下のスタースキーマとディメンションモデリングは、開発に時間がかかりすぎ、クエリがより複雑になります。

すべての結合を処理する必要がある場合は、スタースキーマをお勧めします。

（それ以外の場合、ジスルフィド結合は他の結合ほど一般的ではないため、必要ありません。この作業の場合、ジスルフィド結合の数は30,000に近く、3NFでは十分に速いかもしれませんが、私はそれを経由して処理することにしました正規化されていないテーブルなので、ここには示されていません。）

すべての共有結合の予想される総数は、三次構造の原子数の少なくとも2倍であり、その場合、3NFは非常に遅いため、スタースキーマ形式を使用した非正規化が必要です。 そのスキーマでは、一部のテーブルに2つの外部キーチェックがあります。これは、2つのアトム間に結合が作成されるため、各アトムには独自のprimary_structure_id 、 atom_name_id 、およびresidue_idが必要です。

d_atom_nameディメンションテーブルに入力する方法は2つあります。データから、および別のソースから、前述の化学成分ディクショナリからです。 その形式はPDB形式に似ていますRESIDUE行とCONECT行のみが役立ちます。 これは、 RESIDUEの最初の列に3文字の残基コードが含まれ、 CONECTに原子の名前とその接続が含まれているためです。これらは原子名でもあります。 したがって、このファイルから、すべての原子名を解析してデータベースに含めることができますが、データベースにリストされていない原子名が含まれる可能性を考慮に入れることをお勧めします。

 RESIDUE PRO 17 CONECT N 3 CA CD H CONECT CA 4 NC CB HA CONECT C 3 CA O OXT CONECT O 1 C CONECT CB 4 CA CG HB2 HB3 CONECT CG 4 CB CD HG2 HG3 CONECT CD 4 N CG HD2 HD3 CONECT OXT 2 C HXT CONECT H 1 N CONECT HA 1 CA CONECT HB2 1 CB CONECT HB3 1 CB CONECT HG2 1 CG CONECT HG3 1 CG CONECT HD2 1 CD CONECT HD3 1 CD CONECT HXT 1 OXT END HET PRO 17 HETNAM PRO PROLINE FORMUL PRO C5 H9 N1 O2

このプロジェクトでは、コーディングの速度は、実行の速度やストレージの消費よりも重要です。 正規化しないことにしました。結局のところ、私たちの目標は、イントロで説明した列を含むテーブルを生成することです。

このパートでは、最も重要なテーブルについてのみ説明します。

主なテーブルは次のとおりです。

• proteins ：実験名とコードの表。
• ps ： sequence 、 chain_id 、およびcodeを含む一次構造テーブル。
• ts ：生データから抽出され、 ATOMレコード形式に変換された3次/4次構造を含むテーブル。 これはステージングテーブルとして使用され、抽出後に削除できます。 リガンドは除外されます。
• sources ：実験データが得られた生物のリスト。
• tax_names 、 taxonomy_path 、 taxonomy_paths ：NCBI分類データベースからのリンネ式分類名。 sourcesにリストされている生物から分類パスを取得するために使用されます。
• nr ：NRセットから抽出されたNCBI非冗長タンパク質のリスト。
• pdb_ssbond ：特定のPDBファイル内のジスルフィド結合のリスト。

データのフィルタリングと処理

データは、RCSBPDBリポジトリのスナップショットから取得されます。

各ファイルは、Perlスクリプトを使用してPostgreSQLデータベースの単一のテーブルraw_pdbにインポートされます。 スクリプトは、チャンクごとに10,000挿入のトランザクションを使用します。

raw_pdbの構造は次のとおりです。

インポートスクリプトは次のようになります。

 #!/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 <table> <file>" if not defined($ARGV[0]); $recordCount = 0; $table = $ARGV[0]; $fName = $ARGV[1]; open F, "zcat $fName|"; while (<F>) { 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;

行がインポートされた後、以下で定義する関数を使用して行が解析されます。

raw_pdbデータから、対応するレコードを解析することにより、テーブルts 、 ps 、 proteins 、 sources 、 sources_organela 、およびss_bondを生成します。

psテーブルには、 chain 、 length 、 sequenceの3つの重要な列があります。 タンパク質配列は、各鎖および残基配列順に並べられた各残基のC-アルファ原子を使用して生成され、最初の挿入と最初の代替位置のみを取ります。 chainはTS.chain列から取得され、 lengthは単にsequence文字列の事前計算された長さです。 この記事は単鎖と鎖内接続のみを分析することを目的としているため、ここでの分析では多鎖タンパク質はスキップされます。

SSBONDレコード内では、すべてのジスルフィド結合がpdb_ssbondテーブルに格納されます。このテーブルはpdb_ssbond_extendedテーブルから継承されます。 pdb_ssbondは次のようになります。

これらのインデックスも追加しました。

 "pdb_ssbond_pkey" PRIMARY KEY, btree (id) "ndxcode1" btree (code, chain_id1, res_seq1) "ndxcode2" btree (code, chain_id2, res_seq2)

カットオフ前のジスルフィド結合の分布はガウス分布であると想定されているため（たとえば、KSテストでテストせずに）、同じタンパク質のシステイン間の各距離で標準偏差を計算して、許容される結合長の範囲を見つけて比較しました。カットオフに。 カットオフは、計算された平均プラスマイナス3標準偏差と同じでした。 SSBOND行のPDBファイルに含まれていなかったが、 ATOMレコード間の距離を計算することによって自分自身を挿入した、より多くの可能性のあるジスルフィド結合を導入するために範囲を拡張しました。 ssbondsに選択された範囲は、 ssbondsです。これは、平均（2.05）プラスマイナス4標準偏差です。

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

TSテーブルにはすべての原子の座標が含まれていますが、 " SG "という名前の硫黄を含むシステインのみが使用されます。 そのため、検索するレコードの数を減らしてプロセスを高速化するために、 " SG "硫黄原子のみを含む別のステージングテーブルが作成されます。 硫黄のみを選択した場合、組み合わせの数はすべての原子の場合よりもはるかに少なくなります。前者は194,574で、後者は122,761,100です。 結合されたこのテーブル内で、距離はユークリッド距離を使用して計算され、結果はpdb_ssbondテーブルにインポートされますが、距離が以前に計算された定義済みの長さの間にある場合のみです。 このスピードアップを行う理由は、プロセス全体を再度実行する時間を短縮し、チェックの目的で1日以内に維持するためです。

ジスルフィド結合を洗浄するには、次のアルゴリズムを使用します。

• 接続の両側が同じアミノ酸を指している場合は削除します
• 長さが1.6175344456264と2.48801947642267の間にない結合を削除します
• 挿入を削除します
• 交互の原子位置によって引き起こされた結合を削除しますが、最初の位置は残します

このためのコード（最初の引数としてpdb_ssbondテーブルを使用）は次のとおりです。

 CREATE OR REPLACE FUNCTION pdb_ssbond_clean2( clean_icodes boolean, clean_altloc boolean, mark_reactive boolean, mark_consecutive boolean) RETURNS void AS $$ declare _res integer; BEGIN delete from 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 where chain_id1 = chain_id2 and res_seq1 = res_seq2 and i_code1 = i_code2; update pdb_ssbond set is_consecutive = true , is_reactive = false; delete from pdb_ssbond where not dist between 1.6175344456264 and 2.48801947642267; if clean_icodes then delete from pdb_ssbond where i_code1 not in ('', ' ', 'A') or i_code2 not in ('', ' ', 'A') ; end if; if clean_altloc then delete from 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 < a.ser_num ) ; delete from pdb_ssbond a where exists ( select 1 from pdb_ssbond b where b.code = a.code and b.chain_id2 = a.chain_id2 and b.res_seq2 = a.res_seq2 and b.i_code2 = a.i_code2 and b.ser_num < a.ser_num ) ; end if; if mark_reactive then update pdb_ssbond set is_reactive = false ; update pdb_ssbond set is_reactive = true where exists ( select 1 from pdb_ssbond b where b.code = pdb_ssbond.code and b.chain_id1 = pdb_ssbond.chain_id1 and b.res_seq1 = pdb_ssbond.res_seq1 and b.i_code1 = pdb_ssbond.i_code1 and b.ser_num < pdb_ssbond.ser_num ) ; update pdb_ssbond set is_reactive = true where exists ( select 1 from pdb_ssbond b where b.code = pdb_ssbond.code and b.chain_id2 = pdb_ssbond.chain_id2 and b.res_seq2 = pdb_ssbond.res_seq2 and b.i_code2 = pdb_ssbond.i_code2 and b.ser_num < pdb_ssbond.ser_num ) ; update pdb_ssbond set is_reactive = true where (code, chain_id1, res_seq1, i_code1) in ( select code, chain_id1 as c, res_seq1 as r, i_code1 as i from pdb_ssbond intersect select code, chain_id2 as c, res_seq2 as r, i_code2 as i from pdb_ssbond ) ; update pdb_ssbond set is_reactive = true where (code, chain_id2, res_seq2, i_code2) in ( select code, chain_id1 as c, res_seq1 as r, i_code1 as i from pdb_ssbond intersect select code, chain_id2 as c, res_seq2 as r, i_code2 as i from pdb_ssbond ); end if; if mark_consecutive then update pdb_ssbond set is_consecutive = false; update pdb_ssbond set is_consecutive = true where not exists ( select 1 from pdb_ssbond a where a.code = pdb_ssbond.code and ( (a.chain_id1 = pdb_ssbond.chain_id1 and a.chain_id2 = pdb_ssbond.chain_id2) or (a.chain_id1 = pdb_ssbond.chain_id2 and a.chain_id2 = pdb_ssbond.chain_id1) ) and ( (a.res_seq1 between pdb_ssbond.res_seq1 and pdb_ssbond.res_seq2) or (a.res_seq2 between pdb_ssbond.res_seq1 and pdb_ssbond.res_seq2) or (pdb_ssbond.res_seq1 between a.res_seq1 and a.res_seq2) or (pdb_ssbond.res_seq2 between a.res_seq1 and a.res_seq2) ) and not ( a.code = pdb_ssbond.code and a.chain_id1 = pdb_ssbond.chain_id1 and a.chain_id2 = pdb_ssbond.chain_id2 and a.res_seq1 = pdb_ssbond.res_seq1 and a.res_seq2 = pdb_ssbond.res_seq2 ) ); end if; END; $$ LANGUAGE plpgsql;

これにより、非冗長なタンパク質のセットがnrテーブルにインポートされ、 psおよびproteinsテーブルに結合され、セットにマークが付けられます（ set7 、 set40 、およびset80 ）。 最後に、タンパク質の量に応じて、1セットのみが分析されます。 PDBとNRの間の不一致のチェーンは、分析から削除されます。

ジスルフィド結合のないタンパク質は、どのセットにも属さないタンパク質とともに、研究から除外されます。 データはDSSPで処理され、解像度に問題があるか、処理するには原子が多すぎるこれらのファイルも除外されます。 異なる鎖を持つ接続の数を数えることによってssbondテーブルから簡単に計算されますが、鎖間の接続が保持されなかったため、分析の結果として単一の鎖を持つタンパク質のみが使用されます。

残りのタンパク質については、結合の総数と重複する結合の数が更新され、これは各セットに対して行われます。

ソース生物は、 SOURCEレコードから抽出されます。 それが未知、合成、設計、人工、またはハイブリッドである場合、それは研究から破棄されます。 低解像度のタンパク質も、側鎖が見えない場合にのみ除外されます。

SOURCEレコードは、分類行を含むsourcesテーブルに格納されます。 場合によっては、分類法が欠落しているか、正しくありません。 このような場合、専門家の手動修正が必要です。

NCBIからダウンロードしたソースと分類法から、クラスは各一次構造に割り当てられます。 クラスを割り当てることができない場合、タンパク質は分析リストから削除されます。 生物学データベースが使用されていることを知っているので、生物学者はすべてのソースレコードとNCBI分類クラスの追加チェックを実行することをお勧めします。そうしないと、異なるドメイン間の分類に問題が生じる可能性があります。 ソースセルラーの場所を分類IDと照合するために、ソーステーブルのデータがテーブルsources_organelaに抽出され、すべてのデータがコード、タグ、および値として書き込まれます。 その形式を以下に示します。

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

使用する分類法アーカイブは、7つのダンプファイルを含むZIPファイルです。 これらのファイルの中で、最も重要な2つはnames.dmpとmerged.dmpです。 ドキュメントで詳しく説明されているように、両方のファイルはCSVタブパイプで区切られたファイルです。

• ファイルmerged.dmpには、以前の分類IDのリストと、それぞれがマージされた現在の分類IDが含まれています。
• names.dmpには、次の重要な列がこの順序で含まれています。
  • tax_id ：分類ID
  • name_txt ：種の名前、および該当する場合は一意の名前（種が複数の名前で見つかる場合、これは明確にするのに役立ちます）
• division.dmpには、クラスとして使用するトップレベルドメインの名前が含まれています。
• nodes.dmpは、分類IDを使用した生物の階層構造を含むテーブルです。
  • これには、 names.dmpにある外部キーである親分類IDが含まれています。
  • また、関連するトップドメインデータを正しく保存するために重要な部門IDも含まれています。

このすべてのデータと手動による修正（正しい生命の領域の設定）により、 taxonomy_pathテーブルの構造を作成することができました。 データのサンプリングは次のようになります。

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

分析の前に、バイアスを回避するために、シーケンスの同一性のレベルをチェックする必要があります。 Although the NR set contains sequences which are already compared between each other, an extra check is always recommended.

For each disulfide bond's prior statistical analysis, data is marked if it is reactive or overlapping. After marking overlaps, that information automatically reveals how many consecutive and non-consecutive bonds are inside each protein, and that data is stored in the proteins table from which all protein complexes are excluded in final result.

Each disulfide bond is marked also for its association to sets, by checking both bond sides to see if they belong to the same NR set. Where that is not the case, the disulfide bond is skipped for that set analysis.

To analyze the quantity of bonds by their variance, each length has to be put in a specific class. In this case, only five classes are chosen as written in the function below.

 CREATE OR REPLACE FUNCTION len2class(len integer) 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$ LANGUAGE plpgsql VOLATILE COST 100;

Most of the protein sizes are less than 400 amino acids, so length classification is done by splitting lengths into five ranges, each taking 100 amino acids, except the last one which takes the rest. Below you can see how to use this function to extract data for analysis:

 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) AS expression_system_cellular_location, ( SELECT s.val FROM sources_organela s WHERE s.code = p.code::bpchar AND s.tag::text = 'CELLULAR_LOCATION'::text) AS cellular_location, 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 AND p.chain_cnt = 1 AND p.is_set7 = true AND p.reactive_cnt7 = 0 ORDER BY p.code;

An example result with corrected titles and some manually added columns is shown below.

処理仲介者としてのPostgreSQL

この作業では、データのフェッチから分析までの処理方法を示しました。 科学データを扱う場合、正規化が必要な場合と不要な場合があります。 別の分析に再利用されない少量のデータを処理する場合は、処理が十分に高速な場合は、データを非正規化したままにしておくだけで十分です。

これがすべて1つのバイオインフォマティクスデータベースで行われた理由は、PostgreSQLが多くの言語を統合できるためです。 これには、データベース内で直接統計分析を行うことができるRが含まれます。これは、バイオインフォマティクスツールに関する今後の記事の主題です。

Toptalの同僚であるStefanFuchsとAldoZelenの貴重な相談に特に感謝します。