17

次の形式のデータが与えられた

myDat = structure(list(Score = c(1.84, 2.24, 3.8, 2.3, 3.8, 4.55, 1.13, 
2.49, 3.74, 2.84, 3.3, 4.82, 1.74, 2.89, 3.39, 2.08, 3.99, 4.07, 
1.93, 2.39, 3.63, 2.55, 3.09, 4.76), Subject = c(1L, 1L, 1L, 
2L, 2L, 2L, 3L, 3L, 3L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 6L, 6L, 7L, 
7L, 7L, 8L, 8L, 8L), Condition = c(0L, 0L, 0L, 1L, 1L, 1L, 0L, 
0L, 0L, 1L, 1L, 1L, 0L, 0L, 0L, 1L, 1L, 1L, 0L, 0L, 0L, 1L, 1L, 
1L), Time = c(1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 
1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L)), .Names = c("Score", 
"Subject", "Condition", "Time"), class = "data.frame", row.names = c(NA, 
-24L))

主題、条件、時間の関数としてスコアをモデル化したいと思います。各（人間の）被験者のスコアは、変数Timeで示されるように3回測定されたため、測定を繰り返しました。

サブジェクト効果をランダムに適合させたランダム効果モデルをRに組み込むにはどうすればよいですか？

補遺：これらのデータをどのように生成したかを尋ねられました。ご想像のとおり、1日が長いため、データは偽物です。スコアは時間とランダムノイズの合計であり、条件1にあると、スコアにポイントが追加されます。これは、典型的なPsychのセットアップとして有益です。練習（時間）とスコアを上げる薬（条件== 1）で人々のスコアが良くなるタスクがあります。

この議論の目的のために、いくつかのより現実的なデータがあります。これで、シミュレートされた参加者は、スコアに追加されるランダムな「スキル」レベルを持ちます。また、ファクターは文字列になりました。

myDat = structure(list(Score = c(1.62, 2.18, 2.3, 3.46, 3.85, 4.7, 1.41, 
2.21, 3.32, 2.73, 3.34, 3.27, 2.14, 2.73, 2.74, 3.39, 3.59, 4.01, 
1.81, 1.83, 3.22, 3.64, 3.51, 4.26), Subject = structure(c(1L, 
1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 6L, 
6L, 7L, 7L, 7L, 8L, 8L, 8L), .Label = c("A", "B", "C", "D", "E", 
"F", "G", "H"), class = "factor"), Condition = structure(c(1L, 
1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 
2L, 1L, 1L, 1L, 2L, 2L, 2L), .Label = c("No", "Yes"), class = "factor"), 
    Time = structure(c(1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 
    2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L), .Label = c("1PM", 
    "2PM", "3PM"), class = "factor")), .Names = c("Score", "Subject", 
"Condition", "Time"), class = "data.frame", row.names = c(NA, 
-24L))

それを参照してください：

library(ggplot2)
qplot(Time, Score, data = myDat, geom = "line", group = Subject, colour = factor(Condition))
4

3 に答える 3

7

nlmeライブラリを使用しています...

述べた質問に答えると、次のコードを使用してランダムなインターセプト混合効果モデルを作成できます。

> library(nlme)
> m1 <- lme(Score ~ Condition + Time + Condition*Time,
+ data = myDat, random = ~ 1 | Subject)
> summary(m1)
Linear mixed-effects model fit by REML
 Data: myDat 
       AIC      BIC    logLik
  31.69207 37.66646 -9.846036

Random effects:
 Formula: ~1 | Subject
         (Intercept)  Residual
StdDev: 5.214638e-06 0.3151035

Fixed effects: Score ~ Condition + Time + Condition * Time 
                   Value Std.Error DF  t-value p-value
(Intercept)    0.6208333 0.2406643 14 2.579666  0.0218
Condition      0.7841667 0.3403507  6 2.303996  0.0608
Time           0.9900000 0.1114059 14 8.886423  0.0000
Condition:Time 0.0637500 0.1575517 14 0.404629  0.6919
 Correlation: 
               (Intr) Condtn Time  
Condition      -0.707              
Time           -0.926  0.655       
Condition:Time  0.655 -0.926 -0.707

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-1.5748794 -0.6704147  0.2069426  0.7467785  1.5153752 

Number of Observations: 24
Number of Groups: 8

切片の分散は基本的に0であり、対象内の効果がないことを示しているため、このモデルは時間間の関係をうまく捉えていません。ランダム切片モデルが、反復測定設計に必要なタイプのモデルになることはめったにありません。ランダム切片モデルは、すべての時点間の相関が等しいことを前提としています。つまり、時間1と時間2の相関関係は、時間1と時間3の相関関係と同じです。通常の状況（おそらく偽のデータを生成する状況ではない）では、後者は前者よりも少ないと予想されます。通常、自己回帰構造がより良い方法です。

