統計とか,研究とか,日常についての備忘録的なもの。この広告は、90日以上更新していないブログに表示しています。潜在変数モデルを統合的に扱うことができる統計分析ソフトであるMplus(https://www.statmodel.com/)はかなり広範な分析を実行可能である。基本的な使い方については,インターネット上に多くの例が示されている。以下の清水先生のサイトなどに日本語の解説がある。norimune.netまた,マニュアルやMuthenらのサイトの掲示版も充実しており,実際的な応用に関しては困る場面は少ないだろう。その一方でMplusには簡単にモンテカルロ・シミュレーションを実行できるということはあまり知られていないように感じられる。自分でモンテカルロ・シミュレーションのプログラムを作成することもできるが,Mplus内でできるモンテカルロ・シミュレーションが実行できれば解析的に導出しなくても機械の力によって近似的な解を求めることができる。これは,例えば複雑なサンプリングプロセスを実行する場合のサンプルサイズの決定や検定力の分析に利用できる。https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3737006/)今回はMplusを用いたモンテカルロ・シミュレーションの方法について,簡単な入門例を示す。(今後もう少し複雑な例を示していきたい。)以下,2因子を仮定した確認的因子分析モデルのモンテカルロ・シミュレーションのコードである。それぞれのコマンドについて簡単に解説する。tech9;ちなみに,Mplusのマニュアル(ver.8の場合)のチャプター12がモンテカルロ・・シミュレーションの解説に当てられている。misspesification状況でのパラメタ推定のバイアスやcoverageについて検証が可能となる。「title:」には,いつもどおり,分析についてのコメントを付記する。この部分は分析の際に無視されるので,メモを書いておくとよい。「montecarlo:」はモンテカルロ・シミュレーションを実行するオプションである。y1-y6;」は分析に利用する変数名を記述する。「-」は多数の変数を指定する際に便利な省略企業なので積極的に利用したい。200;」にはサンプルサイズを示す。100;」には,シミュレーションの繰り返しに回数を示す。1234;」は乱数のシードであり,指定しなければデフォルトは0である。100;」は生成したサンプルをどれだけ保存するかを示す。今回はすべての反復の結果を保存する。fa_monte.dat;」で生成したデータを保存する。hoge*.dat」とするとhoge1.datからhoge100.datというアウトプットを作成できる。モンテカルロシミュレーションの結果を逐一再現したい場合にはこのようなデータを保存しておくとよい。MONTECARLO;」と指定することで先に保存したデータを利用できる。これ以外にも,母集団でのクラス数や階層構造,欠測のパタンなどきめ細やかな状況の指定が可能となっている。population:」は母集団でのモデルの指定を行う。y1-y3*.6;」は因子f1が因子負荷量.6でy1~y3に負荷していることを示している。真値の指定は「*」か「@」によって行う。どちらでもかまわないようだ。基本的モデルの読み方は通常のMplusの分析の状況と変わらない。「model:」コマンドでは,実際に推定に利用するモデルを指定する。ここでは,「*」と「@」が通常の分析のように,初期値と固定パラメタの意味になる。注意が必要なのは,ここで指定した値が,MSEやcoverageの計算に使われる母集団値としてみなされる。MSEの計算などでは,populationコマンドの中のものが利用されるわけではない点に注意が必要である(後のoutpu例を参照)。y1-y3*.6;」は因子f1がy1~y3に負荷していることを示していて,MSは.6からどれくらいの平均自乗誤差があるのかが計算される。「output:」では表示するoutputの内容を指定できる。以下,出力される結果について,簡単に説明する。Dev)が記述されている。PercentilesのExpected部分では理論的なパーセンタイル値が示されており,その値よりも極端な値が出る確率が,Proportionsに示されている。22.561)と計算すると0.9899991という値が得られる。22.561の正規分布を想定していることがわかる。Fitではカイ二乗分布を想定している。)Observedは観測されてた,結果に基づいた値とそれ以上の値を取った割合を示している。次にパラメタ推定値の結果の一部を示す。1.000「Population」の列は,modelコマンドの「*」や「@」で指定した値であることに注意が必要である。Average」は100回の推定結果の平均値である。.5963)/.6を計算すると,バイアスの大きさの割合が計算できる(マニュアルp.472より)。Dev.」は100回の推定値の標準偏差である。Average」は100回の推定のそれぞれで得られた標準誤差の平均値である。各回の推定値)^2をシミュレーション回数で割ったもので,平均自乗誤差である。Cover」は95%信頼区間が100回のシミュレーションのうち,Population列の値を含んでいた割合を示している。(ゆえに,.95に近いほうが,95%信頼区間がその役割を適切に表しているといえる。Coeff」は推定された係数が0かどうかの検定を行った場合,何割が有意になるのかを示している。これがMplusを用いたモンテカルロ・シミュレーションの簡単な実行方法と結果の読み方である。resultsの部分でこれらの結果について情報を得ることができる。GENERATION」の節に記載されている。基本的には,多次元正規分布から連続値を生成し,その後カテゴリカルな変数や打ち切りの変数があればそれに応じて値をカットしたり丸めるなどの処理を行っている。データをRなどで作成し,それを用いてモンテカルロ・シミュレーションを行うこともできるようだ。ANALYSIS」の節で説明がある。その他のコマンドや詳細については気が向いたときに改めて記事にすることして,今日はこのあたりでおしまい。putについて月別アーカイブはてなブログをはじめよう!kazzstatさんは、はてなブログを使っています。ブログを報告するスターをつけました引用をストックしましたストック一覧を見る閉じる引用するにはまずログインしてくださいログイン閉じる引用をストックできませんでした。再度お試しください閉じる限定公開記事のため引用できません。