はじめに
一般公開されている、コロナウイルス（SARS-CoV-2）の生物情報の配列データ（biological sequence）の解読を試してみたので体験を共有します。
三回に渡って、生き物の遺伝情報である核酸 (DNA/RNA)の塩基配列やタンパク質を構成するアミノ酸配列の中身を探ってみます。
- Nextstrainでコロナウイルスのゲノム解析を行い解析結果をみてみる
- NCBIの生物情報配列データでコロナウイルスの突起蛋白のアミノ酸配列を解読してみる
- 同じく、コロナウイルスのRNAの塩基配列を解読してみる
三回目は、米国の国立衛生研究所が中心となって運用している生物情報データベース（NCBI）からコロナウイルス（SARS-CoV-2）のゲノムデータをダウンロードして解読してみようと思います。
前回の記事で探ってみた、感染力が増したと疑われる変異株のアミノ酸の変異「N501Y」に対応する、コロナウイルスのRNA塩基配列の変異部分を探ってみたいと思います。
お断り: 筆者はプログラミングが専門で、医療生物は専門ではありません。できるかぎりウラをとったつもりですが、本稿の記載に知識や理解不足から誤りや誤解を生む箇所があるかもしれません。予めご了承ください。ご指摘をいただければさいわいです。
塩基配列の取得
前回解読したコロナウイルスの突起蛋白（S Protein）のアミノ酸配列はNCBIから取得したものです。
これは元になるコロナウイルスの全ゲノム情報から抽出されたものです。
抽出された蛋白データの元になったRNAデータもNCBIに収められています。
NCBIサイトからこれを取得してみます。
NCBI SARS-CoV-2 Resources - NIH
https://www.ncbi.nlm.nih.gov/sars-cov-2/
蛋白データのリソースデータのページを見ると、元になったゲノム配列の登録番号（accession）がわかります。
また、そのゲノム情報のリソースデータを見ると、蛋白データがゲノム中のどこから切り出されたものか（開始終了位置）がわかります。
まず武漢のウイルスの情報を探します。
[基準とする配列を探すリンク]
[武漢のウイルスの塩基配列]
[武漢のウイルスのアミノ酸配列]
続いて豪州で採取されたN501Yを含む変異株の情報を探します。
[登録番号から探すリンク]
[豪州のウイルスのアミノ酸配列]
[豪州のウイルスの塩基配列]
[豪州のウイルスのアミノ酸配列の位置]
以上の情報をまとめると次のようになります。
| Virus | Sequence Type | Accession | File | S Protein Range |
|---|---|---|---|---|---|
| RefSeq | Nucleotide | NC_045512.2 | nuc-refseq.fasta | 21563..25384 |
| RefSeq | Surface Glycoprotein | YP_009724390.1 | pro-refseq.fasta | |
| Au N501Y | Nucleotide | MT745734.1 | nuc-au-200617.fasta | 21525..25346 |
| Au N501Y | Surface Glycoprotein | QLG76817.1 | pro-au-200617.fasta | |
塩基配列のファイル
ファイルの中身をみてみます。
NCBIの核酸の塩基配列はタンパク質のアミノ酸配列と同様、FASTA形式で収められています。
先頭の'>'ではじまるヘッダ行に続き四つの核酸塩基（ATCG）が並んでいます。
RNAではDNAのThymine(T)がUracil(U)に置き換わりますが、NCBIのデータではUではなくTが使われていることがわかります。
$ cat nuc-refseq.fasta
>NC_045512.2 |Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, complete genome
ATTAAAGGTTTATACCTTCCCAGGTAACAAACCAACCAACTTTCGATCTCTTGTAGATCT
GTTCTCTAAACGAACTTTAAAATCTGTGTGGCTGTCACTCGGCTGCATGCTTAGTGCACT
...
TTTAGTAGTGCTATCCCCATGTGATTTTAATAGCTTCTTAGGAGAATGACAAAAAAAAAA
AAAAAAAAAAAAAAAAAAAAAAA
$$ cat nuc-au-200617.fasta
>MT745734.1 |Severe acute respiratory syndrome coronavirus 2 isolate SARS-CoV-2/human/AUS/VIC2150/2020, complete genome
ACTTTCGATCTCTTGTAGATCTGTTCTCTAAACGAACTTTAAAATCTGTGTGGCTGTCAC
TCGGCTGCATGCTTAGTGCACTCACGCAGTATAATTAATAACTAATTACTGTCGTTGACA
...
CACGCGGAGTACGATCGAGTGTACAGTGAACAATGCTAGGGAGAGCTGCCTATATGGAAG
AGCCCTAATGTGTAAAATTAATTTTAGTAGTGCTATCCCCATGTGATTTTAAT
サイズを比べてみます。
$ ls -l nuc-*.fasta
-rw-rw-r-- 1 author author 30430 Feb 7 21:44 nuc-au-200617.fasta
-rw-rw-r-- 1 author author 30500 Jan 18 14:19 nuc-refseq.fasta
$ tail -n +2 nuc-refseq.fasta | tr -d '\n' | wc -c
29903
$ tail -n +2 nuc-au-200617.fasta | tr -d '\n' | wc -c
29813
ファイルサイズはコメントがあるのでもちろんですが、中身のデータ部分のサイズも違います。
豪州の昨年六月の変異株のほうが90個分の塩基が減ってよりコンパクトになっています。
上記のコマンドラインの意味は、tailコマンドで先頭行のコメントを除外、その結果をtrコマンドで改行を除外、wcで文字数を把握。
アミノ酸コード部分の解読
先の表に示したNCBIのリソース情報を見ると、コロナウイルスのRNA全体の塩基配列のうち、第二回で取り上げた突起蛋白（S Protein）を作り出す箇所は、武漢の基準配列の場合は21563～25384の位置に、豪州の変異株の場合は21525～25346の位置にあることがわかります。
タンパク質はアミノ酸でできています。
ゲノムの遺伝情報にはタンパク質のアミノ酸配列を表す部分とそうでない部分があります。
一個のアミノ酸は核酸塩基三つの並びで表されます。
たとえば、アミノ酸のひとつAsparagine(N)はAAUかAAC、Tyrosine(Y)はUAUかUACで表されます。
NCBIのデータに即していえば、UがTに変わります。
| Amino Acid | Amino Code | Nucleobases |
|---|---|---|
| Asparagine | N | AA[TC] |
| Tyrosine | Y | TA[TC] |
それではいよいよ、塩基配列の該当部分（501番目の塩基三組み）が期待通り変化しているか調べてみます。
まず武漢の塩基配列から突起蛋白の501番目のアミノ酸に対応する塩基配列を切り出します。
$ FROM=$(expr 21563 + 500 \* 3)
$ TO=$(expr $FROM + 2)
$ echo $FROM - $TO
23063 - 23065
$ tail -n +2 nuc-refseq.fasta | tr -d '\n' | cut -c $FROM-$TO
AAT
続いて豪州の変異株のほうも調べます。
突起蛋白のアミノ酸配列の開始位置が武漢のものとは異なります。
$ FROM=$(expr 21525 + 500 \* 3)
$ TO=$(expr $FROM + 2)
$ echo $FROM - $TO
23025 - 23027
$ tail -n +2 nuc-au-200617.fasta | tr -d '\n' | cut -c $FROM-$TO
TAT
AATがTATに変化したことで、アミノ酸がAsparagine(N)からTyrosine(Y)に変化したことが確認できました。
tailでコメントを除外、trで改行を除外したあと、電卓のexprでシェル変数のFROMとTOに開始位置と終了位置を設定、cutでその範囲の塩基配列を表示。
まとめ
三回に渡って、コロナウイルス（SARS-CoV-2）のタンパク質のアミノ酸配列や核酸 (RNA)の塩基配列の配列データの解析や解読を行ってみました。
誰でも利用できる一般公開の解析ツールや生物情報データベースを利用することでウイルスの変異を身近に知ることができました。
核酸（RNA/DNA）の塩基配列やタンパク質のアミノ酸配列などの生物情報の交換形式はテキストファイルなので、テキスト処理に優れたUNIX（Linux）のコマンドラインでいとも簡単に扱えることがわかりました。
解読のほか加工や変換もLinuxのテキストコマンド、PythonやPerlなどで作る簡単なスクリプトで実現できそうです。
下記リンクに示したとおり、横浜工文社では、塩基配列やアミノ酸配列の解読や変換に特化したユーティリティを開発しGitHubで公開しています。
参考リンク
横浜工文社の関連ページ
Written 2021/02/08