> m2<-gls(Score ~ Condition + Time + Condition*Time,
+ data = myDat, correlation = corAR1(form = ~ Time | Subject))
> summary(m2)
Generalized least squares fit by REML
  Model: Score ~ Condition + Time + Condition * Time 
  Data: myDat 
       AIC      BIC    logLik
  25.45446 31.42886 -6.727232

Correlation Structure: AR(1)
 Formula: ~Time | Subject 
 Parameter estimate(s):
       Phi 
-0.5957973 

Coefficients:
                   Value Std.Error   t-value p-value
(Intercept)    0.6045402 0.1762743  3.429543  0.0027
Condition      0.8058448 0.2492895  3.232566  0.0042
Time           0.9900000 0.0845312 11.711652  0.0000
Condition:Time 0.0637500 0.1195452  0.533271  0.5997

 Correlation: 
               (Intr) Condtn Time  
Condition      -0.707              
Time           -0.959  0.678       
Condition:Time  0.678 -0.959 -0.707

Standardized residuals:
       Min         Q1        Med         Q3        Max 
-1.6850557 -0.6730898  0.2373639  0.8269703  1.5858942 

Residual standard error: 0.2976964 
Degrees of freedom: 24 total; 20 residual

あなたのデータは、時点間の相関が-.596であることを示していますが、これは奇妙に思えます。通常、少なくとも時点間に正の相関関係があるはずです。このデータはどのように生成されましたか？

補遺：

新しいデータを使用すると、データ生成プロセスがランダム切片モデルと同等であることがわかります（ただし、これは縦断研究では最も現実的ではありません。視覚化により、時間の影響はかなり線形に見えるため、快適に感じるはずです。数値変数として扱います。

> library(nlme)
> m1 <- lme(Score ~ Condition + as.numeric(Time) + Condition*as.numeric(Time),
+ data = myDat, random = ~ 1 | Subject)
> summary(m1)
Linear mixed-effects model fit by REML
 Data: myDat 
       AIC      BIC    logLik
  38.15055 44.12494 -13.07527

Random effects:
 Formula: ~1 | Subject
        (Intercept)  Residual
StdDev:   0.2457355 0.3173421

Fixed effects: Score ~ Condition + as.numeric(Time) + Condition * as.numeric(Time) 
                                  Value Std.Error DF   t-value p-value
(Intercept)                    1.142500 0.2717382 14  4.204415  0.0009
ConditionYes                   1.748333 0.3842958  6  4.549447  0.0039
as.numeric(Time)               0.575000 0.1121974 14  5.124898  0.0002
ConditionYes:as.numeric(Time) -0.197500 0.1586710 14 -1.244714  0.2337
 Correlation: 
                              (Intr) CndtnY as.(T)
ConditionYes                  -0.707              
as.numeric(Time)              -0.826  0.584       
ConditionYes:as.numeric(Time)  0.584 -0.826 -0.707

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-1.44560367 -0.65018585  0.01864079  0.52930925  1.40824838 

Number of Observations: 24
Number of Groups: 8

「yes」条件のスコアが高くなる傾向があることを示す有意な条件効果（約1.7）と、両方のグループが時間の経過とともに上昇することを示す有意な時間効果が見られます。プロットをサポートすると、2つのグループ間で時間の差のある効果（交互作用）は見つかりません。つまり、勾配は同じです。

于 2009-09-04T20:16:35.450 に答える
6

それはあなたの質問に対する答えではありませんが、あなたはあなたのデータのこの視覚化が有益であると思うかもしれません。

library(ggplot2)
qplot(Time, Score, data = myDat, geom = "line", 
  group = Subject, colour = factor(Condition))

データの視覚化

于 2009-09-04T20:54:36.590 に答える
3

（lme4ライブラリを使用）これは、対象の効果をランダムとして、またランダム効果がグループ化される変数に適合します。このモデルでは、ランダム効果は被験者によって変化する切片です。

m <- lmer( Score ~ Condition + Time + (1|Subject), data=myDat )

ランダム効果を確認するには、

ranef(m)

Ian Fellowsが述べたように、データにはランダムな条件と時間のコンポーネントも含まれる場合があります。別のモデルでテストできます。以下の条件では、条件、時間、および切片は、被験者によってランダムに変化することが許可されています。また、それらの相関関係も評価します。

mi <- lmer( Score ~ Condition + Time + (Condition + Time|Subject), data=myDat )

試してみてください

summary(mi)
ranef(mi)

また、切片との相関関係なし、条件と時間の間の相互作用、および他の多数のモデルを使用してこれをテストし、データや理論に最適なモデルを確認することもできます。あなたの質問は少し曖昧ですが、これらのいくつかのコマンドはあなたが始めるのに役立つはずです。

件名はグループ化の要素であるため、他の効果をランダムに当てはめることに注意してください。それは、予測子としても明示的に当てはまるものではありません。

于 2011-08-03T13:40:46.140 に答える