1.5-OPenMxのACEモデル設定

2696 ワード

出典サイト:http://dbtemp.blogspot.com/2011/08/ace-model-specification-in-openmx.html
ACEモデルを使用して、それらの要因が顕著に表型の変異を引き起こしたのかを検出し、これらの要因の作用値を推定することができる.
飽和モデルでは,共分散行列を直接推定するのではなく,三方晶変異に関係する関数として表現する.この3つの変異はそれぞれA,C,Eであった.また,これらの変異は正である必要があるため,Cholesky分解を用いて標準差を分解した(a 2ではなくaを有効に推定した).
従って、OpenMxでは、これら3つの独立した変異源の行列をmxMatrix命令で指定し、mxAlgebra命令で分散成分を算出する.サンプルコードを次に示します
#---------------------------------------------------------------------------------------------------------
#                                                  ACE Model
#   DVM Bishop 13th March 2010; based on p 18, OpenMx manual
#---------------------------------------------------------------------------------------------------------
mytwindata=read.table("mytwinfile")
#read in previously saved data created with Twin Simulate script
myDataMZ=mytwindata[,1:2]  #columns 1-2 are MZ twin1 and twin2
myDataDZ=mytwindata[,3:4]  #columns 3-4 are DZ twin 1 and twin 2
colnames(myDataMZ)=c("twin1","twin2")
colnames(myDataDZ)=colnames(myDataMZ)
mylabels=colnames(myDataMZ)

#   ACE   twinACE   myACEModel
myACEModel 

ここでモデルのフィットが完了し、出力できます.出力はsummaryですべての統計量を出力してもよいし、mxEvalで特定の量を出力してもよい.一般的に、このモデルは最初の飽和モデルと比較してフィット性(goodness-of-fit)を説明します.推定されたパラメータは標準化できます.次は、データを引き出してフォーマットした出力を示すACEモデルの出力結果スクリプトです.
#---------------------------------------------------------------------------------------------------------
#                                 ACE_Model_Output
#  DVM Bishop, 13th March 2010, based on OpenMxUserGuide, p. 18
#---------------------------------------------------------------------------------------------------------
# NB assumes you have run Saturated twin model (    ), and ACE model, 
# and have likelihood and DF from those models in memory
# -2LL for ACE model from previous script
LL_ACE 

ACEモデルの尤度をネストモデルの尤度と比較するために用いることができ,すなわち,AまたはCを除去して新しいモデルの適合性が劣化するかどうかを観測することを試みることができる.ある項目を削除した後、モデルフィットの優位性が著しく変化しない(変数を含むモデルと削除したモデルのchi 2値を比較することによって測定される)場合、一般的には、設計パラメータを保持するより少ない倹約モデルを削除する傾向がある.次に、モデル内のcを削除し、初期値を0に設定し、'free=FALSE'をcの値を固定する.
変更内容は次のとおりです.
mxMatrix(type="Full",
                   nrow=1,
                   ncol=1,
                   free=FALSE, # fix value of this parameter
                   values=0,    # set start value to zero
                   label="c",
                   name="Y"),


出力結果を得た後,AEモデルとACEモデルを直接相対比し,−2 LL値の差とdfの差を得ることができる.一般に,このchi 2値はCパスの大きさと信頼区間を反映することができ,信頼区間が0を含まない場合,Cを除去するとモデル推定が著しく悪化する.次の文で行うこともできます.
twinAEModel 

第1条新しい変数twinAEModelを生成し、その値はACEモデルと同じである.第2のコマンドは、Y行列を再定義し、その値を0に固定する.