Pythonを介したPlotlyの使用により、pdbファイルとして得たDNAの立体データをxyz座標のみのデータに変換し、html出力してブラウザ上で表示させます。

• Plotly …オープンソースのデータ分析・グラフ可視化ツール（公式サイト）
• pdbファイル …PDB（Protein Data Bank）が公開している、X線やNMRを用いた結晶構造解析の結果を記したデータファイル（公式サイト）

今回用いる方法では原子半径を考慮せずに全て均一なサイズの点として表現するため、あくまで原子中心間の距離関係がわかるものとしてお取り扱いください。

関連記事

• [PythonとPlotlyで富士山の立体図をブラウザ上で表示する]
• [PythonとPlotly_Expressで気象グラフをhtml出力する]
• [PythonとmatplotlibでDNAの立体図（pdb）を表示する]

index

• コード全体
• データの取得
• 座標ファイルの変換
• グラフパラメータの指定
• 出力

コード全体

動作環境：

• python 3.6.8
• plotly 4.5.0
• windows10 64bit

import plotly.graph_objects as go
import numpy as np

data = open("input.pdb", "r")

with open('pdb.txt', 'w') as f:
    lines = data.readlines()
    line_list = []
    for i in range(len(lines)):
        if 'ATOM' in lines[i][0:4]:
            line = lines[i].replace(lines[i][0:31], '').replace(lines[i][54:80], '')
            line_list.append(line)

    for i in range(len(line_list)):
        f.write(line_list[i])

data.close()

pts = np.loadtxt('pdb.txt')
x, y, z = pts.T

fig = go.Figure(data=[go.Scatter3d(x=x, y=y, z=z,
                mode='markers',
                marker=dict(
                    size=6,
                    color=z,             
                    colorscale='Viridis',
                )
                )
]
)

fig.update_layout(scene = dict(
                    xaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),
                    yaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),
                    zaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),),
                )

with open('plotly_dna.html', 'w') as f:
    f.write(fig.to_html(include_plotlyjs='cdn'))

重要な箇所を以下に示します。

データの取得

今回用いたpdbデータは、PDBより取得しています（2020/02/22現在）。操作内容を簡単に示します。

上記URL先へ移動した後、まずは上部の検索ボックスより「DNA」を検索します。

文字入力後、【Go】をクリックすると検索結果が表示されます。いくつか構造の候補が表示されますが、DNAの結晶構造はDNAのみの場合よりも、それと相互作用している酵素や、インターカレーションしている化合物との結合を見るために結晶構造が解析されている場合があります。そのため、出来るだけDNAだけが構造解析されているファイルを取得したいと思います。

基本的な結晶構造は早期に取得されている場合があるため、【Sort】の項目から年代順【↑Release Date: Oldest to Newest】を選択します。
そうすると、自動的に古い年代から新しい年代へ昇順にソートされます。

再び複数の検索結果が得られます。この中から「B-DNA」が解析されているものを選びます。細かい説明は省略しますが、DNAは大きくA, B, Z型に分かれており、塩基間の距離やらせんの向き等が異なります。その中でもB型DNA（B-DNA）は右巻き二重らせんとして一般的な構造をとったDNAです。

「1BNA」がありましたので、こちらのデータを取得します。構造は1981年に取得されています。

構造名もしくはPDBID（1BNA）をクリックして個別のページに移動すると、エントリーのより詳細な情報が記載されています。左にある構造の画像から3dを選択することで立体構造がブラウザ上で確認できます。今回はpdbファイルを取得したいため、ページ右側の【Download files】から【PDB Format】を選択します。

選択後、ファイルのダウンロードが行われます。なお、ファイル名は1bna.pdbとなっていますが、今回のコードのinput file名に合わせるためにinput.pdbと変更させておきます。
得られたファイルを今回のPythonのスクリプトファイルと同じディレクトリに移動させて、以降の作業を行います。

座標ファイルの変換

取得したpdbファイルの内容は以下のようになっています。

