欠測を含む順序カテゴリカル経時データの解析 -GEE プロシジャの有用性 - 駒嵜弘 1 藤原正和 2 ( 1 マルホ株式会社 2 塩野義製薬株式会社 ) Ordinal longitudinal data analysis with missing data -Usefulness of Proc

Size: px

Start display at page:

Download "欠測を含む順序カテゴリカル経時データの解析 -GEE プロシジャの有用性 - 駒嵜弘 1 藤原正和 2 ( 1 マルホ株式会社 2 塩野義製薬株式会社 ) Ordinal longitudinal data analysis with missing data -Usefulness of Proc"

えいじろうだいほうじ
5 years ago
Views:

1 欠測を含む順序カテゴリカル経時データの解析 -GEE プロシジャの有用性 - 駒嵜弘藤原正和 2 ( マルホ株式会社 2 塩野義製薬株式会社 ) Ordinal longitudinal data analysis with missing data -Usefulness of Proc GEE- Hiroshi Komazaki,Masakazu Fujiwara 2 Maruho Co, Ltd., 2 Shionogi & Co., Ltd.

2 2 発表の構成カテゴリカル経時データの発生前発表欠測データの発生多重補完 (REG MCMC FCS) 本発表解析 (CATMOD GEE) 2

3 要旨 : 順序カテゴリカル経時データの解析方法として GEE プロシジャによる方法及び CATMOD プロシジャによる方法をそれぞれ紹介するとともにシミュレーションデータを用いて性能評価を行うキーワード : CATMOD プロシジャ GEE プロシジャ Generalized Estimating Equations (GEE) alternating logistic regression(alr) 3

4 Outline 順序カテゴリカル経時データの解析仕様 ( CATMODプロシジャ GEEプロシジャ ) シミュレーション結果考察 4

5 Background 及び Motivation 順序カテゴリカル経時データに対する解析方法 PROC 特徴 CATMOD 重み付き最小二乗法 ( スライドの式参照 ) GENMOD GEE GEE による推定相関構造は IND しか使用できず ( ただしロバスト分散の使用が可 ) GEE による推定相関構造は EXCH が指定可能 ( スライド 8 4 の式参照 ) SAS 9.4(stat 4.) より PROC GEE は応答変数が順序カテゴリカルデータに対しても実行できるよう DIST=MULTINOMIAL がさらに相関構造も指定できるよう LOGOR= EXCH が追加された PROC CATMOD と PROC GEE を用いてカテゴリカル経時データの解析方法を紹介しさらに性能評価を行う 5

6 順序カテゴリカル経時データにおいて CATMOD プロシジャでの解析方法と GEE プロシジャでの解析方法は Liang and Zeger (986) 及び Prentice(988) の GEE 法に類似 Liang and Zeger (986) 及び Prentice(988) の GEE 法との相違点 CATMOD プロシジャ重み付き最小二乗法分散構造は saturated correlation structure GEE プロシジャ相関パラメータ α の推定方程式を改良 Liang and Zeger (986) 及び Prentice(988) の GEE 法を提示後それぞれのプロシジャでの推定方程式を紹介する 6

7 * 順序カテゴリカルデータを表記する上での留意点本発表にて順序カテゴリカルデータは指示変数 (indicator variable) を用いて表記例 )5 カテゴリの応答変数の場合 Yのカテゴリ指示変数 y ( O c) [,2,3,4 ] = I c ( y ) y 2 y 3 y 4 ( 0 0 0) ( 0 0) ( 0) ( ) ( ) P : 群 i : 被験者 t : 時点 7

