1.緒言

　多施設臨床試験デザインは、適当な期間に十分な症例数を収集するのに有用な方法であるが、参加する施設の規模や特性は様々であり、患者背景や治療の実施方法などが施設間で異なる状況も存在し、ベースラインリスクや治療効果が施設間でばらついてしまうことがある。施設間で治療効果が一様であれば、施設効果を無視した通常の解析方法によって得られる結果の一般化に関して議論ができることから、施設間のばらつきを検討することは重要といえる。

　施設間のばらつきを評価する方法として、施設効果を変量効果とみなした回帰モデルがよく用いられる。生存時間解析の枠組みでは、変量効果を積の形でモデルに取り込んだ比例ハザードモデル(Frailty Model)が既に提案されており、多施設臨床試験データへの適用をはじめ、多くの関連研究が行われている。Frailty Modelにおけるパラメータ推定には、非線形混合効果モデルと同様に、変量効果に関して積分した周辺尤度関数が明示的に表現できないという問題がある。これに対し、EMアルゴリズムなどの推定方法が提案されているが、これらの方法は基準ハザード関数に特定の関数を仮定する必要がある。一方、ラプラス近似尤度などの尤度の近似に基づく方法は、推定の性能が近似の精度に影響を受けるものの、基準ハザード関数に特定の関数を仮定する必要がないという利点がある。

　一般にFrailty Modelを用いた解析では、変量効果が対数正規分布やガンマ分布などの特定の分布に従うことが仮定される。しかし、これらの仮定は制限が強く、必ずしも現実的とはいえない。例えば、治療成績や予後などが極端に良い(又は悪い)施設が存在することにより、真の変量効果の分布が歪んだり裾を引いていたりするとする。このような状況において特定の分布の仮定を用いた解析を行うと、変量効果の分布を誤特定し、治療などの固定効果の推定、変量効果の分布(分散成分など)の推定、そして変量効果の予測に偏りが生じたり精度に影響を及ぼしたりすることが考えられる。非線形混合効果モデルの枠組みでは、変量効果の分布に関する仮定を緩めるために、Seminonparametric Density(SNP分布)を変量効果の分布に仮定することが提案されており、変量効果が2峰性や歪みを持つ状況での有用性が議論されている。SNP分布は、正規分布をはじめ、裾の重い・軽い、右・左への歪み、多峰性のある分布などが表現可能である。

　そこで本論文では、多施設臨床試験から得られた打ち切りを伴う生存時間データに対し、変量効果にSNP分布を仮定した比例ハザードモデルを提案し、ラプラス近似に基づく近似尤度によりパラメータの推定を行う。そして、シミュレーションデータと実際の多施設臨床試験から得られた生存時間データへの適用を通し、提案した方法の性能の評価を行うことを目的とする。

2.Frailty Model

2.1.　比例ハザードモデル