HEADER    DNA                                     26-JAN-81   1BNA              
TITLE     STRUCTURE OF A B-DNA DODECAMER. CONFORMATION AND DYNAMICS             
COMPND    MOL_ID: 1;                                                            
COMPND   2 MOLECULE: DNA (5'-D(*CP*GP*CP*GP*AP*AP*TP*TP*CP*GP*CP*G)-            
COMPND   3 3');                                                                 
COMPND   4 CHAIN: A, B;                                                         
COMPND   5 ENGINEERED: YES                                                      
SOURCE    MOL_ID: 1;                                                            
SOURCE   2 SYNTHETIC: YES                                                       
KEYWDS    B-DNA, DOUBLE HELIX                                                   
EXPDTA    X-RAY DIFFRACTION                                                     
AUTHOR    H.R.DREW,R.M.WING,T.TAKANO,C.BROKA,S.TANAKA,K.ITAKURA,                
    .
    .
    .
ATOM      1  O5'  DC A   1      18.935  34.195  25.617  1.00 64.35           O  
ATOM      2  C5'  DC A   1      19.130  33.921  24.219  1.00 44.69           C  
ATOM      3  C4'  DC A   1      19.961  32.668  24.100  1.00 31.28           C  
ATOM      4  O4'  DC A   1      19.360  31.583  24.852  1.00 37.45           O  
ATOM      5  C3'  DC A   1      20.172  32.122  22.694  1.00 46.72           C  
ATOM      6  O3'  DC A   1      21.350  31.325  22.681  1.00 48.89           O  
ATOM      7  C2'  DC A   1      18.948  31.223  22.647  1.00 30.88           C  
ATOM      8  C1'  DC A   1      19.231  30.482  23.944  1.00 36.58           C  
ATOM      9  N1   DC A   1      18.070  29.661  24.380  1.00 40.51           N  
    .
    .
    .

Header以降には、取得したpdbファイルに関する基本的な情報（日時、タイトル、分子名など）が詳細に記載されています。
また、下部には【ATOM】として各原子の座標と分類が記載されています。
簡単な内容についてはPDBj - PDBレコードの説明（簡易版）を参考にしました。ATOMに関するPDBjによる説明を以下に引用します。

今回の目的は、Plotlyのグラフ上に構造を表示させることであるため、これらの原子座標の情報のみを取得していきます。

data = open("input.pdb", "r")

with open('pdb.txt', 'w') as f:
    lines = data.readlines()
    line_list = []
    for i in range(len(lines)):
        if 'ATOM' in lines[i][0:4]:
            line = lines[i].replace(lines[i][0:31], '').replace(lines[i][54:80], '')
            line_list.append(line)

    for i in range(len(line_list)):
        f.write(line_list[i])

data.close()

lines = data.readlines()として、変数dataとして取得した先ほどのpdbファイルにおける各行をリストとして取得します。
そして、linesの長さの分だけforループを行います。ループの中身としては、【ATOM】が行の先頭にあるとき、すなわちその行が原子座標に関する情報であるときに、その行をあらかじめ空のリストとして宣言しておいたline_listに追加していきます。その際、lines[i].replace(lines[i][0:31], '').replace(lines[i][54:80], '')として、行中の座標以外の情報（元素名など）を消去します。

グラフパラメータの指定

pts = np.loadtxt('pdb.txt')
x, y, z = pts.T

fig = go.Figure(data=[go.Scatter3d(x=x, y=y, z=z,
                mode='markers',
                marker=dict(
                    size=6,
                    color=z,             
                    colorscale='Viridis',
                )
                )
]
)

fig.update_layout(scene = dict(
                    xaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),
                    yaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),
                    zaxis = dict(
                        showbackground=False,
                        showticklabels=False,
                        showline=False,
                        showgrid=False,
                        title_text = ""),),
                )

with open('plotly_dna.html', 'w') as f:
    f.write(fig.to_html(include_plotlyjs='cdn'))

得られたpdb.txtをnumpyで読み込みます。読み込まれたデータは、1列目がx、2列目がy、3列目がzとして読み込まれます。

goとして読み込んだplotly.graph_objectsからFigure()を呼び出します。その中で各種パラメータを指定していきます。
go.Scatter3d()の引数を指定していきます。marker=dict()内でマーカー（点としてプロットされるもの）の設定を行います。size=6でマーカーのサイズを、color=zで色付けの参照値（今回はzの値に合わせて着色）を、colorscale='Viridisとしてカラースケールの種類を指定しています。

また、今回は軸を表示させないシンプルな表示方法を用いたいため、fig.update_layout()で各条件を指定しています。x, y, z等の各axisについて、showbackground=False等を指定することで、軸や軸ラベル、グリッドなどのプロット以外のすべてを非表示にするように設定しています。

最後に、plotly_graph.htmlに出力させています。

出力

これにより得られる出力はhtmlですが、かなり長いため詳細な内容は割愛させていただきます。
得られたhtmlの内容をブラウザ上で表示させると以下のようなデータとしてあらわされます。