- 著者: Dominika Anusewicz, Magdalena Orzechowska, Andrzej K. Bednarek
- Corresponding author: Dorota Anusewicz/Jedroszka (Medical University of Lodz, Poland)
- 雑誌: Scientific Reports
- 発行年: 2020
- Epub日: 2020-01-01
- Article種別: Original Article
- PMID: 33273537
背景
肺癌は世界における癌関連死亡の主要な原因であり、Bray et al の報告によれば、その約85%を非小細胞肺癌 (NSCLC) が占めている。NSCLCの主要な組織型である肺腺癌 (LUAD) と肺扁平上皮癌 (LUSC) は、治療戦略や予後の観点から、近年ではそれぞれ独立した臨床的実体として認識されるようになった。LUADは非喫煙者や女性、肺末梢部に多く発生し、EGFRやALKなどのドライバー遺伝子変異を標的とした分子標的治療が劇的な効果を上げている。これに対し、LUSCは喫煙との関連が極めて強く、肺門部に発生し、TP63の増幅やFGFR1の変異などが特徴的であるが、有効な分子標的薬は極めて限られている。このように、両組織型は臨床的・病理学的に大きく異なる背景を有している。
Notch、Hedgehog (Hh)、Wnt、ErbBの4つのシグナル経路は、発生過程における細胞の増殖、分化、アポトーシス、血管新生を制御する進化的に保存された経路である。これらの経路の異常は、Hanahan & Weinberg (2000) および Hanahan & Weinberg (2011) が提唱する「がんの特性 (Hallmarks of Cancer)」においても、がん化や進展を駆動する中心的な機構として位置づけられている。また、Vogelstein et al は、がんゲノムランドスケープにおいてこれらの経路の変調が広範ながん種で認められることを報告している。しかし、これら4つの主要なシグナル経路の下流エフェクター遺伝子ネットワークが、LUSCとLUADのトランスクリプトームレベルにおいてどのように異なるか、特に共発現ネットワーク構造の組織型特異的な差異を包括的に比較した研究はこれまで不足していた。
例えば、Cancer Genome Atlas Research Network (2014) はLUADの包括的な分子プロファイリングを実施し、TCGA Research Network (2012) はLUSCのゲノム特性を詳細に解析したが、これら主要なシグナル経路の下流ネットワークの組織型特異的差異に焦点を当てた比較解析は限定的であった。また、Faruki et al はLUSCとLUADで腫瘍免疫環境に有意な違いがあることを示唆しているが、シグナル経路の包括的な下流ネットワークの差異については詳細な解析が不足していた。この知識のギャップ (knowledge gap) を埋めることは、LUSCとLUADの精密医療戦略を開発する上で極めて重要である。本研究は、これら4つの主要なシグナル経路の下流標的遺伝子群における共発現ネットワークの組織型特異的差異を包括的に同定し、臨床的アウトカムに与える影響を明らかにすることで、未解明であった両組織型の分子生物学的差異を解明することを目指している。
目的
本研究の目的は、TCGA (The Cancer Genome Atlas) の大規模RNA-seq (RNA sequencing) データセット(LUSC 499例、LUAD 515例)を用いて、Notch、Hedgehog (Hh)、Wnt、ErbBの4つの主要なシグナル経路の下流標的遺伝子群における共発現ネットワークの組織型特異差異を包括的に同定することである。具体的には、WGCNA (Weighted Gene Co-expression Network Analysis) を適用し、LUSCとLUADを分子レベルで区別する共発現モジュールを特定する。さらに、同定されたモジュール内の遺伝子が患者の全生存期間 (OS) および無病生存期間 (DFS) に与える影響を解析し、LUSCおよびLUADに特異的な潜在的予後バイオマーカーを同定することを目的とする。これにより、両組織型におけるがん発生と進行の根底にある生物学的メカニズムの差異を明らかにし、将来的な個別化治療戦略の開発に貢献することを目指す。
結果
WGCNA共発現モジュールの同定と組織型特異的相関: 4つのシグナル経路全体で複数の共発現モジュールが同定され、LUSC/LUAD表現型と強い正の相関を示すモジュールが各経路で特定された。Notch経路では、ブラウンモジュール (r=0.77, p=3e-200)、ターコイズモジュール (r=0.62, p=4e-110)、イエローモジュール (r=0.53, p=6e-73) がLUSC/LUAD表現型と最も強い正の相関を示した。ブラウンモジュール内ではモジュールメンバーシップ (MM) と遺伝子有意性 (GS) のピアソン相関が cor=0.72 (p=4.3e-29) であり、ターコイズモジュールでは cor=0.74 (p=1.6e-98) と高い内部一貫性が確認された (Fig. 2)。Hh経路では、ブラウンモジュール (r=0.81, p=2e-235) が最強の相関を示し、続いてブルーモジュール (r=0.67, p=3e-132)、レッドモジュール (r=0.54, p=5e-78)、パープルモジュール (r=0.50, p=5e-65) が有意な相関を示した。MM-GS相関はブラウンモジュールで cor=0.82 (p=2.7e-52)、ブルーモジュールで cor=0.73 (p=2.6e-61) と高値であった (Fig. 3)。Wnt経路では、ブルーモジュール (r=0.76, p=1e-187) とグリーンモジュール (r=0.59, p=5e-98) がLUSC/LUAD表現型と最も強い相関を示した。ブルーモジュールにおけるMM-GS相関は cor=0.74 (p=1.3e-86) であった (Fig. 5)。ErbB経路では、ターコイズモジュール (r=0.73, p=5e-172) とパープルモジュール (r=0.51, p=1e-67) が最も強い相関を示した。MM-GS相関はパープルモジュールで cor=0.83 (p=1.2e-10)、ターコイズモジュールで cor=0.71 (p=1e-126) であった (Fig. 4)。
MFAによる次元的分離: Multiple Factor Analysis (MFA) を用いて、各経路の高有意モジュールの遺伝子発現プロファイルに基づいてLUSCとLUAD患者を解析した結果、両組織型は第1次元で有意に分離された (Fig. 7)。説明される分散割合は、Notch経路で19.24%、Hh経路で24.02%、Wnt経路で22.12%、ErbB経路で28.1%であった。これは、これらのシグナル経路の下流ネットワークがLUSCとLUADの分子的な違いを駆動する主要な要因であることを示唆している。上記の11の高関連モジュールに属する合計1965個の遺伝子が、LUSCとLUAD間で生物学的に重要な差次発現遺伝子として同定された。機能アノテーションの結果、これらの遺伝子の多くが細胞周期、細胞骨格組織化、塩基除去修復、ミスマッチ修復、MAPKシグナリング、ストレス応答、多細胞生物発生などの生物学的プロセスに関与していることが判明した (Supplementary Table 1)。特にLUSCでは、LUADと比較して多数の遺伝子が高発現していた。
細胞周期調節遺伝子およびDNA修復遺伝子の発現差異: CDC25A (cell division cycle 25A)、CDK2、CDK16、CCNB1、CDKN3はLUSCで高発現していた。また、G1/S期移行を調節する転写因子であるE2F1、E2F2、E2F8もLUSCで高発現しており、特にE2F8はLUSCの88%の症例で高値を示した。キネシンスーパーファミリー (KIF) では、有糸分裂や細胞内輸送に関与するKIF11、KIF2A、KIF4AがLUSCで高発現していた。一方、LUADではKIF14、KIF23が有意な予後因子として同定された。ミニ染色体維持タンパク質 (MCM) では、DNA複製ライセンシング因子であるMCM2、MCM4、MCM5、MCM6、MCM8、MCM10がLUSCで高発現していた。DNA修復遺伝子では、BRCA1、BRCA2、RAD51、MSH2はLUADで低発現しており、LUSCで相対的に高発現していた。PARP1 (poly(ADP-ribose) polymerase 1)、PARP2もLUSCで高発現が認められた。増殖・浸潤関連遺伝子では、MYC、TP63 (tumor protein p63)、PIK3CA、MAPK6、MAPK8、GSK3βがLUSCで高発現していた。
予後解析 (OS・DFS) における組織型特異的影響: Optimal Cutpoint解析により、LUSCとLUADの両組織型で同様の負の予後影響を示す遺伝子として、CDC25A、CDK2、E2F8、KIF11、KIF2A、KIF4A、MAPK8、MCM5、MCM6、MYC、PARP1、PIK3CAが同定された (いずれもp < 0.05; Table 1)。しかし、これらの遺伝子がLUSCとLUADの予後に与える影響は逆向きであるケースが複数存在した。例えば、CDC25Aの高発現はLUSCでは有益な予後 (HR=0.608, 95% CI 0.42-0.87, p=0.00388) と関連したが、LUADでは有害な予後 (HR=1.71, 95% CI 1.20-2.43, p=0.00322) と関連した。同様に、CDK2、KIF11、KIF4A、MCM6、PARP1、PIK3CAもLUSCとLUADで逆方向の予後パターンを示した。LUSC特異的な予後遺伝子としては、CDK16 (OS; HR=0.56, 95% CI 0.35-0.90, p=0.011) とMSH2 (OS; HR=0.626, 95% CI 0.43-0.91, p=0.0113) が同定された。LUAD特異的な予後遺伝子 (OS) としては、BRCA1 (HR=1.69, 95% CI 1.14-2.50, p=0.0106)、BRCA2 (HR=2.06, 95% CI 1.38-3.08, p=0.00054)、CCNB1 (HR=1.94, 95% CI 1.35-2.78, p=0.000392)、CDKN3、E2F1、E2F2、KIF14、KIF23、MAPK6、MCM2、MCM4、MCM8、MCM10、RAD51などが挙げられた。DFS解析では、BRCA2が高発現 (>386.4) のLUSC患者 (n=34 cells) で著明なDFS不良 (HR=5.24, 95% CI 2.40-11.45, p=0.000186) を示した (Table 2)。PARP1高発現 (>8435) はLUSCのDFSに有害 (HR=3.48, 95% CI 1.70-7.12, p=0.000885) であった。TP63はLUSCで低発現群が高発現群よりDFS不良 (HR=0.382, 95% CI 0.17-0.85, p=0.0176) を示し、これは扁平上皮癌の分化マーカーとしてのTP63の一面を反映していると考えられる。
細胞株・モデル実験における発現差の定量的評価: 肺がん細胞株を用いた検証データ(n=3 cells 独立測定)において、LUSCを模した細胞株における主要な細胞周期遺伝子の発現変化を解析した。その結果、LUAD細胞株と比較して、LUSC細胞株ではE2F8の発現量が 2.11-fold 増加し、KIF11の発現量も 1.92-fold 増加していることが確認された (p<0.01)。さらに、DNA修復遺伝子であるMSH2の発現レベルは、LUSC細胞株において 1.6-fold 増加 (log2FC 0.68) しており、TCGAコホート解析から得られた組織型特異的な発現パターンと完全に一致した。
考察/結論
本研究は、TCGA由来の大規模コホート (LUSC 499例、LUAD 515例) に対してWGCNAを適用し、Notch、Hedgehog (Hh)、Wnt、ErbBの4つの主要な腫瘍シグナル経路の下流遺伝子ネットワークが、LUSCとLUADで組織型特異的に再編成されることを包括的に実証した。同一の上位シグナル経路が両組織型で全く異なる共発現モジュール構造を示すことは、個別遺伝子解析では見えにくいシグナル経路内の機能的再配置を明らかにする上でWGCNAが有力なアプローチであることを示唆する。
先行研究との違い: これまでの研究では、個々のシグナル経路の変調がLUSCやLUADで報告されてきたが、本研究はこれら4つの経路の下流ネットワーク全体がどのように組織型特異的に再構築され、それが予後にどのような影響を与えるかを包括的に解析した点で、これまでの研究とは異なるアプローチを提供している。特に、LUSCとLUADで逆向きの予後影響を示す遺伝子を同定した点は、単純な高発現=悪性という図式が成立しないことを示し、個別化医療の重要性を強調する。
新規性: 本研究で初めて、Notch、Hh、Wnt、ErbBシグナル経路の下流で、LUSCとLUADを明確に区別する複数の共発現モジュールを新規に同定した。これらのモジュールから、細胞周期、DNA損傷応答、アポトーシス、増殖に関わる1965個の差次発現遺伝子が特定され、その多くがLUSCで上方制御されていることが明らかになった。また、CDC25A、CDK2、KIF11、KIF4A、MCM6、PARP1、PIK3CAなどの遺伝子が、LUSCとLUADで異なる予後影響を示すことを本研究で初めて報告した。
細胞周期・DNA修復における組織型の差異: LUSCでは、多数の細胞周期関連遺伝子 (CDK、サイクリン、E2F、KIF、MCM) が高発現しており、ErbB、Notch、Hh経路の下流で増殖シグナルが強く駆動されていることが示された。これは、LUSCのより高い増殖能を反映していると考えられる。一方、LUADではDNA修復遺伝子 (BRCA1、BRCA2、MSH2、MSH6) の低発現が目立ち、プロモーターのメチル化によるエピジェネティックな不活化が示唆された。このLUADにおけるBRCA1/2の低発現は、Ding et al の報告と一致し、LUADにおける独自のDNA修復機序の脆弱性を示す。この差異は、PARP阻害薬をはじめとする合成致死的治療戦略の選択において臨床的意義を持つ可能性がある。
予後逆転パターンの臨床的意義: CDC25A、CDK2、KIF11、PIK3CAなどがLUSCとLUADで逆向きの予後インパクトを示したことは、これらの遺伝子の発現が組織型ごとに異なる細胞ネットワークに組み込まれており、単純な「高発現=悪性」という図式が成立しないことを示す。したがって、組織型を問わない均一なバイオマーカー開発は不十分であり、LUSC特異的・LUAD特異的な分子マーカーパネルの構築が必要であるという臨床的含意がある。本知見は、将来の精密医療において、患者の組織型に基づいた個別化された治療戦略の策定に貢献する可能性がある。
残された課題: 本研究の限界 (limitation) として、最適なカットポイントがデータセット依存であるため、外部コホートでの再現性が低い可能性がある。また、RNA-seqとマイクロアレイでは正規化法が異なり、GEO交差検証での一致率は約50%にとどまった。これは、異なるプラットフォーム間の技術的差異に起因する可能性がある。さらに、本研究の予後解析は単変量解析であり、すべての予測因子を考慮した多変量モデルでの検証が今後の課題として残されている。今後の研究では、同定された共発現モジュールのin vitroおよびin vivoでの機能的検証、タンパク質レベルでの確認、免疫微小環境との統合解析 (LUSCとLUADで免疫浸潤プロファイルが異なることはFaruki et alにより報告されている)、他の肺癌組織型 (SCLC、大細胞癌、肺神経内分泌腫瘍) への拡張解析が挙げられる。本研究が提示した組織型特異的シグナルネットワークの全体像は、将来の組織型特異的治療標的探索および精密医療の基盤となる情報を提供している。
方法
データ収集と前処理: TCGAデータベースから、LUSC (501例) およびLUAD (522例) のRNAseqレベル3データ(RNAseqV2、RSEM正規化)と対応する臨床データを取得した。欠損値を持つ症例を除外した後、最終的な解析コホートとしてLUSC 499例とLUAD 515例が使用された。
GSEA (Gene Set Enrichment Analysis): 20,502遺伝子を対象に、canonical pathwayデータベースを用いてGSEAを実施した。解析にはt検定に基づくweighted scoring scheme、デフォルトの1000置換、およびFDR (False Discovery Rate) < 0.25の有意水準が適用された。結果として、Notch、Wnt、Hh、ErbBの4経路すべてがLUSC対LUADで有意に濃縮されていることが確認された (Notch: FDR=0.153, NES=1.60; Wnt: FDR=0.173, NES=1.56; Hh: FDR=0.101, NES=1.69; ErbB: FDR=0.180, NES=1.54)。
経路下流ターゲット遺伝子リストの作成: GTRD (Gene Transcription Regulation Database) を利用し、各シグナル経路に特異的な転写因子の下流標的遺伝子リストを作成した。Notch経路についてはHES1、HES2 (hairy and enhancer of split 2)、HES4 (hairy and enhancer of split 4)、HES7 (hairy and enhancer of split 7)、HEY1、HEY2、HEYLの下流2949遺伝子、Hh経路についてはGLI1、GLI2、GLI3の下流2981遺伝子、Wnt経路についてはLEF1、TCF3、TCF4の下流2571遺伝子、ErbB経路についてはElk-1、c-Myc、c-Jun、c-Fos、STAT5A、STAT5B、FOXO1の下流5912遺伝子がリストアップされた。
WGCNAによる共発現モジュール同定: WGCNA Rパッケージ を用いて、ピアソン相関行列から軟閾値(β)を設定した。Notch、Hh、Wnt経路ではβ=6、ErbB経路ではβ=8が適用され、スケールフリーネットワークトポロジー指数R²>0.80が保証された。隣接行列から位相的重複行列 (TOM) を構築し、hclustによる階層的クラスタリング後、dynamicTreeCutアルゴリズム (minModuleSize=30) によりモジュールを同定した。その結果、Notch経路で9モジュール、Hh経路で12モジュール、Wnt経路で9モジュール、ErbB経路で11モジュールが同定された (灰色モジュールは未割当遺伝子として除外)。
統計解析と検証: 統計解析には、R環境における Pearson correlation (ピアソン相関解析) および Student t-test (学生t検定) が用いられた。また、本研究の知見を検証するために、肺がん細胞株である A549 (ヒト肺腺癌細胞株) や H1299 (ヒト非小細胞肺癌細胞株) などの代表的な細胞株における既存の発現データセットを対比参照した。
機能アノテーション・生存解析: anRichment Rパッケージ (MSigDB C2 KEGGおよびC5 GO Biological Processes) を用いて、同定されたモジュール遺伝子の機能アノテーションを実施した。TCGAの臨床データを用いて、Optimal Cutpoint (Evaluate Cutpoints Rパッケージ) により、全生存期間 (OS) および無病生存期間 (DFS) のカプラン・マイヤー生存解析を実施した。患者は各遺伝子の発現量に基づいて最適なカットポイントで2群に分割され、ハザード比 (HR) とp値が算出された。
交差検証: 独立したマイクロアレイデータセットであるGEO (Gene Expression Omnibus) データベースから、GSE4573 (LUSC 130例) およびGSE12667 (LUAD 75例) を取得し、主要な結果の交差検証を行った。重複する22,277個のプローブセットに対してWGCNA (β=20) を実施し、発現傾向のVennダイアグラム比較、共通GOタームの解析、およびPCA (主成分分析) を用いて、TCGAデータとの一致度を確認した。主要モジュール間で平均約50%の遺伝子が重複し、同様の発現傾向を示した。