Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization

―Part VIII. Survival Analysis―

Abstract

生命科学分野で取得されるデータ集合は、雑多（ヘテロ）な構造になり、ヘテロなデータ構造を扱える理論的な枠組みがもとめられている。本連載では、汎用的なヘテロバイオデータの解析手法である行列・テンソル分解を紹介していく。第8回では、イベント発生の時系列性に関する解析である生存時間解析と行列・テンソル分解の関係性について議論する。

生存時間解析とは

生存時間解析は、あるイベントが発生するまでの時間（time-to-event）を扱う統計手法である。対象は、死亡や疾患の発症・再発、治療の中断などのイベントであり、各個人について「観測時間」と「イベント発生の有無（打ち切りを含む）」が記録される。生存時間解析の目的は、イベント発生までの時間分布を推定するとともに、年齢や臨床指標などの共変量が将来の生存見通し（予後）にどのような影響を与えるかを評価することである（図1）。

図１：生存時間解析が扱うtime-to-eventデータの例（Rのsurvivalパッケージに含まれるlungデータセット）

各患者について、生存時間（time）、イベント状態（status）、および年齢や性別などの共変量が記録される。なお、このデータではstatusが一般的な0＝打ち切り、1＝イベント発生ではなく、1＝打ち切り、2＝死亡と符号化されている点に留意されたい。

このようなデータに対しては、これまでの連載［1, 2, 3, 4, 5, 6, 7］で扱った行列分解手法、特に第6回［6］で扱った自己回帰モデルのような時系列解析手法を思い浮かべる読者もいるかもしれない。しかし、本稿で取り上げる生存時間解析（survival analysis）は、以下の理由により一般的な時系列データ解析での議論とは大きく異なる理論的枠組みで発展していることから、ここでは独立して取り上げることにする。

●打ち切り（censoring）の存在：一部の被験者は、観察期間内にイベントを経験しないまま追跡が終了したり、途中で脱落してイベントの有無が不明となる（図1）。こうした不完全な情報を含んだデータを統計的に扱う必要があるため、専用の解析手法が求められる。

●イベントが一度しか起こらない前提：生存時間解析では、死亡や初回再発といった1人につき一度しか起こらないイベントを対象とする。繰り返し起こる「再発イベント」も存在するが、それには再発イベント解析という別の枠組みがあり、本稿では扱わない。ただし、そのようなデータであっても最初の1回だけに注目するのであれば、通常の生存時間解析の範疇である。

●時間の扱い方：生存時間解析では、各個人の計測開始時点を原点とした上で、イベント発生時刻は原点からの経過時刻として扱われる。イベントが発生するまでの詳細な経過データ（バイタルの連続的な時系列など）は手に入らないことが多い。仮にこのような生存時間データで個人を行、イベント発生時刻を列、として｛0, 1｝の二値行列として表したとしても、各個人について1（イベント発生）があるのは最大1箇所のみであり、それ以外は全て0となるため、極めて疎で特殊な構造の行列となることから、このような行列データに対して行列分解など、多変量解析の手法を適用してもほとんど意味をなさない。

このような扱う問題の特殊性ゆえ生存時間解析は、昨今様々なデータサイエンス分野で行われてきた「従来の問題を最新の機械学習モデル（例：深層学習モデル）で単純に置き換えることで高い性能を示す」アプローチが取りづらく、独自の進化を遂げてきた（ただし、後述するように近年では機械学習的アプローチとの融合も少しずつではあるが進んでいる）。

生存時間解析概論

ここでは、生存時間解析における理論的な基本要素を簡潔にまとめる。生存時間解析では、時間の経過とともにイベント（例：死亡、再発、治療中断など）が発生するため、時間tを変数とする関数が多数登場する。その中でも特に基本となる重要な5つの関数を以下に紹介する［8, 9, 10, 11, 12, 13, 14, 15, 16, 17］。

1．確率密度関数f（t）：イベントが発生する時刻の確率密度を表す関数。ある区間［a, b］におけるイベント発生確率は、f（t）をその区間で積分した面積に対応する。

2．累積分布関数$F(t)={\int }_0^{t}f(u)\mathit{du}$：計測開始時刻0から時刻tまでにfを積分した値。「時刻tまでにイベントが起きる確率」を意味する。

3．生存関数S（t）＝1－F（t）：F（t）とは反対に、時刻tまでにイベントが起きていない確率。計測開始時刻0から∞までにSを積分したものは平均生存時間となる。

4．ハザード関数$h(t)=\frac{f(t)}{S(t)}$：時刻tまでイベントが起きていなかった人が、その瞬間にイベントを起こす「条件付き確率密度」。瞬間的なイベント発生の“勢い”を定量化した尺度である。

5．累積ハザード関数$H(t)={\int }_0^{t}h(u)\mathit{du}$：計測開始時刻0から時刻tまでにhを積分した値。

これら5つの関数の定義域・値域および関係性は図2にまとめた。これらは互いに変換可能であり、どれかがわかると全てわかるようになっている。

図２：生存時間解析における5つの基本関数の関係

確率密度関数、累積分布関数、生存関数、ハザード関数、累積ハザード関数の定義および相互関係を示す。これらは相互に変換可能であり、特にS（t）＝exp（－H（t））は生存時間解析の理論的基盤となる重要な関係式である。

