欠 損 値 があるデータを 扱 うMIプ ロシジャについて 小 野 裕 亮 テクニカルサポートグループ ( 株 )SASインスティチュートジャパン jpnyco@jpn.sas.com Copyright 2000, SAS Institute Inc. All rights reserved. 内 容 : multiple imputation( 多 重 代 入 ; 以 後 MI)を SASで 行 なう 方 法 を 紹 介 プログラミングが 中 心 理 論 実 務 (???) 1) MIを 行 なうためのソフトウェア 必 要 要 件 2) Version6では? 3) MIの 一 般 論 4) PROC MIについて 5) PROC MIANALYZEについて 1
ソフトウェア 必 要 要 件 *MIを 行 なうプロシジャ* *バージョン8.1~ *SAS/STATプロダクト * 評 価 版 (experimental) バージョン6には 存 在 しない 日 本 :バージョン8.1 ( 日 本 語 版 )を Microsoft Windows 版 のみリクエストに 応 じて 出 荷 中 http://www.sas.com/japan/service/v8/index.html 11 月 の 現 状 :マニュアルなどは 英 語, Base SASの 一 部 分 メッセージが 日 本 語 V6における 欠 損 値 (.)の 扱 い ほとんどのプロシジャはリストワイズで 削 除 x1 x2 x3 10.5 19.1 30.9 14.2 17.5 35.5. 20.1 35.4 12.8 20.4 35.2 10.4 15.2 33.8 2
SAS V6における 現 状 A) 単 一 代 入 (single imputation) 何 らかの 値 で 事 前 にとりあえず 穴 を 埋 める 努 力 1)STANDARDプロシジャ(Base SAS) 平 均 Easy! 2)ユーザ 自 身 のマクロプログラム データステップなどによる 涙 ぐましい 努 力 例 )Hot-Deck Imput ation 米 国 SUGI 1998, USA Census B)CORRプロシジャ (Base SAS) デフォルトは ペアワイズで 削 除 C)MIXEDプロシジャ(SAS/STAT) Observed-data likelihood SAS V6における 現 状 D) 打 ち 切 りデータ( 生 存 時 間 モデル トービットモデル) E)IMLによって 自 分 自 身 でプログラミング 例 )ペンシルバニア 大 学 のPaul Allison... 3
余 談 :SASデータセットの 欠 損 値 ( 欠 損 といってもいろいろあるが) SASデータセットの 数 値 変 数 における 欠 損 値 基 本. (ピリオド) 特 殊 欠 損 値.A.B... (ピリオド+アルファベット) 分 析 において この 違 いを 積 極 的 に 使 うものは 特 にない ( 例 外 ) MEANSのCLASSステートメント,FREQの 単 純 集 計 UNIVARIATEの 欠 損 値 数 出 力 TRANSREG,PRINQUALプロシジャ 記 号 :(SASマニュアル) Shafer? Yobs... 観 測 されたデータ Ymis 観 測 されなかったデータ R 欠 損 値 の 位 置 を 表 すインデックス Θ 完 全 データに 対 するモデルパラメータ ξ 欠 損 メカニズムを 決 定 するパラメータ 例 ) X2=β0 + β1 X1 + ε Pr(X1=.)=p1 Pr(X2=.) =p2 モデルパラメータ θ = ( β0, β1, σ) 欠 損 メカニズムのパラメータ ξ=(p1 p2) 4
欠 損 値 いろいろ 1)Missing Completely At Random Yobs, Ymisに 依 存 せず 欠 損 がランダムに 生 じている 2)Missing At Random (ランダムな 欠 損 ) Yobsだけに 依 存 して 欠 損 がランダムに 生 じている Pr(R Yobs, Ymis, ξ) = Pr(R Yobs,ξ) 条 件 や 方 法 によっては 欠 損 が 無 視 できる... 3)nonignorable ( 欠 損 がinformativeな 場 合 ) 欠 損 のメカニズムをモデル 化 する 必 要 Multiple Imputation, Observed-data Likelihoodの の 前 提 : 欠 損 を 無 視 できる? L(θ Yobs)が 最 大 になるようなθ L(θ,ξ Yobs, R) が 最 大 になるようなθ,ξ 1) Missing At Random (or MCAR) 2) θとξがdistinct (ベイズ 流 :θとξが 独 立 ) 分 析 者 は ξには 特 に 関 心 ない 分 析 者 は 完 全 データ 上 のモデルパラメータ θに 関 心 5
喩 え 話 : 2 正 規 変 量 X1 と X2において X2に 欠 損 値 が 存 在 A) MCAR サイコロを 振 って X1,X2の 一 部 分 を 観 測 しない B) MAR (Yobsのみに 依 存 ) X1>1 以 上 の 場 合 には サイコロを 振 ってX2の 値 を 観 測 しない C) Ymisに 依 存 X2 > 1の 場 合 には サイコロを 振 ってX2の 値 を 観 測 し ない 平 均 と 分 散 共 分 散 行 列 の 計 算 研 究 者 の 興 味 は μとσ A) リストワイズの 削 除 でもOK B) リストワイズの 削 除 X2の 平 均 および 分 散, 共 分 散, 相 関 C) リストワイズの 削 除 X2の 平 均 および 分 散, 共 分 散, 相 関 6
補 足 : REG:X2=β0 + β1 X1+ε: 研 究 者 の 興 味 はβ0, β1, σ - A) リストワイズの 削 除 でもOK ( 推 定 効 率, 検 出 力 ) - B) リストワイズの 削 除 でもOK 推 定 効 率, 検 出 力 標 準 化 偏 回 帰 係 数, R2 乗 値 は - C) リストワイズ 削 除 : β0,β1,σの 推 定 値 とも 打 ち 切 りを 考 慮 して 尤 度 関 数 を 定 義 =LIFEREG 補 足 : REG:X1=β0 + β1 X2+ε: 研 究 者 の 興 味 はβ0, β1, σ -A) リストワイズの 削 除 でもOK( 推 定 効 率, 検 出 力 ダウン) -B) リストワイズ 削 除 : β0,β1,σの 推 定 値 とも -C) リストワイズ 削 除 : β0,β1,σの 推 定 値 ともOK( 推 定 効 率, 検 出 力 ダウン 標 準 化 偏 回 帰 係 数 R2 乗 値 は ) 7
V6の の 機 能 :MIXEDプロシジャ Observed-data likelihood を 最 大 化 Observed-data Likelihood を 最 大 化 L(θ,ξ Yobs, R) ではなく L(θ Yobs)を 最 大 化? max( L(θ,ξ Yobs, R)) max ( L(θ Yobs))? 本 当 は max( L(θ Ymis, Yobs))がベスト V6の の 機 能 : MIXEDプロシジャ 多 変 量 分 析 が 目 的 ではないので 8
MI 推 定 を 行 なう 手 順 (2つに 分 類 される) MI 実 行 者 の 任 務 モデルA に 従 う 乱 数 で 穴 埋 めを 実 行 この 処 理 を 複 数 回 実 行 し 穴 埋 めされた 複 数 個 のデータセットを 作 成 分 析 者 の 任 務 各 々の 完 全 化 されたデータセットを モデルB にあ てはめる 穴 埋 めによるバラツキも 考 慮 して 複 数 の 結 果 をまとめる MI 推 定 のSASプログラム 3ステップに 分 けられる MI 実 行 者 1)PROC MIを 用 いる ( 重 要 任 務 :モデルAの 選 択 重 要 任 務 :モデルAの 選 択 ) 分 析 者 2) 既 存 プロシジャでモデルBをあてはめる 3) 結 果 をまとめるためPROC MIANALYZEを 使 用 9
ステップ1: MIプロシジャを 実 行 する ステップ2: 分 析 者 は 既 存 プロシジャを 実 行 する PROC REG DATA=WORK. OUT1 OUTEST=WORK.OUT2 COVOUT ; BY _IMPUTATION_; MODEL X3=X1 X2 ; RUN; 10
ステップ3: 分 析 者 は さらに MIANALYZEプロシジャを 実 行 ステップ1:PROC MI 代 入 モデルを 選 ぶ PROC MIで 用 意 されている 方 法 は 全 部 で3つ MULTINORMAL METHOD= * (1)ロジスティックモデルに 基 づく 方 法 PROPENSITY (propensity score method) (2) 回 帰 モデル 法 REGRESSION Version8.1 (3) 多 変 量 正 規 をMCMCで 発 生 させる 方 法 MCMC ( 制 限 ) 1,2は 欠 損 が 単 調 (monotone)である 時 しか 適 用 できない 11
Monotoneな 欠 損 構 造 T1 T2 T3 T4 X X.. X X X X X... X X X. * 前 が 欠 損 だったら 後 ろも 欠 損 (1)Propensity Score Method METHOD=PROPENSITY (a) 欠 損 値 or 非 欠 損 値 の2 値 を 従 属 変 数 としてロジス ティック (b)ロジスティックモデルの 予 測 値 が 近 い 観 測 値 を 抽 出 し て 欠 損 値 を 穴 埋 め - 欠 損 値 or 非 欠 損 値 の 情 報 しか 用 いていない 付 随 するオプション: METHOD=PROPENSITY( GROUP= n ) - 何 グループにするか? 12
(2) 回 帰 モデルに 基 づく 方 法 METHOD=REGRESSION 回 帰 モデルに 基 づき 穴 埋 めを 行 なっていく 欠 損 構 造 がmonotoneの 時 のみ でも MCMC 法 よりも 高 速 基 本 的 に 多 変 量 正 規 分 布 非 欠 損 値 の 部 分 から βおよびσの 推 定 値 β*, σ*を 乱 数 で 生 成 Y=Xβ* + σ* z で 欠 損 値 部 分 を 埋 める (2) 回 帰 モデルに 基 づく 方 法 METHOD=REGRESSION; 重 要 : 重 要 :Version8.1では METHOD=REGRESSIONは 常 に 間 違 った 結 果 になってい ます Version8.2で このバグ は 修 正 されています 13
(3)MCMC 法 METHOD=MCMC 多 変 量 正 規 分 布 Impute -step と Posterior -stepを 交 互 に 行 なう (Ymis μ,σ,yobs) (μ,σ Y) (Ymis μ,σ,yobs) (μ,σ Y). <オプション> INITIAL=EM( BOOTSTRAP = p MLE (デフォルトは 事 後 分 布 ) CHAIN=SIGLE MULTIPLE 単 鎖 or 複 鎖 PRIOR= 事 前 分 布 (デフォルトは Jeffreys) BITER= 各 連 鎖 のまえに 何 回 回 しておくか? ITER= 単 鎖 のときの 間 隔 PROC MIステートメントのオプション 幾 つのデータ( 完 全 化 されたデータ)を 作 成 するか? NIMPU= k ; 乱 数 系 列 のシード 値 SEED= 代 入 値 の 最 大, 最 小 MAX= MIN= 用 いるデータ 名 DATA= 作 成 する 穴 埋 めデータ 名 OUT= 結 果 を 表 示 しない NOPRINT 14
PROC MIANALYZE 推 定 値 および 推 定 値 間 の 分 散 共 分 散 行 列 を 与 えなけれ ばいけない 入 力 データの 形 式 は2 通 り DATA= _TYPE_= EST & _TYPE_= COV PARMS= 推 定 値 COVB= 分 散 共 分 散 行 列 出 力 される 情 報 : ( 別 紙 ) MI 推 定 値 = 各 推 定 値 の 平 均 欠 損 値 を 埋 めることに 伴 う 変 動 Within-impute Variance (W) Between-impute Variance (B) 欠 損 なしの 時 0 Total Variance = W + B /m + B 欠 損 がMI 推 定 値 の 変 動 にどれだけ 影 響 しているか? - a)fraction of missing information about Q - b)relative increase in variance due to nonresponse 15
今 後 の 予 定 (11/28 現 在 ) バージョン8.2で 変 更 拡 張 を 行 なった マニュアルも V8.1のPROC MIは46ページだが 現 在 V8.2は72ペー ジ バージョン8.2でも MIおよびMIANALYZEは 評 価 版 バージョン8.2の 次 バージョンでは プロダクト 版 となる 予 定 なお カテゴリーデータに 対 する 処 理 は 今 後 の 取 り 組 む べき 課 題 Version8.2に 追 加 される 機 能 (MIプロシジャ) 1. EMステートメント MIを 行 なわずに EMアルゴリズムで(μ,Σ)の 推 定 だ けを 行 ないたい 時 2. TRANSFORMステートメント データを 変 換 3. Monotone MCMC method (monotone missing pattern のデータを 作 成 する) 4. MCMC 法 において 定 常 になったかどうかをチェック するための 自 己 相 関 プロット 16
バージョン8.1のドキュメント ( 英 語 ) - MIプロシジャ (Version8.1) http://www.sas.com/service/techsup/faq/stat_proc/mi proc.html - MIANALYZEプロシジャ (Version8.1) http://www.sas.com/service/techsup/faq/stat_proc/mi analyzeproc.html (PDF 形 式 Adobe 社 が 無 料 提 供 しているAcrobat Readerが 必 要 です) アウトプット 例 : 別 紙 ) 17
18