2022 Volume 21 Issue 1 Pages 20-32
高等生物のゲノム配列にはタンパク質に翻訳されるエクソン領域とそれ以外のイントロン領域が存在し，その複雑な構造の中から機能部位を推定するのは一つの大きな課題である．一方，タンパク質アミノ酸配列にはその機能と密接な関係がある特徴的な配列パターンであるモチーフが保存され，遺伝子配列にもよく保存されていると考えられる．本研究ではコドン縮約表現とその近似表現を用いた，動的計画 (DP) 法によるゲノム機能部位推定システムを開発した．さらに，公共データベースのモチーフ辞書PROSITEのアミノ酸配列モチーフに対応する遺伝子配列を定量的に表現したコドン重み行列を基にした，モチーフコドン縮約表現を提案した．はじめに，HumanのTNNC1ゲノム配列に対してタンパク質のコーディング配列 (CDS) の推定実験を行った．次に，Humanをはじめとするいくつかのモデル生物種のTNNC1とHPCAゲノム配列に対してEF-handモチーフの推定実験を行い，システムの有用性を示した．
In eukaryotic genome sequences, there are exons that are translated into proteins, and introns that are not. It is important task to estimate the functional sites in the genome sequence. In the other hand, it is well known that the amino acid sequence of a protein is closely related to its function. This is especially true for particular structural features called motifs, and they are considered to be well reserved sites in the genome sequence. In this work, we have developed the Dynamic Programming (DP) -based functional site estimation system using the codon reduced representation and its approximation. The motif codon reduced representation has been also proposed based on codon weight matrix which is represented frequency of each nucleotide corresponding to the amino acid sequence in PROSITE motif. PROSITE is public database of motif dictionary. Our system successfully estimated Coding sequence (CDS) region in TNNC1 genome sequence of Human. The experiments were also executed using EF-hand motif in TNNC1 and HPCA genome sequence of several model species including Human. These results show the potential applicability of our approach for the functional sites in genome sequence.
生命はゲノムによって決められており，地球上に共存する無数の生物はすべてゲノムを持っている．ゲノムには様々な生物を構築し，その生命を維持するのに必要な遺伝情報が保存されている．ゲノムは編集や翻訳を経てタンパク質となり，特定の三次元構造が形成されて生物学的機能が発現する．これらのアミノ酸配列，立体構造，機能などは遺伝子の本体であるDNAの遺伝情報が基となっている．すなわち，ゲノムDNAの核酸塩基配列を解析することにより，生命の生体機能や進化の過程を解明することができる．高等生物のゲノム配列では，RNAからmRNAへ編集される際，残る領域をエクソン，除去される領域をイントロンと呼ぶ．そして，最終的な発現系であるタンパク質をコードする領域をコーディング領域 (CDS領域) と呼ぶ．これは生命の機能に直接関係している部分であるといえる．
ポストゲノム計画の進展および構造決定技術の進歩により，核酸やタンパク質といった生体高分子の配列データは急激に増加し続けている．生体高分子の配列データベースとして，核酸塩基配列データベースのGenbankやEMBL，アミノ酸配列データベースのUniProt/Swiss-protが挙げられる [1,2,3]．さらには，研究者によってチェックされた冗長性のない高品質な配列が登録されたNCBIのRefSeqがある [4]．これらのデータベースに蓄積された大量のデータから，配列-機能相関など，分子生物学上の意味のある情報 (知識) を獲得する事が望まれている．現在，60種類以上の真核生物のゲノム配列を決定するプロジェクトが進行中だが，増加し続ける配列データに対する生物学的解釈が追いついていないのが現状である [5]．そのため，ゲノム配列から遺伝子そしてタンパク質コーディング領域を予測するための正確かつ迅速なツールは数多く開発されている．Kroghらは，原核生物である大腸菌遺伝子のHMM (Hidden Markov models) を提案し，学習したモデルより遺伝子構成未知のゲノム配列上で尤度の高い遺伝子群，すなわちタンパク質コーディング領域の推定を行った [6]．また，Uberbacherらは，真核生物のゲノム配列上の特徴的なパターンの頻度を算出し，ニューラルネットに求められた値を適用することによってゲノム配列からタンパク質コーディング領域の推定を行った [7]．現在，Claverie，Guigó, Haussler，BurgeとKarlinなどにより，多くの遺伝子予測プログラムが公開されている [8,9,10,11]．しかし，特に複雑な構造をもつ高等生物では，完全にコーディング領域を推定することは難しく大きな課題となっている [12]．
タンパク質アミノ酸配列には，その機能と密接な関係がある特徴的な配列パターンが存在し，これをモチーフと呼ぶ．そしてそれは遺伝子配列にも保存されていると考えられる．モチーフは共通の配列をアライメントすることにより見出すことができる． モチーフを見出すための配列アライメントは動的計画 (Dynamic Programming, DP) 法を用いて求めることができる．アライメントを高速に求めることができるプログラムにBLASTがあり，マルチプルアライメントを行うためのプログラムとしてClustalWが知られている [13, 14]．また，文献等からアミノ酸配列モチーフの情報を集積したデータベースとしてPROSITEが広く知られている [15]．
Ohtomoらは7つのモデル生物種を対象として，PROSITEのアミノ酸配列モチーフとそれに対応する遺伝子配列とタンパク質立体構造情報を集積した，タンパク質-遺伝子モチーフ辞書システムを提案した [16]．このシステムによって，アミノ酸配列レベルでは捉えきれなかったモチーフの配列特徴を塩基配列レベルで見出すことができるようになった．Katoらはコドンとアミノ酸を1対1の関係に表したコドン縮約表現を提案し，タンパク質-遺伝子モチーフを利用したゲノム配列機能部位推定システムを開発した [17]．しかしTroponin C type 1 (TNNC1) ゲノム配列に対してEF-handモチーフ推定実験を行った結果，TNNC1タンパク質配列にモチーフが保存されているにも関わらず推定できなかった．これについて，ゲノム配列のアノテーションをもとに配列を調査し，途中にイントロンが挿入されたEF-handモチーフが保存されていることを確認した．先行研究における推定システムは，ゲノム配列とコドン縮約表現とのパターンマッチングによって実現されている．そのためイントロンが挿入されるパターンには対応できない．そのようなモチーフの発見には，タンパク質や遺伝子配列のアノテーションを手がかりに手動で見つける必要があった．イントロンを考慮した機能部位の推定は，DP法による配列アライメントによって実現できると考えられる．しかし，コドン縮約表現は全てが1対1に対応づいておらず，そのままDP法に適応することは困難である．イントロンを考慮したタンパク質コーディング領域やタンパク質配列モチーフ部位，ゲノム配列の機能部位を推定できれば，未知のゲノム配列の機能やその保存位置について知ることができる．さらにイントロンが挿入された機能部位を生物種間で比較することで，その進化についての知見が得られると考えられる．
本研究ではコドン縮約表現とその近似表現を導入し，それを用いたDP法による配列アライメントシステムを開発した．さらに，アミノ酸配列モチーフに対応する遺伝子配列の塩基使用率を定量的に表現したコドン重み行列をコドン縮約表現とした，モチーフコドン縮約表現を提案した．これをクエリとしたゲノム配列の機能部位推定の結果から，ゲノムの機能部位解析への応用を試みる．
DNAの塩基配列はアデニン (a)，グアニン (g)，シトシン (c)， チミン (t) の4種類の塩基で構成されており，アミノ酸を指定する塩基配列の3塩基1組 (トリプレット=3連子) の単位をコドンと呼ぶ．コドンの組み合わせは43で64通りあり，1つのアミノ酸に対して複数のコドンが対応している． 4種類の塩基とは別に，IUPACではワイルドカードが定義されている．例えば，プリン塩基であるアデニンまたはグアニンはp (プリン)，ピリミジン塩基であるシトシンまたはチミンはy (ピリミジン)，そして任意の塩基はn (any) と表現される．Table 1にワイルドカード表現の一覧を示す．この表現方法を遺伝暗号表に適用したものをコドン縮約表現と呼ぶ [17]．本文中では原則として，ゲノム配列の塩基やワイルドカードを英小文字，タンパク質のアミノ酸を英大文字の1文字表記で表現して区別する．
コドン縮約表現では，2種類のコドンで表現されるフェニルアラニン (F) は"tty"の1パターン，6種類のコドンで表現されるアルギニン (R) は"cgn"，"agr"の2パターンで表すことができる．コドン縮約表現によって20種類中17種類のアミノ酸とコドンを1対1の関係で表すことができる．一方で，1対1の関係で表すことのできないアミノ酸はロイシン (L)，セリン (S)，アルギニン (R) の3種類存在する．これらのアミノ酸を一意に指定できるよう，コドン縮約表現の近似を行った．例えば，"cgn"，"agr"の2種類のコドン縮約表現で表されるアルギニンは，第1塩基目のcとaをm，第3塩基目のnとrをnとして近似することで，1パターンの"mgn"と表現できる．同様に，ロイシンとセリンもコドン縮約表現の近似を行うことで，それぞれ"ytn"，"wsn"と表現できる．コドン縮約表現とその近似を適用した遺伝暗号表をFigure 1に示す．コドン縮約表現の近似は，例えばアルギニンの"mgn"はセリンの"agr"のように余分なアミノ酸にも対応するため，使用には注意が必要である．なお，本研究では使用しないが，終止コドン"taa"，"tag"，"tga"は近似することにより，"trh"で表すことができる．
Genetic code table of codon reduced representation and its approximation
コドン縮約表現の近似によって，アミノ酸配列の逆翻訳は1つの遺伝子配列として表現することができる．例えば，カルシウム結合に関係があるEF-handモチーフに注目する [18]．HumanのTNNC1タンパク質に保存されているEF-handモチーフ配列 "DEDGSGTVDFDEF" は，5残基目のセリンを逆翻訳すると"tcn"と"agy"となるため2通りの遺伝子配列になる．近似表現を用いると"wsn"となり，遺伝子配列を一意に定めることができる．アミノ酸配列から遺伝子配列へ逆翻訳する過程をFigure 2に示す．
Reverse translation of "DEDGSGTVDFDEF"
PROSITEに登録されているアミノ酸配列モチーフの正規表現をそのままコドン縮約表現を用いて表すと，ほとんどのコドンが自由度の高すぎるコドンとなり，検索対象とは無関係な配列を多数検出してしまうと考えられる．そこで，重み行列とワイルドカード表現を用いて遺伝子レベルのモチーフパターンに対応付けを行った．これをモチーフコドン縮約表現と定義する．これにより，遺伝子レベルでよく保存されているモチーフパターンをクエリ配列として表現することができ，ゲノム配列中に保存されているモチーフ部位を効率よく推定できると期待される．
重み行列とは，タンパク質‐遺伝子モチーフによるアミノ酸配列モチーフ部位に対応した塩基配列を基にコドン単位の塩基の出現頻度を計算したものである [16]．この重み行列からは，アミノ酸配列モチーフの正規表現では見えない塩基の偏りが数値として表現され，定量的に評価することができる． アスパラギン酸は"gat"， "gac"の2種類のコドンが存在するが，"gat"が8回と"gac"が2回出現した場合の重み行列をFigure 3に示す．
Codon weight matrix for the aspartic acid pattern
対応付けの方法として，まず重み行列の出現頻度に対して閾値を設ける．次にモチーフパターンにおける各コドンの重み行列と閾値を比較し，閾値以上ならモチーフの対応する塩基とし，以下なら非対応の塩基として扱う．そして，対応する塩基をワイルドカード表現に当てはめる．例えば，"[DENSTG]"の重み行列がFigure 4左であるとすると，閾値0.25ではFigure 4右のように対応コドン"ray"と表現することができる．
Codon reduced representation based on weight matrix with threshold
配列アライメントとは，異なる配列の類似性を求める手段であり，配列中で同じ並び方をしている文字列を探すことで2つないし3つ以上の配列群を比較する手法である．配列間でより多くの文字が一致，アライメントスコアが最大となるよう，ギャップを挿入しながらアライメントを行う．部分的に類似している領域を発見するために，DP法を用いたローカルアライメントがよく使用される．しかしゲノム配列中の機能部位は1つのエクソン領域に保存されているとは限らないため，イントロン領域とマッチしないよう連続でギャップを挿入する必要がある．連続したギャップの挿入が可能なローカルアライメントとしてSmith-Waterman法が知られている [19]．Gotohは，ギャップの開始と連続の場合で異なるペナルティを与えるアフィンギャップを提案し，Smith-Waterman法に適応した [20]．この方法は広く使われており，ギャップが連続する場合のペナルティを極めて小さい値にすることで長大なギャップの挿入にも対応できると考えた．本研究ではSmith-Waterman法とアフィンギャップを用いたゲノム配列機能部位推定のためのアライメント手法を提案した．具体的には，コドン縮約表現による遺伝子配列を用いたアライメントが可能となるよう，塩基のワイルドカードに対応した置換行列を定義した．
初めにSmith-Waterman法とアフィンギャップによる計算式を以下に示す．
ここで，Miは比較配列Mのi番目の塩基，Qjはクエリ配列Qのj番目の塩基である．
上記の式と置換行列を基にしてスコア行列を作成することで，アライメントのスコア計算を行う．スコア行列の作成をする際に，最大値を取るスコアの計算に用いたセルの方向の記録を行う．スコア行列の計算後，最大スコアのセルから記録した方向に向かってトレースバックすることで，アライメント結果を得ることができる．このとき，アフィンギャップがある領域にトレースが進んだ場合，アフィンギャップが終了するまでギャップを挿入するようなトレースを行う．
例として，比較配列M = "gataaaatc"とクエリ配列Q = "gatatc"によるアライメントを行う．置換行列sをTable 2で定義したものとし，開始ギャップペナルティα=10，伸張ギャップペナルティβ=1とする．スコアの最大値からトレースバックを行う場合，アライメントのトレースする方向が分岐するような場合がある．分岐数が少ない場合はそのすべての分岐をトレースすれば最適なアライメントを得ることができるが，分岐数が膨大な場合，すべての分岐をトレースすることは非効率である．そこで，比較配列側に優先したトレースバックとクエリ配列側に優先したトレースバックを行うことにより，アライメントで表現される範囲の最短と最長のパターンを検出する．2つの配列のアライメントのスコア行列作成ならびにトレースバックした例をFigure 5に示す．比較配列Mとクエリ配列Qのアライメントは2パターン導出され，トレースする方向の優先度によって，同じスコアでも異なったアライメント結果を得る場合がある．
An example of calculating sequence alignment
DP法によるローカルアライメントではスコアが最大値からトレースバックを行っていくことで最適なアライメントを得ることができる．しかし，最大値からのトレースバックを開始するだけでは1つの配列中に複数の機能部位が存在した場合，機能部位を取りこぼす可能性がある．そこで，スコアの最大値から始まるアライメントに機能部位が保存されていると仮定すると，最大値に近いスコアから始まるアライメントも機能部位を保存していると考え，トレースバック開始のスコアを以下のように定義する．トレース開始スコア = 最大スコア*検索範囲(%)
このトレース開始スコア以上であればトレースを行い，そのアライメント結果を機能部位候補として導出する．これにより，1配列中に複数存在する機能部位の導出を期待する．
アフィンギャップによるアライメントとコドン縮約表現を適用した置換行列がゲノム配列の機能部位推定に有用であるかを確認するために，はじめにCDS領域の推定実験を行った．具体的には，アミノ酸質配列全体をコドン縮約表現で逆翻訳した遺伝子配列をクエリとし，ゲノム配列に対してアライメントによるCDS領域を推定した．実験にはRefSeqより，複数のエクソンにまたがってCDS領域が存在するHumanのTNNC1ゲノム配列 (NG_008963.1) とアミノ酸配列データ (NP_003271.1) を使用した．入力する各配列データはFasta形式とした．Fasta形式は先頭に '>' から始まるヘッダー行があり，それ以降に配列情報が記述される．ヘッダーにはNCBIで定義された遺伝子ID (gi) ，配列ID (ref)やタンパク質名が含まれている．ここで，推定実験で使用したデータセット中のゲノム配列の各塩基は大文字で表記されていることに注意が必要である．これはRefSeqのオリジナルの配列情報が大文字表記のためである．逆翻訳した配列，クエリ配列は英小文字で表記している．Figure 6にアミノ酸配列データを示し，Figure 7に縮約表現で逆翻訳した遺伝子配列データを示す．逆翻訳した遺伝子配列のヘッダー行には，逆翻訳したことがわかるよう独自の記述を加えている．
An amino acid sequence of NP_003271.1. The '>' line is header information including IDs and protein name information defined by NCBI (excerpt). Sequence information is described after the header information.
A reverse translated nucleotide sequence of NP_003271.1.
1つの遺伝子に含まれるイントロンの数は1個から40個以上までと大きく異なる．イントロンの長さも50塩基程度のものから，10,000塩基以上のものまである．本実験では，開始ギャップペナルティα=20に対して伸張ギャップペナルティβ= 0，検索範囲を80%とした．これにより，長大なイントロン領域を持つゲノム配列に対しても適切なアライメントが期待される．
Figure 8に推定した出力結果を示す．出力結果の'M:'の行はデータセット配列のアライメント結果を示し，'Q:'はクエリ配列のアライメント結果を示している．その配列間の行には, 対応した文字のマッチング記号を示している． '|' は完全一致， '+' はワイルドカードとの一致， ':' はanyとの一致を示している．空欄の場合は不一致である．最後にアライメントスコアやマッチングした位置情報を示している．推定結果は，アライメントの領域が複数のセグメントに分かれていることを推定した．Figure 8中の'segment'が推定した各セグメントの情報であり，(5027 … 5049) は5,027塩基目から5,049塩基目までが1つのセグメントを示す．推定部位が正しいかを確認するため，ゲノム配列のCDSアノテーション情報と比較した．アノテーション情報をFigure 9に示す．アノテーション情報は先頭にIDや登録日時，定義名等の情報が記載される．それ以降に様々なアノテーション情報が記載される．'CDS'アノテーションには，そのゲノム配列上でのCDS領域の位置やその遺伝子名，翻訳後のタンパク質配列情報等が記述されている．推定結果は，各セグメントの開始位置や終了位置に1塩基のずれが生じている場合もあるが，おおよそCDS領域と一致している．アノテーションの最終セグメントは推定結果より末尾3塩基分多いが，これは終始コドンである．以上の結果より，アフィンギャップによるアライメントとコドン縮約表現を適用した置換行列がゲノム配列の機能部位推定に有用であるといえる．
Estimated result of CDS region in TNNC1 genome sequence (excerpt). 'M:' lines show alignment result of dataset sequence. 'Q:' lines show alignment result of query sequence. Matching result between dataset and query sequences is shown matching symbols, '|': exact match, '+': wildcard match, and ':': any match. Correct CDS region is colored red in dataset sequence. Lines after 'segment:' show estimated CDS region. In this result, it is shown that CDS regions have preserved six regions separately.
Annotation of NG_008963.1 (excerpt). The 'CDS' annotation shows gene name and position of CDS region. The 'join' annotation shows position of the CDS region. The '/gene' annotation shows gene name. In the sequence, CDS region of TNNC1 gene is preserved six segments between 5027 and 7767.
先行研究において，HumanのTNNC1ゲノム配列中には3つのEF-handモチーフ中にイントロン領域が挿入されていることが示された [17]．ゲノム配列に対してモチーフコドン縮約表現をクエリとしたCDS領域推定を行えば，イントロン領域にまたがる形で保存されているモチーフ部位の発見が期待される．そこで，モチーフコドン縮約表現を用いてゲノム配列に対する機能部位推定を行った．モチーフコドン縮約表現は，タンパク質-遺伝子モチーフ辞書システムに登録された最もデータ数が充実している，HumanのEF-handモチーフ796部位を対象としたコドン重み行列を基に作成した [16]．PROSITEにおけるEF-handモチーフの正規表現パターンをFigure 10に示す．PROSITEの正規表現はいくつかのルールがある．'[ ]'のパターンはその中のアミノ酸だけを許容し，'{ }'のパターンはその中のアミノ酸以外を許容する．'x'は任意のアミノ酸である．パターンを複数回繰り返したい場合は直後の'()'中にその回数を記述する．重み行列ならびにモチーフコドン縮約表現をFigure 11とFigure 12に示す．重み行列はパターンごとに記述している．モチーフコドン縮約表現の閾値は0.15とした．推定条件は開始ギャップペナルティをα=20，伸張ギャップペナルティをβ=0.001，追跡範囲を80%とした．以後の推定実験は，特に断りがない限り，すべて同様の値を使用した．伸張ギャップペナルティを極めて小さい値としたのは，機能部位と関係のない推定候補であるノイズを削減するためである．
A Regular expression of EF-hand motif in PROSITE (ID: PS00018). A '[ ]' pattern allows for an amino acid in a bracket. A '{ }' pattern does not allow an amino acid in a bracket. An 'x' pattern allows any amino acid. An '()' pattern repeats a previous pattern in the number of bracket times.
A codon weight matrix of EF-hand motif. A codon weight shows each pattern of the regular expression of EF-hand motif.
A motif codon reduced representation of EF-hand motif.
TNNC1のゲノム配列に対して機能部位を推定した結果，2ヶ所のモチーフ対応部位を検出した．検出した2ヶ所のモチーフ部位の推定結果をFigure 13に示す．データセット配列中の赤色と青色の文字はアノテーションのCDS領域を示している．また，ゲノム配列 (NG_008963.1) のアノテーション情報をFigure 14に示す．Figure 14中の'CDS'の赤文字と青文字の領域は，Figure 13中のデータセット配列中の領域とそれぞれ対応している．推定結果とRefSeqのアノテーション情報の比較を行った結果，それぞれ推定できたモチーフはイントロンを挟んだ2つのCDS領域内に保存されていた．
Estimated result of EF-hand motif in the TNNC1 genome sequence (excerpt). Correct CDS region is colored red and blue in the dataset sequences.
Annotation of NG_008963.1 (excerpt). The '/translation' annotation shows a protein sequence. The motif① information is colored red, and the motif② information is colored blue. A green string in the protein sequence shows the EF-hand motif that failed to be estimated.
1つ目のモチーフ対応部位は (6927 ... 6940), (7188 ... 7212) の2つのセグメントに分かれて保存されていることを推定した．推定結果とアノテーションを比較すると，4塩基分イントロン領域の先頭にマッチしてしまい，正確にモチーフ領域とイントロンを分けた推定結果は得られなかった．しかし推定位置に注目すると，イントロンを含んだ場合でもモチーフの開始位置や終了位置は正確に推定できた．2つ目のモチーフ対応部位は (7618 … 7651), (7736 … 7740) の2つのセグメントに分かれて保存されていることを推定した．この推定結果はイントロン領域を含めアノテーション情報と一致しており，正確にモチーフ部位を推定できた．一方もう1つのEF-handモチーフが (7294 … 7298), (7515 … 7548) の2つのセグメントに保存されているが，今回のパラメーターでは推定できなかった．Figure 14中のアミノ酸配列に，該当するモチーフ部分を緑文字で示す．これは対象EF-handモチーフの構成残基が希少なパターンであったためと考えられる．該当モチーフ部位の4残基はアラニンとなっており，そのコドンの2塩基目はシトシンである．Figure 11の4残基目の重み表現に注目すると，2塩基目シトシンの出現頻度は0.067であり縮約表現の閾値以下である．そのため，アライメントスコアが低くなり推定できなかったと考えられる．
次に，イントロンが挿入されているモチーフと挿入されないモチーフが混在して保存されているエントリであり，同じカルシウム結合タンパク質の一種であるhippocalcin (HPCA) に注目した．推定した3ヶ所のモチーフ対応部位をFigure 15に，ゲノム配列 (NG_042176.1) のアノテーション情報をFigure 16に示す．
Estimated result of EF-hand motif in HPCA genome sequence (excerpt). The correct CDS region of a motif③ is colored blue in the dataset sequence.
Annotation of NG_042176.1 (excerpt). The motif① and motif② information is colored red, and the motif③ information is colored blue.
機能推定の結果から，モチーフ対応部位が1つのエクソン領域内に保存されている部位と，2つのエクソン領域にまたがる形で保存されている部位をそれぞれ推定した．アノテーション情報と比較した結果，それぞれの推定部位はモチーフ部位に対応していた．
以上の結果より，モチーフコドン縮約表現は遺伝子レベルで保存されているモチーフをイントロンが挿入された場合でもゲノム配列から正しく推定できた．一方で出現頻度が低いコドンを含むモチーフは正しく推定できなかった．モチーフコドン縮約表現の閾値を低くすることで推定できる可能性があるが，ワイルドカード表現が多くなるため機能と関係のない部位を多数誤って推定してしまうことが懸念される．
3.3 他の生物種を対象としたEF-hand モチーフのコドン縮約表現による機能部位推定モチーフは生物種が異なる場合でもそのアミノ酸配列はよく保存され，そしてその遺伝子配列もよく保存されていると考えられる．ならば，他の生物種由来のモチーフコドン縮約表現においても機能部位の推定は可能であると考えられる．そこでHumanのEF-handモチーフコドン縮約表現を用いて，いくつかのモデル生物種での推定を試みた．RefSeqよりMouse, Rat, Cow, Pig, Frog, Zebrafishのゲノム配列 (NC_000080.7，NC_051351.1，NC_037349.1, NC_010455.5, NC_030680.2, NC_007134.7) から，TNNC1のアノテーションが付加された配列部分のみを抽出し推定実験を行った．Pigのゲノム配列のみ相補鎖で登録されているため，事前に元の配列になるよう変換した．
Mouse, Rat, Cow, PigについてはHumanと同様，2つのEF-handモチーフがイントロン領域にまたがる形で保存されていることを推定した．1つ目のモチーフはHumanと同様に正しくモチーフ領域を推定できた．また，イントロンの先頭4塩基とマッチしてしまう結果も同様に見られた．2つ目のモチーフについては，Mouse, Rat, Pigの推定結果にHumanとの共通の差異が見られた．Mouse，Rat，Pigはモチーフの33塩基目がチミンからシトシンに変化していた．シトシンに変化しているためワイルドカード表現のdがマッチせず，推定結果とアノテーションに差異が生じたと考えられる．代表例として，Mouseにおける2つ目のモチーフの推定結果をFigure 17に示す．比較のためHumanの推定結果を同時に示し，変化している部分を図中の緑枠で示す．しかしアノテーションより，モチーフの開始位置と終了位置は正確に推定できており，モチーフ部位は正しく推定できたといえる．Figure 18にMouseのアノテーション情報を示す．Zebrafishは2つ目のモチーフのみ推定できたが，Frogについては1つも推定できなかった．各配列にモチーフが保存されているにも関わらず推定できなかった原因の1つとして，長大なイントロンに対応しきれなかったことが挙げられる．そこでFrogとZebrafishに対して，伸張ギャップを0にして再度推定実験を行った．伸張ギャップペナルティ以外のパラメーターは同様である．結果として，機能部位と関係ない箇所を多く機能部位候補を出力してしまうが，2つ目のモチーフを両方で推定できた．Figure 19に伸張ギャップを0として推定が成功したFrogの例を示す．Humanにおいて推定できなかったモチーフについては，生物種共通で推定できなかった．いずれの生物種もモチーフの4残基目がアラニンとなっているため，Humanと同様に縮約表現の閾値が原因の1つと考えられる．
Comparison of estimated result of the EF-hand motif② in TNNC1 genome sequence between Mouse and Human. In Mouse, unlike in case of Human, 33rd base is changing from thymine to cytosine (see green boxes in the Figure).
Annotation of TNNC1 genome sequence for Mouse (excerpt). The regions are colored blue shows that are preserved the Mouse TNNC1 EF-hand motif②.
Successful estimated result by changing to the extension gap is zero in the TNNC1 genome sequence for Frog (excerpt). Correct CDS region is colored blue in the dataset sequence.
HPCAについても同様に，各モデル生物種のゲノム配列 (NC_000070.7, NC_051340.1, NC_037329.1, NC_010448.4, NC_030678.2, NC_007130.7) から該当部分のみを抽出した配列に対して推定実験を行った．なおPig以外は相補鎖として登録されているため事前に変換した．
Mouse, RatについてはHuman同様に3つのモチーフを推定することができた．Cow, Pig, Frog, Zebrafishについては1つ目のモチーフは正しく推定できた．しかしCow, Pigは2つ目を推定することができず，Frogは3つ目を推定することができなかった．Zebrafishは2つ目と3つ目両方を推定することができなかった．Mouse, Rat, Cow, Pigの3つ目のモチーフは，Humanと同様15塩基目から16塩基目にイントロンが挿入されることを推定した．Mouse, Rat, Cow, Pigにおける3つ目のモチーフ，イントロンが挿入されたモチーフの推定結果をFigure 20に示す．それぞれのアノテーション情報から，モチーフ部位に対応していることが確認できた．FrogとZebrafishにおいて，3つ目のモチーフに長大なイントロンが挿入されていることがアノテーションより確認できた．伸張ギャップペナルティを0にして再度実験を行ったが，どちらも正しく推定することはできなかった．推定結果を確認すると，イントロン前のモチーフ部分にはマッチしているが，イントロン後は高いスコアでマッチする箇所まで延々とギャップを挿入する結果が複数得られた．伸張ギャップペナルティを0にすることは長大なイントロンに対応できる一方，延々とギャップ挿入を可能とするデメリットを持つ．対象配列によっては，このような結果が得られてしまうことに注意が必要である．
Successful estimated results of the motif③ in HPCA genome sequence for Mouse, Rat, Cow and Pig. Correct CDS region is colored blue in the dataset sequences of each species.
機能部位推定実験の結果，EF-handモチーフはタンパク質やその保存位置によって，挿入されているイントロンの長さは異なることを確認した．また，それらは数十塩基から数千塩基と生物種によって大きく異なる一方，その開始位置は生物種共通であった．TNNC1について，1つ目のモチーフは10塩基目から11塩基目の間，2つ目のモチーフは34塩基目から35塩基目の間にイントロンが挿入される．HPCAについて，3つ目のモチーフは16塩基目から17塩基目の間にイントロンが挿入される．脊椎動物のトロポニンCの遺伝子は5つのイントロンを持ち，第1と第4イントロン以外は生物種が異なる場合でも同様の位置に挿入される [21]．さらに一部のカルシウム結合タンパク質のゲノム配列において，ヘリックスループヘリックスのループ部分，PROSITEのEF-handモチーフ部分にイントロンが挿入されることが示されている [22, 23]．TNNC1とHPCAゲノム配列では生物種によらずEF-handモチーフが保存され，それらのイントロンの挿入位置は固定であることが推定結果とアノテーションから得られた．これらの推定結果は先行研究からも生物学的に妥当な結果であると推測できる．イントロンが挿入されたEF-handがTNNC1ゲノム配列中に保存されていることは先行研究において確認されていた．本研究では7種のモデル生物種のTNNC1とHPCAゲノム配列について，イントロンが挿入されたEF-handモチーフの保存位置に関する情報を得られた．さらに，それらのイントロンの挿入位置はモデル生物種共通であった．これらのイントロンが挿入されたEF-handモチーフに関する情報は，既存のデータベースには掲載されていない，本研究で得られた新たな知見である．
本研究ではゲノム配列の機能部位推定のために，コドン縮約表現とその近似表現を用いたDP法による配列アライメントシステムを開発した．さらに，アミノ酸配列モチーフに対応する遺伝子配列の塩基使用率を定量的に表現したコドン重み行列を基にコドン縮約表現を適用したモチーフコドン縮約表現を提案した．ゲノム配列からモチーフ部位を推定するために2つの点を工夫した．1つは，アライメントの開始位置に幅を持たせることで複数のモチーフパターンが存在する配列に対応した．また，クエリ側へのギャップを抑制することでイントロン領域を考慮したアライメントを行えるようにした．タンパク質アミノ酸配列をコドン縮約表現で逆翻訳した遺伝子配列をクエリとして，TNNC1ゲノム配列のCDS領域推定実験を行った．その結果，複数のイントロンが挿入されている場合でも正しくCDS領域を推定できた．次に，EF-handのモチーフコドン縮約表現によるクエリ配列を用いて，Humanをはじめとした複数のモデル生物種のTNNC1とHPCAゲノム配列に対して機能部位推定を行った．その結果，複数のエクソン領域にまたがるモチーフ部位を推定できた．イントロンを考慮したゲノム配列の機能部位推定は，従来システムでは困難であったが本システムによって容易に推定可能となった．これにより，ゲノム配列の解析への足掛かりとして，モチーフコドン縮約表現によるゲノムの機能部位推定の有用性を示した．
しかしまだ完全に機能部位を推定することは難しい．ギャップペナルティや追跡範囲によって様々な推定結果が得られるため，ユーザーは正しい機能部位であるかを都度検証する必要がある．また，モチーフコドン縮約表現は重み行列と閾値によって表現が異なるため，閾値によって機能部位を取りこぼす可能性もある．今後は，推定結果から正しい機能部位であるかを自動的に導き出すことや適切なモチーフコドン縮約表現の閾値設定を決めることが望まれる．また，本研究で使用した推定システムを以下のページで公開している．
https://sunflower.kuicr.kyoto-u.ac.jp/~ohtomo/GenomeSeqFunctionalSiteEstimationSystem.html