ハザード関数h（t）は、生存時間解析においておそらくもっとも重要で、かつ直感的に意味を理解しづらい関数かもしれない。その意味は「これまでイベントが起こっていなかった個体が、次の瞬間にイベントを起こす可能性の高さ（密度）」である。確率の比として定義されているため、値は0以上の任意の非負実数を取り得るため、この尺度自体は確率値ではないことに留意していただきたい。

特に重要な関係式として、以下がある。

  

$$S(t)=\exp (-H(t))$$（1）

これは以下のような式変形で求められる。まずF（t）とS（t）の定義より、以下のf（t）とS（t）の関係式が得られる。

  

$$f(t)=\frac{d}{\mathit{dt}}F(t)=-\frac{d}{\mathit{dt}}S(t)$$

これを以下のように、h（t）の定義の式に代入する。

  

$$h(t)=\frac{f(t)}{S(t)}=-\frac{S^{\text{'}}(t)}{S(t)}=-\frac{d}{\mathit{dt}}\log (S(t))$$

  

$${\int }_0^{t}h(u)\mathit{du}=-\log (S(t))$$

  

$$-{\int }_0^{t}h(u)\mathit{du}=\log (S(t))$$

  

$$S(t)=\exp \left(-{\int }_0^{t}h(u)\mathit{du}\right)=\exp \left(-H(t)\right)$$

この式（1）は、後述するKaplan-Meier法やNelson-Aalen法における生存時間関数から累積ハザードへの変換や、Cox比例ハザードモデルの生存関数の導出など、生存時間解析においての様々な場面で利用される。

Kaplan-Meier法

実際にデータから上記の関数の具体的な形を推定する場合、計測時間は離散値であり、またデータには打ち切りが含まれることから、様々な工夫がなされる。例えば、データから直接生存関数をノンパラメトリックに推定する方法として、以下のKaplan-Meier法（Kaplan-Meier estimator）がある［8, 9, 10, 11, 12, 13, 14, 15, 16, 17］。Kaplan-Meier法では以下の式に従い、各イベント時刻tごとに、生存確率を逐次的に更新していく。

  

$$\widehat{S(t)}=\prod _{t_j\le t}\left(1-\frac{d_j}{n_j}\right)$$（2）

ただし、t_jは時刻t以前のイベント発生時刻、d_jは時刻t_jにイベントが発生した人数（通常1）、n_jは時刻t_j直前までイベントが発生していない人数である。打ち切りがあった場合は、その分n_jに計上されないため、各時刻ごとに手に入るデータだけが曲線の推定に利用される形となる。この値をステップ関数として描画したものが、Kaplan-Meier曲線である（図3）。ところどころに見られる縦線（ヒゲ）は、打ち切りが発生した時刻を表す（図3）。式（2）を見てわかるように、時間が過ぎていくに従い分母n_jは小さくなっていくため、推定値の不確実性（ばらつき）は増加する。そのため、Greenwoodの公式などを用いて、95％信頼区間を推定する場合もある。Kaplan-Meier法で得られて生存関数は、式（1）により直ちに累積ハザード関数にも変換可能である。

  

$$\widehat{H(t)}=-\log (\widehat{S(t)})$$（3）

図３：生存時間データの層別化と生存関数の推定

群（例：治療群・対照群）ごとに生存関数を推定したKaplan-Meier曲線の例を示す。群間差の検定にはログランク検定が用いられる。

Nelson-Aalen法

Kaplan-Meier法に類似した累積ハザード関数推定手法として、Nelson-Aalen法がある。テイラー近似（e^x≈1＋x）から、

  

$$\exp \left(-\frac{d_j}{n_j}\right)\mathrm{\approx }1-\frac{d_j}{n_j}$$

が得られるため、式（2）に代入することで、Nelson-Aalen法による生存関数は、以下のようにして得られる。

  

$$\widehat{S(t)}=\prod _{t_j\le t}\exp \left(-\frac{d_j}{n_j}\right)$$（4）

同様に、式（3）に式（4）を代入することで、

  

$$\widehat{H(t)}=-\log \left(\prod _{t_j\le t}\exp \left(-\frac{d_j}{n_j}\right)\right)=-\sum _{t_j\le t}\log \left(\exp \left(-\frac{d_j}{n_j}\right)\right)$$

となり、以下のようなシンプルな式が得られるため、これを利用しデータから直接累積ハザード関数をノンパラメトリックに推定する。

  

$$\widehat{H(t)}=\sum _{t_j\le t}\frac{d_j}{n_j}$$（5）

サンプルサイズが小さい場合は、Nelson-Aalen法がKaplan-Meier法よりも推定精度が高いと言われている［9, 10］。

ログランク検定

試薬投与群とプラセボ群など、ある基準に従ってデータを2群に層別化し、それぞれについてKaplan-Meier曲線を描くことを考える（図3）。両群の曲線が大きく乖離している場合、それら群を分けた基準がイベントの発生に影響している可能性が示唆される。この差が偶然によるものかどうかを統計学的に検定する方法としてログランク検定（log-rank test）がある［8, 9, 10, 11, 12, 13, 14, 15, 16, 17］。

ログランク検定では、イベントが発生した各時刻t_jに置いて、群ごとの観測イベント数と期待イベント数を比較する。各時刻におけるリスク集合の人数を群1をn_1j、群2をn_2j、総数をn_j＝n_1j＋n_2jとし、観測イベント数を群1をO_1j、群2をO_2j、総数をd_j＝O_1j＋O_2jとすると、帰無仮説（両群でハザードが等しい）の下で群1の期待イベント数は

  

$$E_{1j}=d_j\frac{n_{1j}}{n_j}$$

と表される。各時刻におけるO_1j－E_1jを全時刻で合計し、分散で標準化した統計量

  

$$Z=\frac{{\sum }_j(O_{1j}-E_{1j})}{\sqrt{{\sum }_j{V}_{1j}}}$$

は漸近的に標準正規分布に従い、その二乗は自由度1のχ²分布に従う。ここで、2群比較では独立に変動できる成分は1つであるため自由度は1となっている。

Cox比例ハザードモデル

生存時間解析において、患者ごとの属性情報（年齢、性別、ステージなど）と生存との関係性を解析したい場合、これら共変量を説明変数として組み込んだ、何らかの回帰モデルが必要となる。これまでに様々な回帰モデルが提案されているが、最も広く用いられているのがCox比例ハザードモデル（Cox proportional hazards model、Cox回帰［8, 9, 10, 11, 12, 13, 14, 15, 16, 17］）である（図4）。ここではその基本的な考え方を紹介する。

図４：Cox比例ハザードモデルによる共変量効果の評価と層別化の概念図

生存時間データに、年齢やステージなどの臨床情報、あるいは遺伝子発現などのオミックスデータといった共変量を加え、それらが生存にどのように影響するかを評価する枠組みを示す。各患者に対して推定されたリスクスコアに基づき、高リスク群・低リスク群などへの層別化が可能となる。

ここでは、$x\in {\mathrm{ℝ}}^p$を個人に紐づいたp次元の共変量ベクトルとする。これは性別や年齢、治療法、がんのステージなど、時間に依存しない定性的・定量的特徴を想定する。年齢や体重、血圧のように時間とともに変化し得る変数であっても、観察期間中でほぼ一定であるとみなせる場合は、時間非依存の共変量とみなしてxの項目に含めてもよい。Cox回帰では、ある個人の時刻tにおけるハザード関数を、次のように定式化する。

  

$$h(t|x)=h_0(t)\text{・}\exp (x^T\beta )$$（6）

ただし、h（t│x）は共変量xを持つ個人の時刻tにおけるハザード関数、h₀（t）はベースラインハザード関数と呼ばれるもので共変量には依存せず時間tにのみ依存する関数、βは共変量ごとの回帰係数ベクトルであり、最適化による推定対象である。

Cox回帰は比例ハザード仮定（proportional hazards assumption）に基づいている。これは、任意の2つの個人i，jにおいてそのハザード比は時間に依存せず一定であるとする仮定である。

  

$$\frac{h(t|x_i)}{h(t|x_j)}=\exp \left({(x_i-x_j)}^T\beta \right)$$

この比率が時間tの関数とはならない（時間tに依存しない）事が比例ハザードの意味である。これは共変量の効果がどの時刻においても、ベースハザード関数に一定の値が掛け算的（multiplicative）に作用することを意味する。

Cox回帰の大きな利点は、ベースラインハザード関数h₀（t）の具体的な形を仮定することなく、回帰係数βを推定できる点にある。これは、尤度関数を工夫することで、h₀（t）を含まない形に変形できるためである。このような定式化はセミパラメトリックなモデル（semi-parametric model）と呼ばれる。観測データがn個あり、それぞれの個人について観測された生存時間をt_i、イベント発生の有無を示す指標をσ_i∈｛0, 1｝、共変量ベクトルをx_iとする。ここでσ_i＝1のときにイベントが発生しており、σ_i＝0のときは打ち切りである。このとき、Cox回帰では以下のような部分尤度（partial likelihood）を構築する。

  

$$L_{\text{partial}}(\beta )=\prod _{i\mathrm{:}{\sigma }_i=1}\frac{\exp (x_i^T\beta )}{{\sum }_{j\in R(t_i)}\exp (x_j^T\beta )}$$（7）

ここでR（t_i）はリスク集合（時点t_iにおけるまだイベントを経験していない全個人）を表す。この尤度は、イベントが発生した時点において、そのイベントが誰に起こったかの確率を定式化したものである。重要なのは、ベースラインハザードh₀（t）はすべての項に共通なので、分母・分子で打ち消し合い、部分尤度関数に含まれないという点である。したがって、Cox回帰では、h₀（t）を明示的に推定することなくβのみを推定できる。

この部分尤度の対数を取ってマイナスをかけることで、Cox回帰は以下の損失関数（Coxロス）の最小化問題として定式化される。

  

$$l_{\text{partial}}(\beta )=-\sum _{i\mathrm{:}{\sigma }_i=1}\left[x_i^T\beta -\log \left(\sum _{j\in R(t_i)}\exp \left(x_j^T\beta \right)\right)\right]$$（8）

例えば、具体的な例として、3人の個人（全て必ずイベント発生する、打ち切りはなし）において、観測された共変量をそれぞれx₁, x₂, x₃∈ℝ^pとし、イベント発生もこの順番で起きたとする（t₁＜t₂＜t₃）。このときのリスク集合は以下のようになる。

  

$$R(t_1)=\{1,2,3\}$$

  

$$R(t_2)=\{1,2\}$$

  

$$R(t_3)=\{1\}$$

よって、部分尤度は、式（7）より、

  

$$L_{\text{partial}}(\beta )=\frac{\exp (x_1^T\beta )}{\exp (x_1^T\beta )+\exp (x_2^T\beta )+\exp (x_3^T\beta )}\text{・}\frac{\exp (x_2^T\beta )}{\exp (x_2^T\beta )+\exp (x_3^T\beta )}\text{・}\frac{\exp (x_3^T\beta )}{\exp (x_3^T\beta )}$$

最後の項は1になるため、実質的には前半2つのみが尤度に寄与する。同様に、Coxロスは式（8）より、

  

$$\begin{array}{c}{l}_{\text{partial}}(\beta )=-[x_1^T\beta -\log \left(\exp (x_1^T\beta )+\exp (x_2^T\beta )+\exp (x_3^T\beta )\right)+x_2^T\beta \\ -\log \left(\exp (x_2^T\beta )+\exp (x_3^T\beta )\right)]\end{array}$$

となる。この関数は連続かつ滑らかで、勾配が明示的に計算できるため、勾配法やニュートン-ラフソン法などの一般的な最適化手法で最適化可能である。これをβに関して微分すると、

  

$${\nabla }_{\beta }l(\beta )=-\left[x_1-\frac{x_1\exp (x_1^T\beta )+x_2\exp (x_2^T\beta )+x_3\exp (x_3^T\beta )}{\exp (x_1^T\beta )+\exp (x_2^T\beta )+\exp (x_3^T\beta )}+x_2-\frac{x_2\exp (x_2^T\beta )+x_3\exp (x_3^T\beta )}{\exp (x_2^T\beta )+\exp (x_3^T\beta )}\right]$$

というように勾配が得られる。より一般的にn個の観測に対して、Coxロスの勾配は、次のような一般形で記述できる。

  

$${\nabla }_{\beta }l(\beta )=-\sum _{i\mathrm{:}{\sigma }_i=1}\left[x_i-\frac{{\sum }_{j\in R(t_i)}{x}_j\exp (x_j^T\beta )}{{\sum }_{j\in R(t_i)}\exp (x_j^T\beta )}\right]=-\sum _{i\mathrm{:}{\sigma }_i=1}\left[x_i-{𝔼}_{\beta }^{(t_i)}[x]\right]$$（9）

ここで${𝔼}_{\beta }^{(t_i)}[x]$は、時点t_iにおけるリスク集合R（t_i）上での指数加重平均である。言い換えれば、Cox回帰では「イベントが実際に起きた人の共変量が、その時点での平均的なリスク構造と比べてどれだけリスクが高かったか」を積み重ねることで回帰係数βを推定している。

Cox回帰における$x_i^T\beta$は予後指標（Prognostic Index）、またはリスクスコア（Risk Score）などと呼ばれ、主成分スコアのようにp個の共変量に重みを付けて線形結合したスカラー値であり、個人iごとに総合的なリスク指標として割り当てられる。この値が高いほど、イベント発生のハザードが高いと解釈されるため、リスクスコアの中央値を用いて被験者を高リスク群（スコアが中央値より高い群）と低リスク群（スコアが中央値以下の群）に層別化することができる（図4）。試薬投与群とプラセボ群のように事前にわかっているラベルを利用するのではなく、共変量と生存時間との関係をモデル化することで、各個体のリスクを連続的に評価し、その推定リスクに基づいて層別化を行う。以降は、各群に対してKaplan-Meier曲線を描き、ログランク検定によってその差の統計的有意性を検定するという流れが広く用いられている。なお、生存時間データに基づく予測モデルの性能評価には、C統計量（concordance index）やBrierスコアなどの指標が利用される［8, 9, 10, 11, 12, 13, 14, 15, 16, 17］。

Cox回帰は、近年ではオミックス解析においても広く用いられている。特に、2005年に米国で開始された大規模がんゲノムプロジェクトであるThe Cancer Genome Atlas Program（TCGA）では、がん種ごとにRNA-Seqデータと生存・死亡ラベルがセットで公開されており、上記の層別化→生存曲線描画→検定という一連の手順を、さまざまながん種にわたって適用することが可能である［14, 15, 18］。これにより、層別化基準（例：ある遺伝子または遺伝子セットの発現の有無）と、がん患者の生存率との関係性を網羅的に探索できるようになっている。このような解析を簡便に実施するためのパッケージや、Webアプリも複数公開されており、専門的な統計やプログラミングの知識がないユーザーであっても、手軽に生存時間解析を実行できる体制が整いつつある［18］。

Cox回帰×機械学習

生存時間データの打ち切りの扱いの煩雑さによるものだと思われるが、生存時間解析分野はCox回帰への依存が大きい。そのため、1972年に提唱されたCox回帰が、現在もなお多変量な生存時間解析における事実上の標準的手法となっている。しかし、リスクスコアの線形性や比例ハザード仮定は、実データでは必ずしも適切とは言えない。このため、2000年代以降、機械学習との融合による拡張が試みられている。

これらの拡張は、大きく以下の5つに整理できる。第一に、リスクスコアの非線形化である。式（8）に示したCoxロスの構造を保ちつつ、リスクスコアx^Tβをニューラルネットワークなどの非線形関数f_θ（x）に置き換えるアプローチである。Faraggi, Dらの論文［19］で、初めてこのアプローチが考案され、その後Cox-nnet［20］、DeepSurv［21］、SALMON［22］といった深層学習系モデルが相次いで開発され、比例ハザード仮定を維持したまま表現力を高める。

第二に、比例ハザード仮定そのものの緩和である。時間と共変量の相互作用を明示的にモデル化し、ハザード比が時間に依存し得る構造へと拡張する手法である［23, 24］。これはCox回帰の枠組みを拡張する方向に位置付けられる。

第三に、目的関数を変更するアプローチである。部分尤度を用いず、イベントの発生順序を直接学習するランキング型手法（例：Survival SVM［25］）などが含まれる。この場合、Cox回帰には依存しない。

第四に、ツリーやアンサンブルに基づく手法である。Random Survival Forest （RSF）［26］では、決定木の分割基準としてログランク統計量を用いるが、モデル自体はCox回帰には依存しない。

第五に、完全尤度型あるいは離散時間型のモデルである。時間を区間化し、ハザードや累積発生確率を直接推定する深層学習モデル（例：DeepHit［27］）などが該当する。これらは比例ハザード仮定を必要としない。

Cox回帰×行列分解

オミックスデータをCox回帰の共変量として利用する場合、遺伝子数やSNP数が患者数を大きく上回るという高次元性の問題が生じる。このため、そのまま多変量Cox回帰に投入することは困難であり、何らかの変数削減や正則化が必要となる。実際には、

・L1正則化を用いたスパースなCox回帰（Cox-Lasso［28, 29］）

・発現変動遺伝子や疾患と相関するSNPなどに事前に絞り込む方法（比較的緩い有意水準での変数選択［30, 31］）

・各変数ごとの単変量Cox回帰を行い、有意なものを候補とする方法［32］

・全変数でクラスタリングを行い、クラスタ中心を利用することで変数を減らす方法［33, 34, 35］

などの戦略が広く用いられている。これらはいずれも高次元性に対処するための実務的な工夫であり、解析目的やデータ規模に応じて選択される。本稿のテーマである行列分解についても、この流れで生存時間解析との統合が模索されている。たとえば、共変量行列Xが遺伝子発現量行列などオミックスデータであった場合、連動して動くグループがあるため、低ランク性が仮定できる。このような場合に、共変量行列をX≈UV^Tと分解した上で、スコア行列Uを、Cox回帰の共変量として用いる手法が提案されている。例えば、第1回［1］で紹介した主成分分析（Principal Component Analysis; PCA）［36, 37, 38］や非負値行列分解（Non-negative Matrix Factorization; NMF）［39, 40］の利用が挙げられる。

生存時間データへの行列分解の適用においては、共変量行列に含まれる患者の予後予測に関係しないパターンまで取り出されうるため、事前に緩めの有意水準のもとでの単変量Cox回帰を実行し、ある程度患者の予後に関係しそうな遺伝子セットを絞り込んでおいてから次元圧縮する場合もあれば、Cox回帰に基づくリスク情報を教師ラベルとして、より生存に関するような共変量パターンをXから抽出する教師あり次元圧縮（第2回［2］）もある。後者の方向性では、第2回［2］で紹介したPartial Least Squares（PLS）［36, 37, 38, 41, 42, 43, 44, 45, 46］が利用されている。また、より複雑な構造の共変量行列セットから患者の予後に関する共変量パターンを抽出するという方向性で、第2回［2］で紹介したCanonical Correlation Analysis（CCA）［47］や、Collective Matrix Factorization（CMF）［48］が利用されたケースがある。

なお、シングルセルおよび空間オミックスデータでは、患者ごとに多数の細胞やスポットが測定されるため、生存時間解析に直接用いることはできず、何らかの工夫が必要となる。よくある方法としては、細胞あるいはスポット単位で算出された情報を患者単位の要約値へと集約する方法である。代表的なものとして、細胞型の構成比を推定するCell-type Deconvolution［46, 49, 50］や、マーカー遺伝子の平均発現量［51, 52］、特定の遺伝子群の発現を平均や重み付き和としてまとめた遺伝子シグネチャースコア［53, 54］などが用いられる。このような解析においても、行列分解がしばしば利用される。

Cox回帰×テンソル分解

行列分解と同様に、テンソル分解を生存時間解析と統合した例も報告されている。共変量データがマルチオミックス［55］、あるいは患者×表現型×時間［56］、患者×バイオマーカー×時間［57］のような横断データといったケースでは、データ構造のオーダーが行列を超えてテンソルとなる。このような共変量テンソルに対しても、前処理としてテンソル分解を適用し、抽出された因子行列にCox回帰を適用するアプローチが多く用いられる一方で、テンソルの再構築誤差とCoxロスをEnd-to-Endで同時に最適化するモデリング［55］も提案されている。さらに、テンソル共変量に対する回帰係数自体を高階テンソルとして扱い、その係数テンソルに低ランク性を仮定するテンソルCox回帰モデル［57］も報告されている。これは第6回［6］で紹介した、テンソル回帰の一種と言える。

おわりに

今回は、Kaplan-Meier法、ログランク検定、Cox回帰といった、生存時間解析の基礎的な手法を紹介した。生存時間データ特有の問題（打ち切りや時間構造の不均一性）により、従来は機械学習モデルの適用が限定的であったが、近年はCox回帰を基盤とする拡張モデルがいくつか提案されている。行列・テンソル分解も同様に、共変量行列を低次元化することで、高リスク群と低リスク群の分離に寄与するようなパターンを抽出する目的で用いられる。より高次元かつ多モードな共変量データに対しては、テンソル分解を用いた解析手法も報告されている。

次回は、ユーザーが解析的に目的関数の勾配（パラメータに対する偏微分）を導出することなく、ツール側が自動的に勾配を計算することで、複雑なモデリングや最適化を可能にする技術である自動微分と、行列・テンソル分解との関係について議論する。

References

• [1]  Tsuyuzaki, K. (2021/4/23). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part I. Matrix Factorization―. JSBi Bioinformatics Review, 1(2), 18-25. 10.11234/jsbibr.2021.1
• [2]  Tsuyuzaki, K. (2021/10/5). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part II. Simultaneous Matrix Factorization―. JSBi Bioinformatics Review, 2(1), 15-29. 10.11234/jsbibr.2021.6
• [3]  Tsuyuzaki, K. (2022/6/2). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part III. Tensor Factorization―. JSBi Bioinformatics Review, 3(1), 20-33. 10.11234/jsbibr.2022.1
• [4]  Tsuyuzaki, K. (2022/11/1). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part IV. Tensor Factorization―. JSBi Bioinformatics Review, 3(2), 33-46. 10.11234/jsbibr.2022.2
• [5]  Tsuyuzaki, K. (2024/6/3). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part V. Random Walk―. JSBi Bioinformatics Review, 5(1), 1-15. 10.11234/jsbibr.2024.1
• [6]  Tsuyuzaki, K., Matsumoto, H. (2025/6/11). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part VI. Sequential Data―. JSBi Bioinformatics Review, 6(1), 11-28. 10.11234/jsbibr.2025.1
• [7]  Tsuyuzaki, K. (2026/1/6). Mathematics for Heterogeneous Biological Data Fusion Analysis with Matrix-Tensor Factorization―Part VII. Matrix Balancing and Optimal Transport―. JSBi Bioinformatics Review, 6(2), 60-68. 10.11234/jsbibr.2025.5
• [8]  David G. Kleinbaum, Mitchel Klein, 神田英一郎，藤井朋子，古川敏人．（2015/3/27）．エモリー大学クラインバウム教授の生存時間解析、サイエンティスト社
• [9]  杉本知之．（2021/9/7）．生存時間解析（統計解析スタンダード）、朝倉書店
• [10]  辻谷將明、池亀和博．（2024/8/5）．生存時間データ解析：Cox比例ハザードモデルからディープラーニングまで、共立出版
• [11]  朝倉こう子、濱崎俊光．（2015/11）．医学データの統計解析の基本（第5回）生存時間データの解析、日本DDS学会、30(5), 474-484.
• [12]  武冨奈菜美，山本和嬉．（2023/3）．生存時間解析・信頼性解析のための統計モデル，日本統計学会誌，52(2), 69-112. 10.11329/jjssj.52.69
• [13]  河野裕治，青柳陽一郎．（2022/11/18）．生存時間解析―時間とイベント発生の関係―，59(11), 1118-1124
• [14]  Zhao, Z. & Zobolas, J. & Zucknick, M. & Aittokallio, T. (2024/3/4). Tutorial on survival modeling with applications to omics data. Bioinformatics. 40(3), btae132. 10.1093/bioinformatics/btae132
• [15]  Tran, D. & Nguyen, H. & Pham, V.D. & Nguyen, P. & Nguyen Luu, H. & Minh Phan, L. & Blair DeStefano, C. & Jim Yeung, S.C., & Nguyen, T. (2025/3/4). A comprehensive review of cancer survival prediction using multi-omics integration and clinical variables. Briefings in Bioinformatics. 26(2), bbaf150. 10.1093/bib/bbaf150
• [16]  Wang, P. Li, Y., Reddy, C.K. (2019/2/27). Machine Learning for Survival Analysis: A Survey, ACM Computing Surveys, 51(6): 110, 1-36, 10.1145/3214306
• [17]  Pereira, B. de B., Rao, C.R. (2005). Survival Analysis Neural Networks, SBRN, 3(2), 50-60. 10.21528/lmln-vol3-no2-art1
• [18]  カラマンシー．（2025）．TCGAデータを用いた生存時間解析，LabCode. Available from: https://zenn.dev/labcode/books/9d0b4746f5983d
• [19]  Faraggi, D. & Simon, R. (1995). A neural network model for survival data. Stat Med. 14(1), 73-82. 10.1002/sim.4780140108
• [20]  Ching, T. & Zhu, X. & Garmire, L.X. (2018). Cox-nnet: an artificial neural network method for prognosis prediction of high-throughput omics data. PLoS Computational Biology. 14(4), e1006076. 10.1371/journal.pcbi.1006076
• [21]  Katzman, J.L. & Shaham, U. & Cloninger, A. & Bates, J. & Jiang, T. & Kluger, Y. (2018). DeepSurv: personalized treatment recommender system using a Cox proportional hazards deep neural network. BMC Medical Research Methodology. 18, 24. 10.1186/s12874-018-0482-1
• [22]  Huang, Z. & Zhang, X. & Xiang, S. & Johnson, T.S. & Helm, B. & Yu, C.Y. & Zhang, J. & Salama, P. & Rizkalla, M. & Han, Z. & Huang, K. (2019). SALMON: survival analysis learning with multi-omics neural networks on breast cancer. Frontiers in Genetics. 8(10), 166. 10.3389/fgene.2019.00166
• [23]  Kvamme, H. & Borgan, Ø. (2021). Continuous and discrete-time survival prediction with neural networks. Lifetime Data Anal 27, 710-736. 10.1007/s10985-021-09532-6
• [24]  Yuki, S. & Mizukoshi, C. & Abe, K. & Shimamura, T. (2026). scTREND: An annotation-free single-cell time-resolved and condition-dependent hazard model. bioRxiv. 10.64898/2026.01.26.701686
• [25]  Van, Belle, V. & Pelckmans, K. & Van Huffel, S. & Suykens, J.A.K. (2011). Support vector methods for survival analysis: a comparison between ranking and regression approaches. Artif Intell Med. 53(2), 107-118. 10.1016/j.artmed.2011.06.006
• [26]  Ishwaran, H. & Kogalur, U.B. & Blackstone, E.H. & Lauer, M.S. (2008). Random survival forests. Ann Appl Stat. 2(3) 841-860. 10.1214/08-AOAS169.
• [27]  Lee, C. & Zame, W. & Yoon, J. & van der Schaar, M. (2018). DeepHit: A Deep Learning Approach to Survival Analysis With Competing Risks. Proceedings of the AAAI Conference on Artificial Intelligence, 32(1). 10.1609/aaai.v32i1.11842
• [28]  Tibshirani, R. (1997). The lasso method for variable selection in the Cox model. Stat Med. 16, 385-395. 10.1002/(sici)1097-0258(19970228)16:4<385::aid-sim380>3.0.co;2-3.
• [29]  Simon, N. & Friedman, J. & Hastie, T. & Tibshirani, R. (2011). Regularization paths for Cox’s proportional hazards model via coordinate descent. J Stat Softw. 39(5), 1-13. 10.18637/jss.v039.i05
• [30]  Kim, H. & Bredel, M. (2013). Feature selection and survival modeling in The Cancer Genome Atlas. Int J Nanomedicine. 8(Suppl 1), 57-62. 10.2147/IJN.S40733
• [31]  Suresh, K. & Psoter, K.J. (2023). Post hoc survival analyses using RNAseq data: handle with care. Am J Physiol Lung Cell Mol Physiol. 324(1), L1-L4. 10.1152/ajplung.00037.2022
• [32]  Gentles, A.J. & Newman, A.M. & Liu, C.L. & Bratman, S.V. & Feng. W. & Kim, D. & Nair, V.S. & Xu, Y. & Khuong, A. & Hoang, C.D. & Diehn, M. & West, R.B. Plevritis, S.K. & Alizadeh, A.A. (2015). The prognostic landscape of genes and infiltrating immune cells across human cancers. Nature Medicine. 21(8), 938-945. 10.1038/nm.3909
• [33]  Chaudhary, K. & Poirion, O.B. & Lu, L. & Garmire, L.X. (2018). Deep learning-based multi-omics integration robustly predicts survival in liver cancer. Clin Cancer Res. 24(6), 1248-1259. 10.1158/1078-0432.CCR-17-0853
• [34]  Zhang, L. & Lv, C. & Jin, Y. & Cheng, G. & Fu, Y. & Yuan, D. (2018). Deep learning-based multi-omics data integration reveals two prognostic subtypes in high-risk neuroblastoma. Front Genet. 9, 477. 10.3389/fgene.2018.00477
• [35]  Ramazzotti, D. & Lal, A. & Wang, B. & Batzoglou, S. & Sidow, A. (2018). Multi-omic tumor data reveal diversity of molecular mechanisms that correlate with survival. Nat Commun. 9, 4453. 10.1038/s41467-018-06921-8
• [36]  Witten, D.M. & Tibshirani, R. (2010). Survival analysis with high-dimensional covariates. Stat Methods Med Res. 19(1), 29-51. 10.1177/0962280209105024
• [37]  Bair, E. & Tibshirani, R. (2004). Semi-supervised methods to predict patient survival from gene expression data. PLoS Biol. 2(4), E108. 10.1371/journal.pbio.0020108
• [38]  Bair, E. & Hastie, T. & Paul, D. & Tibshirani, R. (2006). Prediction by Supervised Principal Components. Journal of the American Statistical Association, 101(473), 119-137. 10.1198/016214505000000628
• [39]  Lai Y, Hayashida M, Akutsu T. (2013/4/23). Survival analysis by penalized regression and matrix factorization. The Scientific World Journal, 632030. 10.1155/2013/632030
• [40]  Huang Z. (2021). Integrative analysis of multimodal biomedical data with machine learning. Purdue University ProQuest Dissertations & Theses.
• [41]  Yamaguchi, K. & Abdelbaky, S. & Yu, L. & Oakes, C.C. & Abruzzo, L.V. & Coombes, K.R. (2025). PLASMA: Partial LeAst Squares for Multiomics Analysis. Cancers, 17(2), 287. 10.3390/cancers17020287
• [42]  Vahabi, N. & McDonough, C.W. & Desai, A.A. & Cavallari, L.H. & Duarte, J.D. & Michailidis, G. (2021). Cox-sMBPLS: An Algorithm for Disease Survival Prediction and Multi-Omics Module Discovery Incorporating Cis-Regulatory Quantitative Effects. Frontiers in Genetics. 12, 701405. 10.3389/fgene.2021.701405
• [43]  Nygård, S. & Borgan, Ø. & Lingjærde, O.C. & Storvold, H.L. (2008). Partial least squares Cox regression for genome-wide data. Lifetime Data Anal 14, 179-195. 10.1007/s10985-007-9076-7
• [44]  Sadiq, M. & Alnagar, D.K.F. & Abdulrahman, A.T. & Alharbi, R. (2022). The Partial Least Squares Spline Model for Public Health Surveillance Data, Computational and Mathematical Methods in Medicine, 8774742, 7, 10.1155/2022/8774742
• [45]  Nguyen, T.S. & Rojo, J. (2009). Dimension reduction of microarray data in the presence of a censored survival response: a simulation study. Statistical Applications in Genetics and Molecular Biology. 8(1), 4. 10.2202/1544-6115.1395
• [46]  Ge, Q. & Wang, Z. & Wang, Y. & Chu, T. & Wang, Z. & Zhu, W. & CZhang, Y. & Chen, Y. & Ye, D. & Xu, W. & Wang, Z. & Zbudurexiti, M. (2026). Integrative single-cell and spatial transcriptomics with explainable AI reveal lethal prognostic axis in prostate cancer. npj Digital Medicine. 9, 119. 10.1038/s41746-025-02297-4
• [47]  Shao, W. & Cheng, J. & Sun, L. & Han, Z. & Feng, Q. & Zhang, D. & Huang, K. (2018). Ordinal Multi-modal Feature Selection for Survival Analysis of Early-Stage Renal Cancer. MICCAI 2018. Lecture Notes in Computer Science, 11071. 10.1007/978-3-030-00934-2_72
• [48]  Žitnik, M., & Zupan, B. (2014). Survival regression by data fusion. Systems Biomedicine, 2(3), 47-53. 10.1080/21628130.2015.1016702
• [49]  Mizukoshi, C. & Kojima, Y. & Hayashi, S. & Abe, K. & Kasugai, D. & Shimamura, T. (2026). scSurv: a deep generative model for single-cell survival analysis, Bioinformatics, 42(1), btaf671, 10.1093/bioinformatics/btaf671
• [50]  Netskar, H. & Pfefferle, A. & Goodridge, J.P. & Sohlberg, E. & Dufva, O. & Teichmann, S.A. & Brownlie, D. & Michaelsson, J. & Marquardt, N. & Clancy, T. & Horowitz, A. & Malmberg, K.-J. (2024). Pan-cancer profiling of tumor-infiltrating natural killer cells through transcriptional reference mapping. Nature Immunology, 25, 1445-1459. 10.1038/s41590-024-01884-z
• [51]  Zheng, H. & Li, Y. & Zhao, Y. & Jiang, A. (2023). Single-cell and bulk RNA sequencing identifies T cell marker genes score to predict the prognosis of pancreatic ductal adenocarcinoma. Scientific Reports 13(3684). 10.1038/s41598-023-30972-7
• [52]  Xie L. & Kong, H. & Yu, J. & Sun, M. & Lu, S. & Zhang, Y. & Hu, J. & Du, F. & Lian, Q. & Xin, H. & Zhou, J. & Wang, X. & Powell, C. A. & Hirsch, F.R. & Bai, C. & Song, Y. & Yin, J. & Yang, D. (2024). Spatial transcriptomics reveals heterogeneity of histological subtypes between lepidic and acinar lung adenocarcinoma. Clinical and Translational Medicine. 14(2), e1573. 10.1002/ctm2.1573
• [53]  Su, M. & Qiao, K.-Y. & Xie, X.-L. & Zhu, X.-Y. & Gao, F.-L. & Li, C.-J. & Zhao, D.-Q. (2021) Development of a Prognostic Signature Based on Single-Cell RNA Sequencing Data of Immune Cells in Intrahepatic Cholangiocarcinoma. Frontiers in Genetics. 11(615680). 10.3389/fgene.2020.615680
• [54]  Wang, X. & Venet, D. & Lifrange, F. & Larsimont, D. & Rediti, M. & Stenbeck, L. & Dupont, F. & Rouas, G. & Garcia, A.J. & Craciun, L. & Buisseret, L. & Ignatiadis, M. & Carausu, M. & Bhalla, N. & Masarapu, Y. & Villacampa, E.G. & Franzén, L. & Saarenpää, S. & Kvastad, L. & Thrane, K. & Lundeberg, J. & Rothé, F. & Sotiriou, C. (2024). Spatial transcriptomics reveals substantial heterogeneity in triple-negative breast cancer with potential clinical implications. Nature Communications 15(10232). 10.1038/s41467-024-54145-w
• [55]  Zhang, J.Z. & Xu, W. & Hu, P. (2022). Tightly integrated multiomics-based deep tensor survival model for time-to-event prediction. Bioinformatics. 38(12), 3259-3266. 10.1093/bioinformatics/btac286
• [56]  Zhao, J. & Zhang, Y. & Schlueter, D.J. & Wu, P. & Eric Kerchberger, V. & Trent Rosenbloom, S. & Wells, Q.S. & Feng, Q. & Denny, J.C. & Wei, W.Q. (2019). Detecting time-evolving phenotypic topics via tensor factorization on electronic health records: Cardiovascular disease case study. J Biomed Inform. 98, 103270. 10.1016/j.jbi.2019.103270
• [57]  Park, S.H. & Zhou, R. & Zhang, X. & Li, L. & Liu, L. (2024). Tensor landmark analysis with application to ADNI data. Stat Med. 43, 1-15. 10.1002/sta4.70014

略語リスト

・C統計量

Concordance index

・TCGA

The Cancer Genome Atlas Program

・PCA

Principal Component Analysis

・NMF

Non-negative Matrix Factorization

・PLS

Patrial Least Squares （PLS）

・CCA

Canonical Correlation Analysis （CCA）

・CMF

Collective Matrix Factorization （CMF）

著者略歴

露崎 弘毅 2015年、東京理科大学生命創薬科学科博士後期課程終了。博士（薬科学）。同年より、理化学研究所（所属2）に在籍し、シングルセルオミックスのデータ解析や解析ツール開発に従事。現在は千葉大学（所属1）で医療データの解析の傍ら、R, Python, Julia言語のパッケージングを業務として行なっている。パッケージングに特化したハッカソンBio “Pack” athonを主催。 ホームページ：https://researchmap.jp/kokitsuyuzaki