Jak symulować powtarzane pomiary wyników na wielu odmianach w R?

9

@whuber zademonstrował, jak symulować wyniki na wielu odmianach (y1, y2, i y3) za jeden punkt czasowy.

Jak wiemy, dane podłużne często pojawiają się w badaniach medycznych. Moje pytanie brzmi: w jaki sposób symulować wyniki wielokrotnych wyników wielowymiarowych w R? Na przykład wielokrotnie mierzymyy1, y2, i y3 w 5 różnych punktach czasowych dla dwóch różnych grup leczenia.

Tu.2
źródło

Odpowiedzi:

2

Aby wygenerować normalne dane wielowymiarowe o określonej strukturze korelacji, musisz zbudować macierz kowariancji wariancji i obliczyć jej rozkład Cholesky'ego za pomocą cholfunkcji. Iloczyn rozkładu Choleskiego pożądanej macierzy vcov i niezależnych losowych normalnych wektorów obserwacji da losowe normalne dane z tą wariancją macierzy kowariancji.

v <- matrix(c(2,.3,.3,2), 2)
cv <- chol(v)

o <- replicate(1000, {
  y <- cv %*% matrix(rnorm(100),2)

  v1 <- var(y[1,])
  v2 <- var(y[2,])
  v3 <- cov(y[1,], y[2,])

  return(c(v1,v2,v3))
})

## MCMC means should estimate components of v
rowMeans(o)
AdamO
źródło
2

Użyj funkcji rmvnorm (), wymaga 3 argumentów: macierzy kowariancji wariancji, średnich i liczby wierszy.

Sigma będzie miała 3 * 5 = 15 wierszy i kolumn. Jeden na każdą obserwację każdej zmiennej. Istnieje wiele sposobów ustawiania tych 15 ^ 2 parametrów (ar, dwustronna symetria, nieustrukturyzowany ...). Jednak wypełniając tę ​​macierz, należy pamiętać o założeniach, zwłaszcza gdy ustawisz korelację / kowariancję na zero lub gdy ustawisz dwie wariancje na równe. Na początek macierz sigma może wyglądać mniej więcej tak:

 sigma=matrix(c(
    #y1             y2             y3 
    3 ,.5, 0, 0, 0, 0, 0, 0, 0, 0,.5,.2, 0, 0, 0,
    .5, 3,.5, 0, 0, 0, 0, 0, 0, 0,.2,.5,.2, 0, 0,
    0 ,.5, 3,.5, 0, 0, 0, 0, 0, 0, 0,.2,.5,.2, 0,
    0 , 0,.5, 3,.5, 0, 0, 0, 0, 0, 0, 0,.2,.5,.2,
    0 , 0, 0,.5, 3, 0, 0, 0, 0, 0, 0, 0, 0,.2,.5,
    0 ,0 ,0 ,0 , 0, 3,.5, 0, 0, 0, 0, 0, 0, 0, 0,
    0 ,0 ,0 ,0 ,0 ,.5, 3,.5, 0, 0, 0, 0, 0, 0, 0,
    0 ,0 ,0 ,0 ,0 ,0 ,.5, 3,.5, 0, 0, 0, 0, 0, 0,
    0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3,.5, 0, 0, 0, 0, 0,
    0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3, 0, 0, 0, 0, 0,
    .5,.2,0 ,0 ,0 ,0 ,0 ,0 ,0 , 0, 3,.5, 0, 0, 0,
    .2,.5,.2,0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3,.5, 0, 0,
    0 ,.2,.5,.2,0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3,.5, 0,
    0 ,0 ,.2,.5,.2,0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3,.5,
    0 ,0 ,0 ,.2,.5,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,.5, 3

    ),15,15)

Zatem sigma [1,12] wynosi .2, a to oznacza, że ​​kowariancja między pierwszą obserwacją Y1 a drugą obserwacją Y3 wynosi .2, zależnie od wszystkich pozostałych 13 zmiennych. Nie wszystkie rzędy po przekątnej muszą mieć tę samą liczbę: to moje uproszczone założenie. Czasami ma to sens, a czasem nie. Zasadniczo oznacza to, że korelacja między trzecią obserwacją a czwartą jest taka sama jak korelacja między pierwszą a drugą.

Potrzebujesz także środków. To może być tak proste jak

 meanTreat=c(1:5,51:55,101:105)
 meanControl=c(1,1,1,1,1,50,50,50,50,50,100,100,100,100,100)

Tutaj pierwsze 5 jest średnimi dla 5 obserwacji Y1, ... ostatnich 5 to obserwacje Y3

uzyskaj 2000 obserwacji swoich danych za pomocą:

sampleT=rmvnorm(1000,meanTreat,sigma)
sampleC=rmvnorm(1000,meanControl,sigma)
 sample=data.frame(cbind(sampleT,sampleC) )
  sample$group=c(rep("Treat",1000),rep("Control",1000) )

colnames(sample)=c("Y11","Y12","Y13","Y14","Y15",
                   "Y21","Y22","Y23","Y24","Y25",
                   "Y31","Y32","Y33","Y34","Y35")

Gdzie Y11 jest pierwszą obserwacją Y1, ..., Y15 jest piątą obsadą Y1 ...

Seth
źródło
1
Aby utworzyć macierze jak w pierwszym przykładzie, Seth, spróbuj tego: n <- 3*5; sigma <- diag(1, nrow=n, ncol=n); sigma[rbind(cbind(1:n-1,1:n),cbind(1:n,1:n-1))] <- 1/2. Podobne podejście wygeneruje drugi przykład. Mają jednak wspólny problem: utraciłeś kowariancje wśródyw każdym okresie - te matryce nie odzwierciedlają struktury powtarzanych miar.
whuber
@ whuber twoja składnia jest pomocna, ale różni się od tego, co napisałem. Myślę, że różnica jest trochę ważna. Myślę o tym, co napisałem jako AR (1) i masz wpis w korelacji krzyżowej między ostatnią obserwacją jednej zmiennej a pierwszą obserwacją następnej zmiennej. Innymi słowy, myślę, że sigma [5,6] powinna wynosić 0.
Seth
Ach, teraz widzę, co robisz: w pierwszym przykładzie tworzysz trzy serie AR (1). Przegapiłem to, ponieważ uważam, że OP martwi się również o korelacje między seriami: to właśnie rozumie się przez „wyniki na wielu odmianach”. Z tego punktu widzenia wydaje się, że chcesz zobaczyć te macierze korelacji jako takie5 przez 5 blokuj macierze przy każdym wpisie a 3 przez 3 macierz, a nie jako 3 przez 3 macierz blokowa z 5 przez 5Bloki.
whuber
Myślałem, że mój drugi sigma jest prostym przykładem dopuszczenia dodatniej wariancji między Y1 a Y3. Zredagowałem nieco odpowiedź, aby wyjaśnić, że macierz można skonfigurować w zależności od procesu generowania danych. Jest zdecydowanie wiele sposobów na skórowanie tego kota.
Seth
W porządku, ale twoje podejście stwarza trudności, ponieważ łączenie korelacji wielowymiarowej między yiz modelem AR. Na przykład, czy zdawałeś sobie sprawę, że druga macierz nie jest pozytywnie określona? („Wariancja” c (-102, 177, -204, 177, -102, 0, 0, 0, 0, 0, 102, -177, 204, -177, 102) jest ujemna.)
whuber