8 順序カテゴリカルデータのGEE( Liang and Zeger (986) Prentice(988) ) ( 主効果 β の推定方程式 ) U ( β ) D V ( α) ( Y ) = 0 = = 0 0k β X : カテゴリ k の切片 : デザイン行列 ( 群時点群時点 ) Y : y k D 各群 (=0,) 各症例各時点 (t=,2,3,4) の指示変数例 : 症例 i の指示変数 Y の期待値 β φ y = { y y y y } ) = ln ( k ) = β 2 k k = L( k 0k 5 累積ロジットモデル X β V (α ) 各群 (=0,) の各時点各カテゴリ間の分散共分散行列 8

9 ln β β X = 0 ln ロジットモデル累積ロジットモデルカテゴリ数 :5 の場合 β β X k k k = 0 5 ) ( 2 ln ln ln ln 各カテゴリにて共変量は同じと仮定比例オッズモデル β

10 順序カテゴリカルデータのGEE ( 相関パラメータ α の推定方程式 ) U ( α ) E W ( Z η ( α)) = 0 = = 0 Z η E 相関係数 η α ( ) Z Z の期待値 ( tk )( t k ) α η = ( y k k )( y k k [ ( ) ( )] 2 k k k k ex( f α ) i( tk )( t k ) ( α) = i( tk ) ( t k ) ex( f α ) i( tk ) f ( t k ) ( i( tk t k X と X の関数 ) )( ) : ) Miller(993) W Z i( tk )( t k ) の分散共分散行列 0

11 順序カテゴリカルデータの GEE 推定量 ( ) = = = 0 ) ( 0 ) ( t t t t t t X G X X G X φ β β CATMOD プロシジャの推定量分散共分散構造が saturated correlation structure* のとき CATMOD の推定量は上式の推定量の赤線部分に対応している (Miller,993) *saturated correlation structure 群間時点間カテゴリ間で異なる相関パラメータ α を指定 β ( ) = t t t t V n G φ α φ ) ( β

12 PROC GEEではHeagerty and Zeger(996) の方法を使用 Liang and Zeger (986) Prentice(988) のGEE β の推定方程式 ( スライド 8) 相関パラメータ α の推定方程式 ( スライド0) Heagerty and Zeger(996) の GEE β の推定方程式 ( スライド 8) 相関パラメータ α の推定方程式 ( スライド4) ξ, c2 ) ( y y ) = E i( t, t )( c itc it c 2 2

13 GEE プロシジャの順序カテゴリカルデータの GEE 法は Heagerty and Zeger(996) の手法に基づいている相関パラメータ算出の際 alternating logistic regression(alr) を使用 β 主効果の推定方程式は Liang and Zeger (986) Prentice(988) の GEE と同じ α 相関パラメータは時点間カテゴリ間の条件付き期待値を利用して推定条件付き期待値の回帰パラメータが対数オッズ比となる ( スライド 7 参照 ) ξ, c2 ) ( ) y y = E i( t, t )( c itc it c 2 3

14 ALRによる相関パラメータ α の推定方程式 U ξ = = 0 α t ( ) * α V ( Y ξ ) = 0 y * ( y,, y ) = 4 ( ) n i 4 各群 (=0,) の各時点 (t=,2,3,4) の指示変数各データを ( 全時点 - 当該時点 ) 個分 ( カテゴリ数 -) 個分重複生成している ( 詳細は次ページ参照 ) ξ, c2 ) ( y y ) = E ( t, t )( c tc t c 2 時点 t でカテゴリ c2 が得られた条件のもとで時点 t でカテゴリ c が得られる期待値 (t >t) V i = diag [ ξ ξ )] i ( i 4

15 例 : 説明の簡略化のため時点数 t:3 カテゴリ数 k:3 を想定以下の 2 時点のデータが得られたとする ( 時点 3 のデータは説明変数で使用し応答変数には用いないため省略 ) 時点時点 2 yi y i 2 カテゴリ 2 指示変数 (,0) (,) n i : 被験者 iの時点数 y 各時点各カテゴリごとの条件付き期待値を求めるにあたり応答変数は以下の通りデータを重複発生する必要がある * i n n 2 i ( y,, y, y,, y,, y ) = i ( k ) i ( k ) i2 ( k ) i2 ( k ) i( n ) ( k ) y = = [ yi 2 yi 2 yi2 2 ] [(,,0,0) (,,0,0) (,,,) ] * i i i 5

16 例 : 各時点各カテゴリごとの条件付き期待値の構造は以下の通り ξ = = i [ ξ ξ ξ ] 2 [ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ] 2() 3 2(2) 23 2(2) 従って本例の場合残差 y = 2(22) 3() 3(2) = [ yi 2 yi 2 yi2 2 ] [(,,0,0) (,,0,0) (,,, )] * i ξ = = [ ξ ξ ξ ] 3(2) 3(22) 時点 2 でが得られた条件のもと時点でが得られる期待値 i 2 23() は 2 個のデータで構成される [ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ξ ] 2() 3 2(2) 23 2(2) y * 2(22) ξ 3() 3(2) 3(2) 3(22) 23() 23(2) 23(2) 23(2) 23(2) 23(22) 23(22) 6

17 ξ ( y y ) = E ( t, t )( c tc t c, c2 ) 2 のモデル [ E( Y Y )] log it = itc it c γ i( tt )( c c ) Y 対数オッズ比二項目は切片 it c log µ µ itc itc E Y µ it c ( ) itc Y it c2 E( Y Y ) 2 itc it c 2 γ i( tt )( c c ) = log ρ ( ( ) t i( tt )( cc 2 ) OR Y = = itc, Y it c Z 2 i( tt )( cc 2 ) α log ρi( tt )( c ) c2 2 時点間カテゴリ間の対数オッズ比より相関パラメータ ρ i t Z i = を推定相関係数 ( tt )( c c2 ) が求まる ( のとき相関構造はexchangeable) α 数式の詳細は Heagerty and Zeger(996) Lisitz(99) SAS HELP 参照 7

18 シミュレーションの流れカテゴリカル経時データの発生欠測データの発生多重補完 (REG MCMC FCS) 解析 (CATMOD GEE) 8

19 シミュレーションデータ Donneau(205) 及び Ibrahim and Suliadi(20) を参考にシミュレーションデータを作成 [ Pr( Y k x, t) ] logit = β0 β x β I( t = j) β x I( t j) i y x t ij i k x i t tx i = j : 被験者 : 時点 (,2,3,4) : 応答変数のカテゴリ (,2,3,4) : 応答変数 : 投与群 (0: プラセボ : 実薬 ) : 時点 (,2,3,4) シミュレーション回数 :000 回 k β 0 β 02 β 03 β 04 β x β t β xt パラメータ真値

20 シミュレーションデータ各時点各カテゴリの周辺分布プラセボ群 VISIT K % 実薬群 % VISIT K

21 シミュレーションデータ観測確率モデル logit( ) Y t = 2 0. t 時点前のデータが欠測の有無に影響各群最終時点までに 40% 脱落 2

22 CATMOD プロシジャの指定方法データセットの構造時点 ~ 時点 4 のカテゴリのパターン sim Imutation Number arm 2 _3 _4 count : : : : : : : : 各時点の応答変数パターンを度数で格納度数 0 のパターンは代わりに 0.0 を格納 ( 値の妥当性は未検証 ) 分散共分散行列のランク落ち累積ロジスティックモデルの準完全分離を防ぐため ERROR: Logits cannot be comuted. One of the marginal robabilities is equal to zero. ERROR: The covariance matrix of the resonse functions is singular in oulation. 22

23 CATMOD プロシジャの指定方法 roc catmod data=in_data ; resonse clogits; oulation arm; weight count ; model _*_2*_3*_4 = ( , , , , ~~~~~~~~~~~~~~~~~~~ ~~~~~~~~~~~~~~~~~~~ , , , , , , , , ~~~~~~~~~~~~~~~~~~~ ~~~~~~~~~~~~~~~~~~~ 累積ロジットモデルの指定投与群の指定度数変数の指定列 ~4 : ( 応答変数のカテゴリ数 -) 個の切片列 5 : 投与群列 6~8 : 時点列 9~ : 投与群時点の交互作用 , , , ) ( = 'intercet', 5='arm',6 7 8='VISIT', 9 0 ='arm*visit') ; run; 23

24 GEE プロシジャの指定方法データセットの構造 sim Imutation Number arm ID VISIT AVAL : : : : : : 応答変数の順序カテゴリカルデータ (AVAL:,2,3,4,5) 被験者ごと VISIT ごとに応答変数データを格納 24

25 GEE プロシジャの指定方法 roc gee data=in_data; class arm ID; model AVAL_M =visit arm visit*arm / dist=multinomial; reeated subject=id / logor=exch; run; 相関構造に exchangeable を指定累積ロジットモデルの指定 arm : 投与群 (0 or ) visit : 時点 (,2,3,4) Visit*arm : 投与群時点の交互作用 25

26 結果 : モデルを正しく特定 ver CATMOD 真値推定値 (MSE) BIAS RB αエラー MCMC (0.078) REG (0.090) FCS (0.09) 補完なし (0.092) 完全データ (0.059) GEE 000 ( µ ) 真値推定値 (MSE) BIAS RB αエラー MCMC (0.086) REG (0.095) FCS (0.095) 補完なし (0.094) 完全データ (0.059) i 000 µ µ 000 i µ 26

27 結果 2: 相関構造誤特定 ver CATMOD 真の相関構造真値推定値 (MSE) BIAS RB αエラー MCMC (0.080) REG (0.085) FCS (0.086) 補完なし (0.090) 完全データ (0.058) GEE 真値推定値 (MSE) BIAS RB αエラー MCMC (0.094) REG (0.098) FCS (0.097) 補完なし (0.08) 完全データ (0.065)

28 結果 3: 欠測メカニズムがMNAR ver CATMOD 真値推定値 (MSE) BIAS RB αエラー MCMC (0.20) REG (0.46) FCS (0.47) 補完なし (0.37) GEE 真値推定値 (MSE) BIAS RB αエラー MCMC (0.37) REG (0.58) FCS (0.57) 補完なし (0.58)

29 結果のまとめモデルを正しく特定できたとき性能 : GEE CATMOD GEEでは補完によりバイアス増加相関構造を誤特定したとき性能 : GEE > CATMOD GEEの頑健性を確認できた GEEでは補完によりバイアス減少欠測メカニズムを誤特定したとき性能はGEE> CATMOD 相関構造の誤特定よりはバイアス小 MSEはいずれも増加 GEEでは補完によりバイアス増加 α エラーの増大 (>>2.5) は特に確認できなかった 29

30 発表のまとめカテゴリカル経時データの解析方法の紹介 (CATMOD GEE) GEEプロシジャでは相関パラメータαを対数オッズ比より推定性能 : GEE>CATMOD GEEでは補完した方がバイアスが増大することもあった累積ロジットモデルのとき Weighted-GEEによる解析はSAS のPROC GEEでは実装されていないため今後の導入が期待される 30

31 参考文献 Donneau, A. F., Mauer, M., Molenberghs, G., & Albert, A. (205). A simulation study comaring multile imutation methods for incomlete longitudinal ordinal data. Communications in Statistics-Simulation and Comutation, 44(5), Heagerty, P. J., & Zeger, S. L. (996). Marginal regression models for clustered ordinal measurements. Journal of the American Statistical Association, 9(435), Ibrahim, N. A., & Suliadi, S. (20). Generating correlated discrete ordinal data using R and SAS IML. Comuter methods and rograms in biomedicine, 04(3), e22-e32. Liang, K. Y., & Zeger, S. L. (986). Longitudinal data analysis using generalized linear models. Biometrika, 73(), Lisitz, S. R., Laird, N. M., & Harrington, D. P. (99). Generalized estimating equations for correlated binary data: using the odds ratio as a measure of association. Biometrika, 78(), Miller, M. E., Davis, C. S., & Landis, J. R. (993). The analysis of longitudinal olytomous data: Generalized estimating equations and connections with weighted least squares. Biometrics, Prentice, R. L. (988). Correlated binary regression with covariates secific to each binary observation. Biometrics,

スライド 1

スライド 1 欠測のあるデータにおける population-averaged 及び subject-specific アプローチの性能評価多田圭佑サノフィ株式会社研究開発部門医薬開発本部統計解析プログラミング部統計解析室土居正明駒嵜弘 Performance evaluation of population-averaged and subject-specific approach with missing