2017-06-01 6 views
1

コンテキスト含まれるようにぎざぎざにglmer(二項)を翻訳:相関ランダム効果(時間)

を私は個人が0-4からの評価を与えられている12項目のリスク評価(4が最高リスクであることを持っています)。リスクアセスメントは、個人ごとに複数回行うことができます(最大= 19、ほとんどの場合、5回以下の測定値しかありません)。 リスクのベースラインレベルは個人によって異なるので、私はランダムインターセプトモデルを探していますが、リスクの動的性質、すなわち「時間」をランダム係数として追加する必要もあります。

結果はバイナリです:

  • さらに問題測定レベルで行われ、私は基本的に1以上の範囲内で発生したダイナミックどのような変更で探していますことを意味する(FO.bin) 12項目のうち、どのように彼らが測定間の期間にさらなる犯行を犯した個人に貢献したか

私が本質的に探しているのは、个々の人(同じ特性を共有している人)の評価履歴、文脈的な要因、時間の経過とともに変化する要因に基づいて、将来的には劇的なものになるでしょう。

目標:

私は時変(レベル1)を追加することによって、私の 'ベーシック' モデルに追加したいと時不変(レベル2)の予測:

  • 時間には、法令違反、拘禁中の拘禁時間などの刑事司法プロセスに関する疑似変数が含まれています。これらは不変

  • 時間は、初犯の時点で女性、非白人であること、および年齢などの連続予測因子であるとしてダミー変数を含ん評価の間の期間に発生した「イベント」として反映されています 私はこれをlmer4を使ってOKに設定しました。インタラクションや相互作用がある場所など、レベル1とレベル2の予測子を追加すると、興味深い結果が得られます。しかし、拡張されたモデルの複雑さは、収束しないことについての警告メッセージを含むあらゆる種類の警告メッセージを投げている。したがって私は、Rjagsを使用してベイジアンフレームワークに切り替えることが適切であると感じて、自分の発見についてより自信を持って感じることができます。

問題:

は、基本的には、翻訳の一つです。

# the number of Risk Assessments = 552 
N <-nrow(data)                                                             

# number of Individuals (individual previously specified) = 88 
J <- length(unique(Individual))                                            

# the 12 items (previously specified) 
Z <- cbind(item1, item2, item3, item4, ... item12)                         

# number of columns = number of predictors, will increase as model enhanced 
K <- ncol(Z)                                                               

## Store all data needed for the model in a list 

jags.data1 <- list(y = FO.bin, Individual =Individual, time=time, Z=Z, N=N, J=J, K=K)                    

model1 <- function() { 
    for (i in 1:N) { 
    y[i] ~ dbern(p[i]) 
    logit(p[i]) <- a[Individual[i]] + b*time[i] 
  } 
  
  for (j in 1:J) { 
    a[j] ~ dnorm(a.hat[j],tau.a) 
    a.hat[j]<-mu.a + inprod(g[],Z[j,]) 
  } 
  b ~ dnorm(0,.0001) 
  tau.a<-pow(sigma.a,-2) 
  sigma.a ~ dunif(0,100) 
  
  mu.a ~ dnorm (0,.0001) 
  for(k in 1:K) { 
    g[k]~dnorm(0,.0001) 
  } 
} 

write.model(model1, "Model_1.bug") 

を見て:これは私の試みは、バグモデルにこれを翻訳することです

Basic_Model1 <- glmer(BinaryResponse ~ item1 + item2 + item3 + ... + item12 + time + (1+time|individual), data=data, family=binomial) 

:これはlme4使用して、時間とリスク評価における12の項目に基づいて、私の「基本」モデルであります出力は、私の直感では、私は時間のために、様々な係数を追加していませんでしたということで、私がこれまでに行っていること

私のように時間を反映するために、私のバグコードを微調整するにはどうすればよい
Basic_Model2 <- glmer(BinaryResponse ~ item1 + item2 + item3 + ... + item12 + time + (1|individual), data=data, family=binomial) 

の唯一同等です変化する係数、すなわちBasic_Model1?

私が見つけた例に基づいて、私はUループを追加する必要があることを知っているので、U [j]を監視することができます。時間を伴うlogitステートメントは、私は木の木を見ることができないポイントになった!

私よりも多くの専門知識を持つ人が私に正しい方向を向けることを望んでいます。最終的には、レベル1とレベル2の予測変数を追加してモデルを拡張しようとしています。 lme4を使ってこれらを調べたところ、相互作用のクロスレベルのやりとりを指定する必要があることが予想されるため、私はこのように拡張するのに十分柔軟なアプローチを探しています。私は非常にコーディングに新しいので、私と優しいしてください!

答えて

0

そのような場合は、時間の間に自己回帰ガウスモデル(CAR)を使用できます。あなたのタグはwinbugs(またはopenbugs)なので、関数car.normalを以下のように使用することができます。このコードは、あなたのデータセットに適合させる必要があります!

データ

y

は、列のラインと時間で観察とマトリックスであるべきです。 iの時間が同じでない場合は、 NAの値を入力してください。
また、一時的なプロセスのパラメータを定義する必要があります。これは、重みのある近傍の行列です。私は申し訳ありませんが、私は完全にそれを作成する方法を覚えていない...オーダー1の自己回帰のために、このようなものでなければなりません:ぎざぎざとあなたの情報については

jags.data1 <- list(
    # Number of time points 
    sumNumNeigh.tm = 14, 
    # Adjacency but I do not remember how it works 
    adj.tm = c(2, 1, 3, 2, 4, 3, 5, 4, 6, 5, 7, 6, 8, 7), 
    # Number of neighbours to look at for order 1 
    num.tm = c(1, 2, 2, 2, 2, 2, 2, 1), 
    # Matrix of data ind/time 
    y = FO.bin, 
    # You other parameters 
    Individual =Individual, Z=Z, N=N, J=J, K=K)  

モデル

model1 <- function() { 
    for (i in 1:N) { 
     for (t in 1:T) { 
    y[i,t] ~ dbern(p[i,t]) 
    # logit(p[i]) <- a[Individual[i]] + b*time[i] 
    logit(p[i,t]) <- a[Individual[i]] + b*U[t] 
    }} 

    # intrinsic CAR prior on temporal random effects 
    U[1:T] ~ car.normal(adj.tm[], weights.tm[], num.tm[], prec.nu) 
    for(k in 1:sumNumNeigh.tm) {weights.tm[k] <- 1 } 
    # prior on precison of temporal random effects 
    prec.nu ~ dgamma(0.5, 0.0005)  
    # conditional variance of temporal random effects 
    sigma2.nu <- 1/prec.nu         


    for (j in 1:J) { 
    a[j] ~ dnorm(a.hat[j],tau.a) 
    a.hat[j]<-mu.a + inprod(g[],Z[j,]) 
    } 
    b ~ dnorm(0,.0001) 
    tau.a<-pow(sigma.a,-2) 
    sigma.a ~ dunif(0,100) 

    mu.a ~ dnorm (0,.0001) 
    for(k in 1:K) { 
    g[k]~dnorm(0,.0001) 
    } 
} 

dmnormを使用してCARモデルをコーディングする必要があります。

+0

私はどのように乗っているかを報告します。あなたの説明はすべて非常にはっきりしています。プロセスが同時にコーディングを改善するのを助けているので、CARモデルのコーディング方法を調べるのに問題はありません! – Helen

関連する問題