総説
エピゲノムワイド関連研究（EWAS）による形質関連DNAメチル化サイトの探索
2021 年 2 巻 1 号 p. 58-75
詳細
2021 年 2 巻 1 号 p. 58-75
近年の測定技術の発展は、ゲノム情報だけでなく様々なエピゲノム情報の取得を可能にした。中でもDNAメチル化は、他のエピゲノム情報の担い手と比べて安定的に遺伝子発現制御を行う仕組みである。ヒトのDNAメチル化は、疾患の有無や、これまでの生活習慣によって変動する。このことから、DNAメチル化は種々の疾患の病態解明のカギになると考えられ、さらに診断・発症リスク評価等のバイオマーカーとしての活用も期待されている。DNAメチル化アレイによって、ヒトゲノムのDNAメチル化プロファイルを比較的安価に取得できるようになり、大規模なエピゲノムデータの解析が可能になった。本稿では、エピゲノムワイド関連研究（Epigenome-wide association study, EWAS）に焦点を当て、DNAメチル化のデータ形式や、頻用される解析方法及び解析上の注意事項について概説する。また、近年提案されたアプローチや残されている課題についても概説する。
測定技術の発展によるDNA配列決定の低コスト化は、大規模ヒトゲノムデータの蓄積と分析を可能にした。結果、ヒトゲノムのDNA配列の生物学的機能・意義の解明に貢献し、遺伝学分野に大きな進展をもたらしている。しかし同時に、ヒトゲノムデータだけでは全てを明らかにできないことも判明した。例えばヒトゲノムの場合、約30億の塩基配列に遺伝情報が蓄積されている。この塩基配列は、同一人物由来ならば脳、肺、心臓等のどの細胞からDNAを採取しても全く同じ配列になる（体細胞変異による違いはここでは考慮しない）。各臓器の細胞では、30億塩基対の遺伝情報が全て機能しているのだろうか？そのようなことはなく、細胞の役割や状況に応じて機能する遺伝情報は異なっている。このような仕組みが生命にあることは古くから考えられており、それを研究する学問はエピジェネティクス（Epigenetics）と呼ばれる。エピジェネティクスという言葉は1942年にC. H. Waddingtonによって提唱された［1］。後成説を意味するEpigenesisと、遺伝学を意味するGeneticsを合体した複合語である。エピジェネティクスの定義は、研究の発展に伴い時代に応じて変わりつつある。Waddingtonは「発生過程で遺伝子型が表現型を決めるようになる仕組み」［2］と呼んだが、近年は「DNAの塩基配列の変化をともなわずに、染色体における変化によって生じる、安定的に受け継がれうる表現型」という定義が定着している［3］。
エピゲノムは、ゲノムDNAとゲノムDNAが巻き付いたヒストンに付与されたエピジェネティック修飾の総体である［4］。エピゲノムの中でDNAメチル化（DNA methylation）は、エピジェネティクスの分子基盤の担い手である。以下にDNAメチル化とは何かを概説する。図1に示す通り、DNAは4種類の塩基（アデニン（A）、グアニン（G）、シトシン（C）、チミン（T））のいずれかを持つヌクレオチドが多数つながってできた分子で、塩基配列によって遺伝情報を記録している。DNAメチル化は、このDNA中の塩基配列のシトシン（C）の5位という場所にメチル基（-CH3、図中ではMeと記載）が付与することをいう。DNAメチル化は、ヒトを始めとする真核生物、原核生物の両方に存在する（ただしシトシン以外の塩基にメチル化が起きている例もある）。DNAメチル化は、シトシンならどれでも起きるわけではない。ヒトの場合、図1に示す通りシトシン（C）―グアニン（G）という二塩基配列のシトシン（C）に起きる。それ以外のシトシンには通常DNAメチル化は起きない。この二塩基配列の出現する場所は、シトシン（C）とグアニン（G）の間にリン酸（p）があることを踏まえてCpGサイトと呼ばれる。すなわちヒトの場合、DNAメチル化はCpGサイトのシトシン（C）で起きると言える。CpGサイトはDNAメチル化サイトとも呼ばれる。DNAメチル化サイトのシトシンにメチル基が付与されることをDNAメチル化（又は単純にメチル化、DNA methylation）と呼び、メチル基が外されることをDNA脱メチル化（又は単純に脱メチル化、DNA demethylation）と呼ぶ。DNAは2本鎖構造を取っているので、片方の鎖にCpGサイトがあれば、相補配列であるもう片方の鎖にもCpGサイトが存在する。図1のように、通常、片方の鎖がDNAメチル化されている場合、もう片方の鎖もDNAメチル化されている（フルメチル化と呼ぶ）。CpGサイトはヒトゲノムの様々な場所に存在するが、すべてのCpGサイトがDNAメチル化されているとは限らない。図2のように、サイト1はメチル化されていないが、サイト2やサイト3はメチル化されているということがあり得る。哺乳類のゲノムDNAでは、CpGサイトのおよそ60%～90%がメチル化されている［5］。どのメチル化サイトがメチル化されているかは、細胞の種類や状態によって異なる。ゲノム中のメチル化パターンによって、DNAの働きを制御しているのである。DNAメチル化は、エピゲノム情報の中では比較的安定しており、DNA複製後もDNAメチル化のパターンは維持される。DNA複製後、片方の鎖しかDNAメチル化されていない状態（ヘミメチル化）になるが、その後メチル化パターンがもう片方の鎖にもコピーされる。このDNAメチル化パターンを維持する仕組みについては書籍2～4を参考にしてほしい［2, 3, 4］。
A, C, G, T: DNAの塩基、Me: メチル基（-CH3）
DNAメチル化がどのようにしてDNAの働きを制御しているか、その一例を図3に示した。通常、遺伝子の発現はDNA中の遺伝子の上流にある調節領域（プロモーター）に転写因子が結合して、遺伝子の転写を開始する。しかし、プロモーターにあるDNAメチル化サイトがメチル化されていると、転写因子の結合が阻害される。その結果、遺伝子の発現が抑制される。このため、プロモーター上にあるDNAメチル化は遺伝子発現をオフにするスイッチであると認識されている。また、プロモーター以外のDNAメチル化は異なる働きを持つ。DNAメチル化のより詳細な役割は、書籍［4］を参考にしていただきたい。
ヒトゲノムの塩基配列は一生変わらないが、DNAメチル化はそうではない。DNAメチル化はエピゲノム修飾の中でも比較的安定しており、長期間変化されずに維持される可能性があるものの、様々な要因によってメチル化や脱メチル化が起きる。例えば、加齢や喫煙等の生活習慣によって一部のDNAメチル化サイトの状態が変わることがわかっている。DNAメチル化が疾患と関係し時に治療薬のターゲットになることもわかっている。特にガンの発症と深いかかわりがあり、DNAメチル化を阻害する薬剤がガンの治療薬として開発されている［6, 7］。もし疾患と関係のあるDNAメチル化サイトを見つけられれば、疾患発症機構の解明にもつながると期待できる。また、DNAメチル化を疾患の診断や発症を予測するバイオマーカーとして利用する試みもある（詳細は9章で述べる）。このような理由から、ヒトを対象に疾患や生活習慣と関連するDNAメチル化サイトの探索が進められている。
疾患や生活習慣と関連があるDNAメチル化サイトをヒトゲノム全体から探索するアプローチとして、エピゲノムワイド関連研究（epigenome-wide association study, EWAS）が行われている。この方法は、一塩基多型（single nucleotide polymorphism, SNP）を対象としたゲノムワイド関連研究（genome-wide association study, GWAS）と同様、関連を評価したい疾患や生活習慣を形質とし、各DNAメチル化サイトと形質の関連を評価し、有意に関連するDNAメチル化サイトをゲノム全体から探索する方法である。本来、エピゲノム情報はDNAメチル化以外にヒストン修飾なども含むため、DNAメチロームワイド関連研究と呼ぶべきところかもしれないが、EWASの呼び名が一般的である。図4は、2020年までの各年に出版されたEWAS論文の数である（後述するEWAS Atlasの情報に基づく）。2009年に初めてEWAS論文が報告されてからEWAS論文数は急速に増加し、その後も報告され続け、2020年末の時点で累計878報が出版されてきた。様々な形質に対してEWASが行われており、これまでに肥満度（body mass index, BMI）［8］や2型糖尿病［9］等の疾患、喫煙などの生活習慣［10］、血中のタンパク質濃度［11］等のEWASが行われてきた。日本人を対象としたEWASの論文も発表されている［12, 13, 14, 15, 16, 17］。尚、後述するEWASのメタアナリシスも本稿ではEWASとして扱う。以降、DNAメチル化アレイデータを用いたEWASの流れを概説し、各工程で用いられる統計解析やバイオインフォマティクス手法を概説する。
（EWAS Atlasのデータを基に作成）
ゲノム中のDNAメチル化状態を測定するため、様々な実験方法が開発されてきた。各種測定方法の概要は書籍［4］や総説［18］で説明されているので、本稿では簡単に紹介する。ヒトゲノム中のDNAメチル化情報を取得する方法は大きく2種類に分かれる。一塩基レベルのDNAメチル化状態を観測する方法と、ゲノム中の領域単位のDNAメチル化状態を観測する方法である。EWASでは、形質と関連する特定のDNAメチル化サイトを同定することが目的なため、前者のアプローチで測定したデータを用いる。一塩基単位のDNAメチル化状態をヒトゲノム全体に亘って取得する測定方法として、バイサルファイト変換を利用したショートリードシークエンサーやマイクロアレイによる測定が用いられている。前者はゲノム全体を測定することから、全ゲノムバイサルファイトシーケンシング（whole genome bisulfite sequencing, WGBS）と呼ぶ。WGBSを使えばヒトゲノム中の全てのDNAメチル化状態を一度に測定できる（ただし、ショートリードシークエンサーで塩基配列決定できる領域に限る）。次にバイサルファイト変換を利用したマイクロアレイは、基盤上にプローブと呼ばれるDNA配列をハイブリダイズしたチップを使って測定する方法である。DNAメチル化を測定するためのマイクロアレイをDNAメチル化アレイまたはメチル化アレイと呼ぶ。EWASではIllumina社製のDNAメチル化アレイが主流である。DNAメチル化アレイは、測定するDNAメチル化サイトが決まっており、表1に示す通りこれまでにヒトを対象とした3種類のDNAメチル化アレイがIlluminaから販売されている。これらは測定可能なメチル化サイト数が異なっており、最新のDNAメチル化アレイMethylationEPICならばヒトゲノム中の85万カ所のDNAメチル化サイトのメチル化状態を測定できる。2021年8月現在、MethylationEPIC以外のアレイは販売終了している。しかしながら8章で説明する通り、過去の研究でHumanMethylation27やHumanMethylatoin450を用いて取得された大量のデータが公開されているため、これらの製品も把握しておくと良い。また、近年Illumina社からマウスを対象としたDNAメチル化アレイ（Infinium® Mouse Methylation）の販売が開始されている。これによってマウスのDNAメチル化アレイデータもヒトと同様に取得できるようになった。
| 製品名 | 測定箇所 |
| Infinium HumanMethylation27 BeadChip | 2万7千 |
| Infinium HumanMethylation450 BeadChip | 48万5千 |
| Infinium MethylationEPIC BeadChip | 86万5千 |
下側の製品程新しい。MethylationEPIC以外は販売終了。
WGBSとDNAメチル化アレイを比較した場合、WGBSの方がヒトゲノム全体のDNAメチル化状態を網羅的に測定できるというメリットがあるがその分高コストである。逆に、DNAメチル化アレイの方は測定できるDNAメチル化サイトに限りはあるが、低コストで測定でき、大量のサンプルを必要とする研究に向いている。多数のサンプルを必要とするEWASでは、DNAメチル化アレイのデータを用いていることがほとんどである。本総説はIllumina社のDNAメチル化アレイを用いる前提で話を進める。
この章では、DNAメチル化データがどのようなデータ型で表現されるのか、生物学的・実験的背景を踏まえて説明する。DNAメチル化データは図5のようなp行×n列のマトリックスで表すことが多い。pは解析対象となるDNAメチル化サイトの数であり、nは解析対象となるサンプルサイズである。この図を見るとわかる通り、通常、DNAメチル化の状態は割合を示す連続値（0～1）で表現する。「シトシンにメチル化が有るか無いかの2値のデータになるのでは？」と不思議に思う読者もいるかもしれない。DNAメチル化データが割合なのは、一つの細胞のDNAメチル化状態を観測しているのではなく、複数の細胞のDNAのメチル化状態を観測していることに起因する。通常、細胞からDNAを抽出する場合、一つの細胞からDNAを抽出するのではなく、幾つもの細胞の混合物からDNAを抽出している。この時、細胞毎にDNAのメチル化の状態は異なっている。図6に示す通り、複数細胞から得たDNAのメチル化程度を観測しているため、観測した全細胞の内、いくつの細胞でDNAメチル化が起きていたかという割合として定量化される。図6のメチル化サイト1の場合、測定した全細胞がDNAメチル化されていないためメチル化割合=0となり、サイト2のように全体の半数でDNAメチル化されていれば0.5 （=50%）、サイト3のように全ての細胞でDNAメチル化されていれば1.0 （=100%）というデータが1サンプルから得られる。DNAメチル化アレイで測定したデータはβ値として数値化される。これが、DNAメチル化割合に相当する。すなわち、DNAメチル化アレイデータはβ値のマトリックスなのである。
DNAメチル化アレイではDNAメチル化の程度を示す指標としてβ値に加えM値という指標が用いられることがある。β値とM値は相互に変換可能で、次式で計算できる［19］。
ここでβは、あるサンプルの特定のDNAメチル化サイトのβ値を意味し、MはM値を意味する。M値はβの底を2とした対数オッズ比であり、両者は底を2としたロジット変換で相互変換していると言える。EWASでは、各DNAメチル化サイトのβ値やM値と、疾患の有無や生活習慣のような形質との関連を統計解析によって評価する。β値は割合データであるため、［0, 1］の範囲しか値を取り得ない。そして低値や高値の値を取るDNAメチル化サイトの場合、β値は不等分散性（heteroscedasticity）を示すため、多くの統計解析手法で仮定されている等分散性の仮定から逸脱することになる。一方、M値は取り得る値の範囲に制限はなく、近似的に等分散性（homoscedasticity）を示すことから、β値をM値に変換することでより適切に統計解析に利用でき得る［19］。β値はDNAメチル化割合そのものを示しており直感的な生物学的解釈が可能である一方、M値は生物学的解釈が困難である。そのため、β値をそのまま統計解析に用いている事例も多い。β値とM値両者の利点を生かす方法として、統計解析の際にはM値を用いて評価し、メチル化割合の平均値や中央値、効果量を表記する際には、β値に戻して評価することもある［19, 20］。
DNAメチル化アレイの生の測定データは、idat形式のファイルとして保存されているが、これは図5のようなマトリックスデータになっていない。idatファイルから図5のようなデータを得るには前処理が必要である。idatファイルは、DNAメチル化サイトごとにDNAメチル化の程度をシグナル値として記録しており、バックグラウンド補正（background correction）や正規化（normalization）の過程を経てβ値やM値を算出し、マトリックスデータに加工する。idatファイルからβ値（M値）のマトリックスデータへの変換は専用のソフトウェアを用いることで実施可能である。DNAメチル化アレイを開発したIllumina社が、変換用のソフトウェアGenomeStudio 2011.1 （https://jp.support.illumina.com/array/array_software/genomestudio/downloads.html）を公開している。最新のGenomeStudioのバージョンはv2だが、こちらはSNPデータの処理用であり、メチル化データの処理には以前のバージョンである2011.1を用いる必要がある。表1に示した通りIllumina社のDNAメチル化アレイは3種存在するが、GenomeStudio 2011.1はどのタイプのデータにも適用可能である。また、idatファイルからの変換用RパッケージがBioconductorで開発されている。代表的なものとして、minfi ［21］、methylumi ［22］、wateRmelon［23］等が挙げられる。これらのパッケージはマニュアルも整備されており利用しやすい。これらのRパッケージはHumanMethylation27やHumanMethylation450のデータ処理用に開発されたものである。その後に開発されたMehylationEPICに対応しているのは、minfi［24］、wateRmelonである。他のパッケージも将来対応する可能性がある。idatファイルからβ値（M値）のマトリックスに変換する過程で、多くのバイオインフォマティクス的手法が開発されている。本稿では詳細を取り扱わないが、各処理の内容が知りたい方は総説［25］やGenomeStudioのマニュアル（https://jp.support.illumina.com/downloads/genomestudio_m_module_v18_ug_(11319130_b).html）をご覧頂きたい。β値のマトリックスデータを取得した後、さらにクオリティコントロールが必要である。取得したデータの中で、低品質と判定されたサンプルやDNAメチル化サイトを統計解析前に除いたり［25］、クロスハイブリダイズの影響［26］、測定プローブ中のSNPの影響［27］も加味して、解析対象とするDNAメチル化サイトを絞り込む必要がある。
DNAメチル化アレイデータはDNAサンプルから測定するため、度々SNPなどの塩基配列データと同じ感覚で捉えてしまうことがあるが、1章で述べた通りDNAメチル化はDNAサンプルを取得する臓器や細胞種、状態によって異なる。本章ではこの点を踏まえて改めてDNAメチル化データの注意事項を説明する。
ヒトは、原則としてどの組織・どの細胞から取得したDNAであっても基本的に同じ塩基配列情報を有する（体細胞変異の影響を無視した場合）。しかしながらDNAメチル化情報は、データを例えば血液（白血球）、肺の細胞、脳の細胞のいずれから取得したのかによって大きく異なる。多くの場合、EWASで用いるDNAメチル化データは、コホート研究等で採取された血液由来のDNAサンプルから測定することが多い。そのため、白血球のDNAメチル化プロファイルでEWASを行っている。そのようなデータを用いて形質と関連するDNAメチル化サイトを探索する場合、なぜ白血球のDNAメチル化状態とその形質が関連するのかを考えて研究を進めなくてはならない。もちろん白血球以外の組織を対象とする場合も有り得るので、何由来のDNAから測定したデータか解析前に都度確認する必要がある。白血球のDNAメチル化データでEWASを行う場合、白血球が深くかかわる自己免疫疾患［28, 29］を形質とすれば、発症機構の解明に結び付けやすく、結果の解釈・考察が行いやすいだろう。白血球と直接関係しない形質でもEWASは実施されている。例えば、アルツハイマーと関連するDNAメチル化サイトを探索したいのなら、疾患とかかわりの深い脳組織のDNAメチル化データでEWASを行うのが理想であろう。実際に脳組織での実施例もある［30］。しかしながら、脳組織の採取は容易ではないため、比較的容易に行える白血球のDNAメチル化データで実施した例もある［31］。この場合、白血球でEWASをなぜ行ったのか説明し、有意な関連を持つDNAメチル化サイトが見つかったならば、なぜ白血球で有意な関連が確認できたのか、検討する必要がある。また疾患の診断や発症リスクを評価するためのバイオマーカーとしてEWASでDNAメチル化サイトを探索するなら、比較的採取し易い白血球を用いることはバイオマーカーの実用化という面では妥当といえる。
SNP等のDNAの塩基配列情報は、生まれてから死ぬまで原則一生変わることは無い。そのためSNPデータと形質の関連を評価する際には、SNPが形質に影響を及ぼしたという因果関係を考えやすい利点があった。一方、DNAメチル化状態は加齢や生活習慣の影響を受けて変化するため、DNAメチル化データにはこの利点が無いことに注意すべきである。前向き研究か、後ろ向き研究か等の疫学的視点で研究デザインを考えてデータを取得し、解析結果を解釈する必要がある。
EWASでは、DNAメチル化サイトごとにβ値やM値といったメチル化状態と形質との関連を評価する。ここでは関連解析の内容とDNAメチル化アレイデータを扱う上での注意事項を概説する。
６-１．回帰モデルによる関連解析DNAメチル化と形質の関連解析には様々な方法があるが、患者群と健常者群のような2群間でDNAメチル化サイトのメチル化状態に差があるか評価する解析がよく行われる。このような群間で比較する解析を、メチル化差解析（differential methylation analysis, DMA）と呼ぶことがある。DMAの最も単純な方法として、β値やM値をt検定やWilcoxonの順位和検定で比較することが考えられる。しかしながら多くの場合、これらの検定方法の使用は適切ではない。DNAメチル化データは生物学的・実験的な交絡因子によって変動することがわかっており、これらの因子の影響を除いて関連を評価することが求められているからである。したがって、交絡因子を考慮した解析手法の適用が望まれる。現在のEWASでは、下記のように種々の交絡因子を加味した線形回帰モデルを考える。
ここで、yにDNAメチル化の程度（β値やM値）、xにDNAメチル化との関連を評価したい形質を配置する。仮に患者群と健常者群の2群比較を行いたい場合、患者群をx=1、健常者群をx=0と数値化する。xは連続値であっても良い。例えば、肥満の程度と関連するDNAメチル化サイトを探索したい場合、xに体格指数（body mass index, BMI）の値を入れることもある。covariatek, celll, batchmは各種共変量であり、詳細は6-2で後述する。w0が切片で、wx, wk, wl, wmが回帰係数である。εは誤差項であり、正規分布に従うと仮定する。EWASではDNAメチル化の程度yと形質xの関連を評価したいので、wxに注目する。推定値
Kupers et al. によって行われた新生児の出生時体重を形質としたEWASメタアナリシスの結果。図は論文で公開されている要約統計量を基に作成。青い点は有意水準に達するP値を示したDNAメチル化サイト。横線：有意水準=1.06×10-7
EWASの場合、後述する様々な交絡因子の影響を受けたDNAメチル化状態が観測されていると考え、（2）式のようにDNAメチル化状態を従属変数y、形質を独立変数xとして採用することが多い。しかし、前向き研究などにおいて、DNAメチル化状態が原因で形質に影響すると考える場合は、DNAメチル化状態を独立変数x、形質を従属変数yとして採用することもある。その場合は、（2）式のような線形回帰モデルだけでなく、形質のデータの型によって一般化線形モデルの枠組みでロジスティック回帰やCoxの比例ハザードモデルを使い分ける。EWASの結果を見る際には、DNAメチル化情報が独立変数なのか従属変数なのか注意して読み取る必要がある。また、通常のEWASはサンプル間に近縁関係が無いことを仮定している。もし、近縁者（related individuals）を含めて関連解析を行う場合、混合モデルを採用することもある。
EWASの有意水準は幾つに設定すると良いだろうか？EWASでは多数のDNAメチル化サイトとの関連を検定する。HumanMethylation450やMethylationEPICの場合、数十万回の検定を行うことになる。5%という有意水準では多数の偽陽性が生じるため、多重検定の補正が必要と考えられる。多くのEWASでは、ボンフェローニ補正［33］かBenjamini-Hochberg法等によるFDR（false discovery rate）で制御する。ボンフェローニ補正は保守的すぎるため、有意水準を幾つにするのが適切か？という検討も進められている。Saffariらはシミュレーションデータの検討結果から、ファミリーワイズエラー率=5%を維持するには、HumanMethylation450の有意水準として2.4×10−7、ゲノム全体を対象とした探索の場合は、有意水準3.6×10−8が適切であると提案している［34］。ManshellらはMethylationEPICにおける有意水準として9×10−8を提案している［35］。
６-２．DNAメチル化アレイデータにおける交絡因子6-1で述べた通り、DNAメチル化アレイデータには様々な交絡因子の影響が入り込んでいるため、これらの影響を考慮して関連解析を行う必要がある。大きく分けて3種類の因子を考える必要がある。一つ目は（1）疫学的交絡因子、もう一つは（2）細胞組成（cell-type composition）、最後が（3）バッチ効果（batch effect）である。以下で順に紹介する。
（1）疫学的交絡因子まず疫学的交絡因子とは、ヒト集団を対象とした時に考慮すべき交絡因子のことである。一般的に性別、年齢、喫煙習慣、地域等がよく取り上げられる。性差［36］、年齢差［37］のあるDNAメチル化サイトの存在はEWAS初期から報告されている。EWASを実施する際には、年齢や性別を共変量として加えるべきであろう。喫煙習慣によって変動するDNAメチル化サイトも多数知られているため［10］、共変量として加えることが望ましい。その他、評価したい形質に応じて、適切な調整因子を加える必要がある。DNAメチル化アレイデータを測定したサンプルに対してSNPアレイデータも取得している場合、集団階層化（population stratification）を補正するためにSNPデータに主成分分析（principal component analysis, PCA）を適用して得られた主成分を共変量として加えることもある［38］。
（2）細胞組成続いて、細胞組成を概説する。EWASに用いるDNAメチル化データは、サンプルの取得し易さや既存試料としてバイオバンキングされている点から、血液由来のデータであることが多い。血液由来のDNAとは白血球由来のDNAを意味する。しかしながら、「このメチル化データは白血球のDNAメチル化プロファイルだ」と捉えるだけでは十分ではない。白血球には、様々な種類が存在する。具体的には顆粒球、単球、CD4+T細胞、CD8+T細胞、B細胞、NK細胞がある（実際には更に細分化できるが省略する）。これら6種の細胞間でDNAメチル化のプロファイルは異なっている［39, 40］。また、サンプル1では単球が全白血球中の5%を占めるが、サンプル2では単球は10%を占めるというように、白血球6種の組成はサンプルによって異なる。つまり、EWASの対象となる白血球から取得したDNAメチル化データは、これらの6種類の白血球のDNAメチル化プロファイルを組成割合で重み付き平均した値である。そのため、この影響を考慮しないで関連解析を行うと、細胞の組成割合が交絡因子となって誤った結果を招く恐れがある。そこで、DNAメチル化アレイデータからサンプル毎の細胞組成を推定し、得られた組成割合を関連解析の共変量として追加する必要がある。オミクスデータから細胞組成を推定する過程を細胞型デコンボリューション（cell-type deconvolution）と呼び、これまでに様々な方法が提案されている［41］。EWASでは、特にJaffeらによって実装されたmifiパッケージのestimateCellCountによる細胞組成の推定値が頻用されている［42］。Jaffeらは、他の研究で取得された6種の白血球それぞれのDNAメチル化アレイデータ［40］をリファレンスとしてHousemanらが開発した方法［39］を用いて細胞組成を計算する方法を構築した。これにより、白血球のメチル化アレイデータを入力情報として、6種の白血球の細胞組成をサンプル毎に推定できる。EWASではこのようにして得た細胞組成を（2）式の共変量に投入して関連解析を行う。
DNAメチル化の測定対象となる細胞は白血球に限らない。脂肪細胞や各種臓器の細胞等様々な細胞が測定対象となになり、これらの細胞も様々な細胞種の混合である。そのため、細胞組成を考慮した解析が必要である。Jaffeらの方法はあくまでも白血球のための方法であり、それ以外の細胞種にそのまま適用することはできない。これらの手法の開発も現在進められている。詳細は総説［41］をご覧頂きたい。
（3）バッチ効果最後にバッチ効果とは、「測定条件間で質的に異なる行為によって生じ、生物学的または科学的変数とは無関係な測定のサブグループのこと」である［43］。言葉だけでは伝わりづらいため、例を挙げて説明する。例えば、全16サンプルに対してDNAメチル化アレイを測定したいとする。この時、“8サンプルを月曜に測定し残りを火曜日に測定した場合”、“8サンプルずつロットの試薬、アレイ、装置等を変えて測定した場合”、“8サンプルずつを別々の実験者が測定した場合”等にバッチ効果が生じ得る。このような何らかの測定条件の違いによって、測定値に偏り（バイアス）をもたらすのがバッチ効果である。DNAメチル化アレイでは、図5のマトリックスデータを生成する際に、正規化（normalization）を適用するが、多くの方法ではバッチ効果までは取り除けない。解決策としては①実験デザインによる対処、②データ取得後の対処の2通りが考えられている。
まず、①実験デザインによる対処は、比較したいサンプル間で測定条件に偏りが出ないように測定するということである。②データ取得後の対処を考えるよりも先にこの対処を行うべきである。例えば、DNAメチル化プロファイルを患者群8例と健常者群8例の2群で比較したい場合を考える。仮に実験手技の都合で1日に8サンプルしか測定できないならば、初日に患者群4例、健常者群4例を測定し、翌日に各群残りの4例ずつを測定すれば良い。こうすれば、測定日の違いがバッチ効果として測定データに入り込んだとしても、患者群と健常者群の比較には影響しない。逆に悪い例として、初日に患者群8例、翌日に健常者群8例を測定した場合、患者群と健常者群でDNAメチル化割合に有意な差を見出したとしても、それがバッチ効果によるものか、疾患の有無によるものか区別できない。測定を行う前に、この研究で比較したい内容を考え、比較に影響しないような実験デザインを組むことが大事である。
次に②データ取得後の対処とは、統計学的手法によりバッチ効果を取り除く方法である。大別すると2通り存在し、メチル化データからバッチ効果を除く専用の前処理方法と、関連解析時にバッチ効果を考慮する方法がある。事前にバッチ効果を除く方法としては、経験的ベイズ法に基づく方法が頻用されている［44］。関連解析時の補正方法としてよく使われるのは、（2）式にバッチ効果となり得る条件を共変量として加える方法である。例えば、測定日の違いやアレイのロット番号・測定単位をカテゴリデータとして加える。この際、線形混合モデルを採用して、バッチ効果を固定効果ではなくランダム効果として入れることも多い。また、具体的なバッチ効果となり得る変数が無い場合は、代理変数法（surrogate variable analysis, SVA）を用いる手法もある［45］。他にもDNAメチル化アレイに搭載されているコントロールプローブのデータに対してPCAを実施し得られた主成分を共変量に用いる方法や、更にDNAメチル化データの各種共変量を投入した回帰モデルを当てはめ、そこから得られた残差に対してPCAを適用して得られた主成分を共変量に用いるアプローチもある［33］。この方法を応用して、サンプルの収集プロトコールの違いによる影響を補正する方法も開発されている［46］。他にも様々な方法が提案されている。詳細は総説［43］をご覧頂きたい。
６-３．他の関連解析方法6-1では関連解析方法として、線形回帰モデルや一般化線形モデルを用いた方法を紹介したが、それ以外の方法も提案されている。
まず、メタアナリシスである。効果量の小さいDNAメチル化を検出するためには、より多くのサンプルサイズが必要である。そこで、複数の研究グループが解析した結果をメタアナリシスで統合することがある。このアプローチはEWASメタアナリシスやエピゲノムワイドメタアナリシス（epigenome-wide meta-analysis）等と呼ばれる。逆分散法（inverse variance-weighted method, IVW）を用いて、固定効果モデル（fixed effect model）で統合することが多い。ただし、統計学的異質性（I2統計量） が高い場合には、ランダム効果モデル（random effect model）を採用することもある。近年発表されている大規模EWAS論文は、一般化線形モデルから得た効果量をメタアナリシスで統合したEWASメタアナリシスによって行われている［8, 32, 38］。図7で紹介した結果も、EWASメタアナリシスによるものである。この論文ではIVWによって24のコホートの結果を統合することで、8,825例という大規模のEWASを実現した。
一般化線形モデルとは異なるアプローチでDNAメチル化と形質の関連を評価する方法も開発されている。DMAの方法として、例えばD3M（detection of differential distributions of methylation levels）［47］が挙げられる。D3Mはメチル化量を分布値データ（distribution-valued data）として表現し、分布の違いを2群間で比較する手法である。他にもDM-BLD（differential methylation detection using a hierarchical Bayesian model exploiting local dependency）［48］は階層ベイズ法を用いて2群間の比較を行う。これらの手法は共変量を組み込んでモデル化することは想定されていないため、事前にβ値（M値）から共変量の影響を回帰モデル等で取り除き、得られた残差で比較する必要がある［47］。
また、細胞組成を考慮した新しい関連解析方法も提案されている。6-2（2）で紹介した方法は、細胞組成を調整することで、白血球の“平均的な”DNAメチル化と形質の関連を評価していた。しかし、白血球の各細胞種におけるDNAメチル化と形質の関連を評価したいこともある。このような問題を解決するため、Takeuchiらは、非線形リッジ回帰を用いた方法を提案している［49］。
これまでに紹介した解析は、DNAメチル化サイト単位の解析であったが、DNAメチル化データを基に形質と関連する生物学的パスウェイを同定したいことも多い。遺伝子発現データ解析の場合、遺伝子セット解析を用いて生物学的パスウェイの解析を行うことが多く、GSEA （Gene Set Enrichment Analysis）［50］やDAVID （the database for annotation, visualization and integrated discovery）［51］が頻用されていた。これらの方法をDNAメチル化データに応用することが期待されるが、直近に遺伝子が存在しないDNAメチル化サイトも多く存在する。そこでWahlらの研究では、DNAメチル化サイトに遺伝子のアノテーションを付与する距離の条件を複数用意して遺伝子セット解析を実施した［8］。他にもDNAメチル化データ用のパスウェイ探索手法の開発が進められており、更なる手法の提案が望まれる。
以上のようにDNAメチル化アレイデータ解析からは未だ多くの課題が見いだされ、方法論の開発が進められている。
EWASを行い、有意水準に達するP値を示したDNAメチル化サイトを検出でき、その関連はこれまでに報告の無い新規のものだったとする。めでたく論文化といきたいが、そこで研究完了ではない。有意なDNAメチル化サイトを検出した後の重要な2つの工程を以下に記載する。
（1）再現性の確認EWASの結果、新規のDNAメチル化サイトが検出された論文を雑誌に投稿すると、レフリーから独立したサンプルで再現性（replication）を確認するように求められることが多い。そこで、独立したサンプルでもDNAメチル化や形質情報の測定を行い、関連解析を実施する必要がある。とはいえ、DNAメチル化アレイは高価であり、サンプルも有限である。再現性を確認したいDNAメチル化サイトは限られているため、ゲノム全体のDNAメチル化情報を取得するのではなく、狙ったDNAメチル化サイトの領域のみを測定するという方法もある。例えば、パイロシークエンシングを使って特定のサイトのDNAメチル化割合を測定すれば、DNAメチル化アレイで測定するよりは安価に済むだろう。サンプルサイズには限りがあるため、再現性を確認するためのサンプルが無い場合は、公開データを活用する方法がある。公開データの説明は8章で行うが、もし再現性の確認に必要なデータが公開されていれば、追加のデータを取得する必要なく研究を進めることができる。再現性確認の際の有意水準には、ボンフェローニ補正を用いることが多い。5%を検証に回したDNAメチル化サイトの数で除算した値を再現性確認の有意水準として採用する。
（2）同定したDNAメチル化サイトの意味付け同定したDNAメチル化サイトの生物学的・医学的意味付けも重要である。すなわち、なぜそのDNAメチル化サイトが形質と関連するのか、調査する必要がある。最低でも、検出したDNAメチル化サイトがヒトゲノム上のどの位置にあるのか、近隣の遺伝子との位置関係はどうなっているのかといったアノテーション情報を調べるべきであろう。Illumina社が公式サイトでDNAメチル化アレイのアノテーション情報を公開しており、これを活用する手がある（https://jp.support.illumina.com/array/downloads.html）。RのBioconductorからアノテーションを取得することも可能である。HumanMethylation450なら、IlluminaHumanMethylation450kanno.ilmn12.hg19, MethylationEPICならIlluminaHumanMethylationEPICanno.ilm10b4.hg19のRパッケージから取得できる。また、SNP等の塩基配列データで活用されるANNOVAR［52］やSnpEff［53］といたアノテーション付与ツールも有効である。
SNPをはじめ、生まれ持ってのDNA配列の違いがDNAメチル化状態に個人差をもたらすことも知られている［54, 55］。例えばヒトゲノム配列で通常ならCpGサイトがある箇所で、シトシン（C）がアデニン（A）に置き換わっていた場合、CpGサイトでなくなるためDNAメチル化が起こらなくなる。また、DNA配列の違いがその近隣にあるメチル化サイトのDNAメチル化状態と関連することもある。DNAメチル化状態と関連のあるSNPはmQTL（methylation quantitative loci） SNPと呼ばれる。仮にEWASで形質と関連するDNAメチル化サイトを見つけた場合、8章で紹介するiMethyl等のmQTLのデータベースを用いてそのサイトと関連するmQTL SNPが無いか調査すると良い。mQTL SNPが見つかった場合、そのSNPがDNAメチル化サイトを介して形質に影響をもたらしているのではないかという可能性も検討する必要がある。mQTL SNPを探索することはDNAメチル化の機能を考察する上で重要である。
ここまでの話を踏まえて、例として筆者らが実施した血中レジスチン濃度のEWASを挙げる［13］。レジスチンは、タンパク質の一種であり主に単球で発現する。血中のレジスチン濃度は生活習慣病のリスクマーカーになると言われている。そこでレジスチン濃度の調節機構解明のため、DNAメチル化に焦点を当てて、白血球由来のメチル化データでEWASを実施した。その結果、レジスチン濃度と関連するDNAメチル化サイト（cg02346997, cg22322184) の同定に成功した。これらのメチル化サイトの機能を調べるためアノテーションを確認したところ、どちらもレジスチンタンパク質をコードする遺伝子（RETN）のプロモーター領域にあることがわかった。そこでこれらのサイトは、図3で紹介した機構によってRETN遺伝子の発現調節を単球中で行い、その結果血中のレジスチンタンパク質濃度に寄与していると考えられる。更に、これらのDNAメチル化サイトと同じプロモーター領域にSNP （rs3219175, rs34861192) が存在する。これらのSNPは血中のレジスチン濃度と関連することが以前よりわかっていたが、どのようにして濃度を調節していたかは不明であった［56］。EWASで同定したDNAメチル化サイトとこれらのSNPの有意な関連を確認できたことから、SNPがDNAメチル化を介してレジスチンの遺伝子発現を制御し、その結果血中のタンパク質濃度に影響したと考察した。以上のように、DNAメチル化データを取得した組織が（単球を含んだ）白血球であることを踏まえつつ、同定したメチル化サイトのアノテーションおよび、SNPの影響も考慮することで、同定したメチル化サイトの生物学的意味を検討することが可能となった。
これまでにEWASの実施に有益なデータベースが開発・公開されている。本章では、3種類の観点でEWASと関係のあるデータベースを紹介する。最初に（1）EWAS論文やその結果を登録するデータベース、次に（2）個体レベルのDNAメチル化アレイデータを登録しているデータベース、最後に（3）日本人を対象としたDNAメチル化データベースを紹介する。これらを活用することで、EWASで検出したDNAメチル化サイトの意味付けや、自身がデータを保有していなくてもDNAメチル化アレイデータを活用した研究が実施可能になる。
（1）EWAS論文やその結果を登録するデータベース図4で紹介した通り、これまでに多数のEWAS論文が発表され、様々な形質関連メチル化サイトが報告されてきた。結果、研究者が新たにEWASを行い、有意に形質と関連するDNAメチル化サイトを見つけられたとしても、それが過去に報告されているかどうかを知るには膨大な既報のEWAS論文を調査する必要が出てきた。このような問題を解決すべく開発されたEWASのデータベースがEWAS Atlasである（https://ngdc.cncb.ac.cn/ewas/atlas/）［57］。このデータベースは、China National Center for Bioinformation （CNCB）傘下のNational Genomics Data Center （NGDC）によって開発・運用されており、現在も更新され続けている。このデータベースには、EWAS論文で報告されている形質とDNAメチル化サイトの有意な関連解析結果が登録されている。EWAS Atlasの用途は多岐に亘るが、例えば上記で記載した通り、「関連解析で有意に形質と関連していたDNAメチル化サイトは、これまでに同様の報告がされているか？」を調べることができる。また、「目的とする形質を対象としたEWAS論文があるか？その形質との関連が報告されたDNAメチル化サイトは何か？」、「目的とするDNAメチル化サイトは何らかの形質との関連が報告されているか？関連を報告している論文は何か？」などを調べることも可能である。類似したデータベースとしてEWASdb （http://www.bioapp.org/ewasdb/）というデータベースも存在する［58］。こちらは、Harbin Medical Universityによって開発・運用されており、EWAS Atlasと同様の使い方ができる。しかし、このデータベースは2018年の途中から更新が行われていないようである。EWAS Catalog （http://www.ewascatalog.org/）というデータベースもある［59］。こちらはUniversity of Bristol のMRC Integrative Epidemiology Unit （MRC-IEU）によって開発・維持されており、現在β版である。こちらもEWAS Atlasと同様の使い方ができる。
（2）個体レベルのDNAメチル化アレイデータを登録しているデータベース次に個体レベルのDNAメチル化アレイデータを登録しているデータベースを紹介する。研究に用いられたデータは、論文掲載後に公開・共有されるべきという意見があり、投稿先の雑誌によっては研究に用いたオミクスデータの公開を求めることがある。DNAメチル化アレイデータもその対象となっており、これまでに多くのデータが公開されてきた。表2にDNAメチル化アレイのデータを公開しているデータベースの一覧を示した。ここに記載した殆どはDNAメチル化アレイ以外のオミクスデータも公開している。GEOやArrayExpressは、マイクロアレイや次世代シークエンサー等で取得したオミクスデータベースを公開するためのレポジトリである。論文に用いられた多数のデータがこれらのデータベースに登録・公開されている。GDC Data portalは、NIHが運営するデータベースで、The Cancer Genome Atlas （TCGA）等で取得されたがん組織・がん細胞から取得したオミクスデータを公開しているデータベースである。ENCODEは、ヒトや他の生物の組織・細胞から様々なエピゲノム情報を取得するプロジェクトで、そこで取得されたデータを公開している。このように様々なレポジトリやプロジェクトが個体レベルのDNAメチル化アレイデータを公開している。しかしながら、それぞれのデータベース毎に生データからの前処理の仕方や、データの管理方針が異なるため、統合した解析が難しいという課題があった。そこで、GEO, ArrayExpress, TCGA, ENCODEのデータベースで登録されたDNAメチル化アレイのデータを収集し、共通の方法で前処理して公開しているデータベースEWAS Data Hub （https://ngdc.cncb.ac.cn/ewas/datahub/index）［60］ が開発された。EWAS Data Hub には、2021年8月現在Illumina HumanMethylation450, MethylationEPICのデータが登録されている。EWAS Data HubはEWAS Atlasと同様、NGDCによって開発・運営されている。このデータベースでは、共通の前処理が行われたDNAメチル化アレイデータに加え、年齢、性別、形質等のデータもあわせて取得できる。
| データベース名 | 運営 | URL |
| Gene Expression Omnibus (GEO) | National Center for Biotechnology Information (NCBI) | https://www.ncbi.nlm.nih.gov/geo/ |
| ArrayExpress | European Bioinformatics Institute (EBI) | https://www.ebi.ac.uk/arrayexpress/ |
| Genomic Data Commons (GDC) Data Portal | National Institute of Health (NIH) | https://portal.gdc.cancer.gov/ |
| Encyclopedia of DNA Elements (ENCODE) | ENCODE Consortium | https://www.encodeproject.org/ |
| EWAS Data Hub | National Genomics Data Center (NGDC) | https://ngdc.cncb.ac.cn/ewas/datahub |
以上のデータベースを使えば、個体レベルのDNAメチル化アレイデータが取得でき、研究者が望む条件で解析できる。EWASメタアナリシスや、EWAS以外の解析にも利用できるだろう。これらのデータベースの中には、DNAメチル化データに加え、同一サンプル由来の他のオミクスデータ（遺伝子発現量等）も公開している研究もあるため、トランスオミクス解析も実施可能である。
個体レベルの公開データを用いる場合、次のような事項に注意する必要がある。これらのデータベースから取得できるDNAメチル化アレイデータは、2章で紹介したマトリックスで公開されていることもあれば、生データ（idatファイル）で公開されていることもある。生データの場合、自分で前処理を行い、マトリックスに変換する必要がある。前処理後のマトリックスデータが公開されているとしても、どのような前処理が適用されているか調べ、その処理内容が自分の行いたい研究に悪影響を与える可能性がないかチェックする必要がある。通常、元の研究で用いられた交絡因子の情報は、DNAメチル化アレイデータとともに、データベース内で公開されている。ただし、様々な理由から交絡因子が公開されていないこともある。また、その研究で使用されていない交絡因子は公開されていないことが多い。公開データを活用する場合、自分が必要とする交絡因子が全て含まれているか、データ取得時に確認する必要がある。
（3）日本人を対象としたDNAメチル化データベース最後に、日本人を対象としたDNAメチル化データベースiMethyl （http://imethyl.iwate-megabank.org/）を紹介する。iMethylはいわて東北メディカル・メガバンク機構によって開発・維持されているデータベースであり、日本人約100名のDNAメチル化情報、遺伝子発現情報、SNP等の多型情報をゲノムブラウザで公開している［61］。iMethylでは個人ごとのデータは公開されていないものの、オミクスデータは人種によって異なるプロファイルを持ち得るため、日本人のデータを調べる場合に有益である。DNAメチル化情報、遺伝子発現情報が、CD4+T細胞、単球、好中球（顆粒球の一種）毎に、WGBSとRNA-seqで取得されており、これらのDNAメチル化情報、遺伝子発現情報、遺伝子多型情報の互いの関連解析結果を要約統計量としてゲノムブラウザ上に公開している。1章で概説した通りDNAメチル化は遺伝子発現を制御しているが、このデータベースを使えば、実際に各細胞種でDNAメチル化とどの遺伝子の発現量が関連するのかを調べることができる。また、7-（2）で述べた通りDNA配列の違いがDNAメチル化状態の個人差に影響をもたらすこともある［54, 55］。iMethylでは、前述の通り日本人におけるDNAメチル化情報と遺伝子多型情報の関連解析結果をmQTLとして登録している。iMethylを使えば、日本人の白血球におけるmQTLの探索も行うことができる。
DNAメチル化アレイデータの用途は単純なEWASに限らない。遺伝子発現データの研究で行われてきたように、DNAメチル化アレイデータを用いたアウトカムを予測するシグネチャの開発研究も進められてきた［62, 63, 64］。開発されたシグネチャは、疾患診断や発症リスク予測のための新たなバイオマーカーになると期待されている。中でも特に注目されているのが、エピジェネティクス年齢（epigenetic clock, エピゲノム年齢と訳されることもある）である。これはDNAメチル化アレイデータから推定した年齢である。ここから得られた年齢の推定値は、実年齢（chronological age）ではなく、生体内の年齢（biological age）である。実年齢とエピジェネティクス年齢の差分を評価することで、サンプルが実年齢より老化しているかどうかを評価することができ、その差分値が疾患発症のリスクマーカーになると期待されている。これまで様々なエピジェネティクス年齢が開発されてきたが、中でもHorvathの開発したシグネチャが注目されている［65］。彼は公開されているDNAメチル化データを基に罰則付き回帰を用いて年齢を予測するモデルを開発し、DNAメチル化データから予測した年齢をエピジェネティクス年齢と定義した。このエピジェネティクス年齢と実年齢の差分は、将来的な循環器疾患発症のリスクマーカーになることを示した［66］。他にも様々なエピジェネティクス年齢が開発されており、疾患発症の予知マーカーになることが報告されている。詳細は総説［67］をご覧頂きたい。このようにDNAメチル化アレイデータには、EWAS以外にも多くの用途がある。今後も新たなシグネチャの開発が期待できる。
本稿では、EWASを主軸にDNAメチル化アレイデータを用いた解析の流れを紹介した。DNAメチル化データと形質の関連解析方法だけでなく、交絡因子の除去や、細胞組成推定を用いた解析等、新たな方法の開発は続いている。GWASと同じ様に今後サンプルサイズの大規模化が進み、要約統計量を活用した解析方法も発展していくだろう。
現在のEWASは、測定コストの面から多数のサンプルを必要とするため、DNAメチル化アレイのデータが主流となっているが、将来的にシークエンシング技術がより低コスト化していくことで、WGBSで測定したデータを用いたEWASが主流になると期待される。事実、既にWGBSデータを用いたEWASが日本から発表されている［68, 69］。WGBSデータを用いたEWASでも、DNAメチル化アレイで用いた方法・ノウハウがそのまま流用できる点もあるだろうが、WGBS特有の統計手法も必要になるだろう。EWASの方法論開発は今後より激化していくことだろう。
本稿の執筆にあたり、名古屋大学医学部附属病院 木下文恵先生には、ご多忙の中、細部まで原稿の内容・表現を確認頂き助言を頂いた、この場をお借りして感謝を申し上げたい。本総説でご紹介したレジスチンのEWASは、北名古屋ゲノム疫学研究（KING Study）、九州大学福岡コホート研究（Fukuoka study）の参加者を対象としており、これらのコホート研究参加者及び関係者の皆様に御礼申し上げます。これらの研究でお世話になった、久留米大学医学部 横田充弘先生、山本健先生、自治医科大学医学部 市原佐保子先生、九州大学大学院医学研究院 大中佳三先生、愛知みずほ大学人間科学部 松原達昭先生、岩倉病院 浅野 展行先生他、多くの共同研究者の先生方にこの場を借りて深謝いたします。
 |
中杤 昌弘 名古屋大学大学院工学研究科化学・生物工学専攻を修了。博士（工学）。現在、名古屋大学大学院医学系研究科 総合保健学専攻 実社会情報健康医療学の准教授として講座を運営。遺伝疫学データやリアルワールドデータの解析研究を進めている。 |