　変量効果を含む比例ハザードモデルにおける、i番目の施設(i=1,…,M)のj番目(j=1,…,ni）の対象者の時点tにおけるハザード関数を以下のように定義する。

但し、λ0(t)は時点tにおける基準ハザード関数、xij＝(xij1,…,xijp)Tは治療を含むp個の共変量からなるベクトル、βkはxijkに対応する固定効果パラメータとする。特に、xi1は治療に関する共変量、β1は治療効果を表すものとする。また、bi=(bi0,bi1)Tは変量効果パラメータからなるベクトルとし、bi0とbi1はそれぞれ施設iのベースラインリスクからの偏り(施設間の患者層の違いなど)、治療効果の全体平均からの偏り(施設と治療効果の交互作用)を示す。但し、施設と治療効果の交互作用がないと想定される場合には、bi0のみのモデルを扱う場合も考えられる。変量効果は通常、対数正規分布やガンマ分布など特定の型の分布に従うことが仮定される。

2.2.　変量効果に対するセミパラメトリックな仮定

　1章で述べたように、施設間差を評価する解析を行う際、変量効果の分布が特定の分布に従っているという証拠は一般的には存在しない。ここでは、変量効果の分布に対する仮定を緩めるために、変量効果は形状や尺度に対してより柔軟なSNP分布に従うと仮定する。SNP分布はデータから形状や尺度が推定されるため、この柔軟な分布を変量効果に仮定することは、セミパラメトリックな方法と呼ぶこともできる。

　まずSNP分布について述べる。q変量のランダムな変数のベクトルuについて、SNP分布は以下のように表される。

但し、α=(α1,…,αq)は|α|=Σqk=1αkをみたす整数αkからなる指標であり、uα=uα11…uαqq、φq(・)はq変量の標準正規分布関数を表す。分布の形状に影響する係数パラメータααは、〓をみたすようにデータから推定される。Kは分布の柔軟性を調整するパラメータであり、あらかじめ定めるものとする。この調整パラメータKは、Kが大きくなるほど、多様な分布が表現できる。先行研究より、K=1又はK=2とした場合でも十分柔軟性のある分布が表現されることが示されている。従って本論文では、調整パラメータKに対し、0、1、2のいずれかをあらかじめ設定することとした。

　次に、SNP分布のFrailty Modelへの応用について述べる。(2)式においてq=2、K=2と設定した場合について具体的に述べる。上記の設定の下で、(2)式は以下のように表される。

但し、α=(α00,…,α20)T、z=(1,u1,…,u22)Tである。ここでB=(Eu[ZZT])1/2を用いて、Bα=(cosθ1,sinθ1cosθ2,sinθ1,sinθ2cosθ3,…,sinθ1sinθ2sinθ3sinθ4sinθ5)と変数変換を行うと、常〓の条件式を満たす。以降、制約のある係数パラメータαの代わりにθ=(θ1,θ2,θ3,θ4,θ5)をデータより推定することとする。更に、R=D1/2(Dは2×2の非負定行列)を用いて、b=Ruと変換すると、(3)式は以下のように表される。

但し、dはDの要素ベクトルである。この(4)式を、(1)式のモデルに含まれる変量効果biが従う確率密度関数として用いることとする。ちなみに(4)式は、K=0とすると、平均(0,0)T、分散Dの2変量正規分布に一致する。

　以上のモデルより、Kを固定した下では、固定効果パラメータβ、変量効果の分布に含まれる分散パラメータd、形状パラメータθがデータから推定される。

3.近似尤度に基づく推定方法

　2章で提案したモデルにおけるパラメータの推定方法について述べる。基準ハザードも含めた周辺尤度は、全ての変量効果パラメータを積分することによって得られる。ここで、基準ハザードに特定の関数形を仮定する必要が生じる。これに対し、本論文では、順位データに基づき基準ハザードに仮定をおかないMarginal Rank Likelihoodを用いることとする。Marginal Rank Likelihoodは、Cox回帰の部分尤度と同様の考え方から導かれ、以下のように表される。

但し、b=(b1T,…,bMT)T、ηij=Σpk=1xijkβk+bi0+xij1bi1、〓、YkI(t)は時点tにおいてリスク集合に含まれるか否かの指示変数である。(5)式は、積分不能のため明示的に表すことができない。そこで、ラプラス近似を応用しMarginal Rank Likelihoodに対する近似尤度を求め、以下の対数尤度関数を得た。

以上より、パラメータ(β,d,θ)の推定値は、(6)の式の対数尤度関数を最大化することによって得られる。

4.シミュレーションによる性能評価

　提案する方法の性能評価を行うために、多施設臨床試験の状況を想定したシミュレーション研究を行った。標本の大きさは、多施設臨床試験の状況をふまえて施設数40、施設内症例数5とし、固定効果(治療効果など)と変量効果が1変量からなる比例ハザードモデルに基づき、データを発生させた。変量効果の分布には、正規性が崩れた場合として、正規分布の混合分布0.8*N(0,0.2)+0.2*2N(0,1)と0.9*N(0,0.2)+0.1*N(0,1.4)、正規性が成立している場合として、正規分布N(0,0.5)とN(0,0.1)を用いた。打ち切り割合は10%、又は30%とした。比較には、提案する方法(K=0,1,2)と正規性を仮定したペナルティー付き尤度に基づく方法(分散推定方法により、ML,REML)を用いた。固定効果の推定の性能評価にはバイアス、平均2乗誤差(MSE)、経験標準誤差、推定標準誤差を用いた。変量効果の分布の推定における性能評価には、調整パラメータの設定間でθとdを直接比較するのが難しいため、分布の期待値、分散、尖度を推定し、比較に用いることとした。

　まず、変量効果の設定に正規分布の混合分布を用いた場合の結果を述べる。固定効果の推定にはバイアス・MSEともに一定の傾向はみられなかった。変量効果の分布に関しては、提案する方法(K=2)が、尖度に関して過小評価の傾向がみられたものの、分散の推定においては他の方法と比べ、真値に近い値を示した。

次に、変量効果の設定に正規分布を用いた場合の結果を述べる。固定効果の推定にはバイアス・MSEともに一定の傾向はみられなかった。提案する方法(K=2)では、分散の推定においては、他の方法と比べ、真値に近い値を示した。しかし、提案する方法(K=2)では尖度に関して過大評価の傾向がみられた。更に、ペナルティー付き尤度に基づく方法は大抵の推定値でバイアスとMSEの観点から良い性能を示していたが、真の変量効果の分布がN(0,0.1)の場合にはラプラス近似による方法と比べ、変量効果の分散を過大評価する傾向がみられた。

5.Z100 Trialデータへの適用

　多施設臨床試験から得られた生存時間データに対し、提案する方法の適用を行った。対象とした試験はZ100の第III相臨床試験(Z100 Trial)である。Z100はがん患者に対して用いられる免疫調整剤である。基礎研究や第II相用量設定試験に基づき、子宮頸がんの放射線療法にZ100の高用量(40μg)を併用したときの有効性を示すため、低用量(0.2μg)を対照としたランダム化試験が行われた。主要評価項目として死亡までの時間が用いられた。本論文では、施設内症例数が1例の施設を除き、計38施設185例(施設内症例数:最小値2,最大値17,中央値4)を解析対象とした。施設を無視したCox回帰の結果では、低用量群の高用量群に対するハザード比は1.41(95%信頼区間:0.94,2.11)であった。

　提案する方法による解析では、固定効果に治療(高用量群/低用量群)と傍結合織浸潤度(高度/軽中等度)、そして施設間差(変量効果)を含めた(1)式と同様の比例バザードモデルを用いた。変量効果は、2変量のSNP分布(調整パラメータ:K=0,1,2)に従うと仮定した。但し、K=0の場合には変量効果に2変量の正規分布を仮定していること同じである。表にZ100 Trialのデータに対し、提案する方法(K=0,1,2)を適用した結果を示す。固定効果パラメータの推定結果は、Kの設定間でほぼ変わらなかった。変量効果の分布に関しては、K=0とK=1の場合でほぼ同様の結果が得られた。K=2では、θの推定値が、正規分布を仮定しているK=0の場合と異なる値を示しており、正規性の仮定から乖離していることが考えられた。次に、図にθとdの推定値を用いて推定された変量効果(b0,b1)の周辺分布を示す(K=1は、K=0と変わらないため省略)。b0の周辺分布はKの設定間でほぼ同様の結果であり、b1ではK=2の方が少しばらつきが小さく推定されているものの、一峰性、対称性を示し、正規分布から大きく乖離してはいなかった。更に、提案するモデル(K=2)の下で、経験ベイズ推定により各施設のベースライン効果と治療効果を予測したところ、ベースライン効果はばらついていたものの、治療効果に関してはそれほど大きなばらつきはみられなかった。

　これらの結果から、Z100 trialのデータでは、変量効果の分布は正規性の仮定からの乖離が小さいことが示唆された。更に、治療と施設に大きな交互作用がみられなかったことから、多少のベースラインリスクにおける施設間差は存在するものの、施設間で治療効果の大きなばらつきはみられず、治療効果が一様であることが示された。

6．考察

　本論文では、多施設臨床試験データに対し施設間のばらつきを評価する方法として、変量効果に柔軟なSNP分布を仮定した比例ハザードモデルを用いることを提案した。この方法により、変量効果に特定の分布型を仮定するという制限を緩めることができた。提案する方法では、変量効果の分布のパラメータ推定値自体の解釈は難しく、推定された変量効果の分布をグラフで視覚的に確認する必要がある。そこで、実際にZ100ランダム化試験データに適用する際、グラフによる確認を行ったところ、変量効果の分布は正規性の仮定からはそれほど乖離していないことが確認された。このグラフを用いた方法は、正規性の前提の確認に有用であると考えられる。

　シミュレーションの結果より、提案する方法は調整パラメータをK=2と設定すると、K=1又はK=0とした場合(K=0の場合は正規分布を仮定していることとなる)に比べ、真の変量効果の分布が正規分布の場合にも、裾を引いた分布(正規性が崩れている)の場合にも、変量効果の分布の分散の推定値が真値に近い値を示した。これらの結果より、セミパラメトリックな分布の仮定がパラメトリックな分布の仮定に比べ、変量効果の分布の推定、特に分散の推定を改善することが示唆された。尖度の推定に伴うバイアスについては、θの推定が初期値の設定に影響をうけていることが考えられる。今回は常にゼロを初期値として用いたが、θのパラメータ数が少なければグリッド探索により初期値を求めるのがよいかもしれない。

　Frailty Modelでは、固定効果パラメータと変量効果に関するパラメータが直交していないため、変量効果の分布の誤特定が固定効果パラメータの推定精度に影響することが考えられる。しかし、シミュレーションによる性能評価の結果からは、変量効果の分散を過大評価、つまり変量効果の分布を誤特定する傾向がある場合にも、固定効果パラメータの推定への影響はみられなかった。先行研究では、治療効果や背景因子などの施設間ではなく施設内でばらつく因子に関しては影響が小さいという報告もあるが、シミュレーションなどによる更なる研究が必要と思われる。

　多施設臨床試験データの他に、クラスター間に相関のある生存時間データの例として、空間データや再発データなどがあげられる。これらのデータでは、クラスター間の相関構造が複雑で、変量効果(クラスター効果)の分布が今回想定したような裾を引いた分布の他に、多峰性や大きな歪みを持つ状況も考えられる。このような状況では、クラスター間の異質性の検討に、特定の分布型を変量効果に仮定したモデルによる解析を用いることは適切でないかもしれない。本論文で提案する方法は、これらの多施設臨床試験データ以外の状況においても応用可能である。

7．結論

　本論文では、多施設臨床試験から得られた打ち切りを伴う生存時間データに対し、変量効果にSNP分布を仮定した比例ハザードモデルを提案し、パラメータの推定にはラプラス近似に基づく近似尤度を用いた。シミュレーションによる性能評価の結果、提案する方法(K=2)が他の正規性を仮定したモデルに比べ、特に変量効果の分布の分散の推定における偏りが小さく、提案する方法の有用性が示唆された。また、Z1OOTrialデータに対して適用した結果、提案する方法(調整パラメータ;K=2)により、変量効果の分布が正規性から大きく乖離していないことが確認された。

図.変量効果の周辺分布の推定結果

(a)b0の周辺分布

(b)b1の周辺分布

表.Z100 Trialデータに対する解析の結果

　本論文は、多施設臨床試験から得られた打ち切りを伴う生存時間データに対し、変量効果に特定の分布型を仮定しない比例ハザードモデルとラプラス近似に基づく近似尤度によりパラメータを推定することを提案したものであり、シミュレーションデータと実際の多施設臨床試験から得られた生存時間データへの適用を通し、提案した方法の性能の評価を行い、下記の結果を得ている。

1.提案する方法の性能評価を行うために、多施設臨床試験の状況を想定したシミュレーション研究を行った。極端な生存時間を示す施設(外れ値)がある状況として、変量効果の分布に、正規分布の混合分布を仮定した。変量効果の分布の推定に関しては、提案する方法(調整パラメータ；K=2)が、尖度に関して過小評価の傾向がみられたものの、分散の推定においては他の方法と比べ、真値に近い値を示した。この結果より、提案する方法は調整パラメータをK=2と設定すると、K=1又はK=0とした場合(K=0の場合は正規分布を仮定)に比べ、裾を引いた分布(正規性が崩れている)の場合には、変量効果の分布の分散の推定値が真値に近い値となることが示された。

2.また、正規性が成立している場合として、変量効果に正規分布N(0,05)とN(0,0.1)を仮定した。その結果、固定効果の推定にはバイアス・平均2乗誤差ともに一定の傾向はみられなかった。提案する方法(K=2)では、分布の分散の推定においては、他の方法と比べ、真値に近い値を示した。従って、提案する方法では、真の変量効果の分布が正規分布の場合にも、変量効果の分布の分散の推定値が真値に近い値となることから、セミパラメトリックな分布の仮定がパラメトリックな分布の仮定に比べ、変量効果の分布の推定、特に分散の推定を改善することが示唆された。

3.シミュレーションによる解析では、固定効果の推定に関しては、正規性が崩れた場合、正規性が成立していない場合のどちらの場合でも、バイアス・平均2乗誤差ともに一定の傾向はみられず、変量効果を誤特定下場合の固定効果の推定の影響は認められなかった。この結果は、先行研究では、治療効果や背景因子などの施設間ではなく施設内でばらつく因子に関しては影響が小さいという報告と矛盾しないものではあるが、理論的には固定効果パラメータと変量効果に関するパラメータが直交していないことが影響することが考えられる。従って、シミュレーションなどによる更なる研究が必要と思われる。

4.実際の多施設臨床試験(Z100の第III相臨床試験:Z100 Trial)から得られた生存時間データに対し、提案する方法の適用を行った。但し、変量効果は、2変量のSNP分布(調整パラメータ:K=0,1,2)に従うと仮定した。固定効果パラメータの推定結果は、Kの設定間でほぼ変わらなかった。変量効果の分布に関しては、ベースラインリスク、治療との交互作用を示す変量効果どちらの周辺分布も、一峰性、対称性を示し、正規分布から大きく乖離してはいなかった。更に、提案するモデル(K=2)の下で、各施設のベースライン効果と治療効果を予測したところ、ベースライン効果はばらついていたものの、治療効果に関してはそれほど大きなばらつきはみられなかった。これらの結果から、Z100 Trialのデータでは、変量効果の分布は正規性の仮定からの乖離が小さいことが示唆された。また、治療と施設に大きな交互作用がみられなかったことから、多少のベースラインリスクにおける施設間差は存在するものの、施設間で治療効果の大きなばらつきはみられず、治療効果が一様であることが示された。

5.多施設臨床試験から得られた打ち切りを伴う生存時間データの解析において、施設間のばらつきを評価する方法として用いられる、変量効果を積の形でモデルに取り込んだ比例ハザードモデルにおいて、変量効果に特定の分布型でなく、データからセミパラメトリックに推定されるSeminonparametric Density(SNP分布)を仮定した。この方法により、変量効果に特定の分布型を仮定するという制限を緩めることができた。この方法は、真の変量効果の分布が歪んだり裾を引いていたりするような状況において有用と考えられる。

　以上、本論文は変量効果を伴う比例ハザードモデルにおいて、変量効果の分布型の仮定をゆるめた統計的方法を提案した初めての研究であり、提案された方法論は、相関を伴う生存時間データを解析するための有用な方法であると考えられ、学位の授与に値するものと考えられる。