플라이어 그룹에 대한 샘플 공분산 행렬 계산

이 예제에서는 http://gettinggeneticsdone.blogspot.com/2009/11/split-apply-and-combine-in-r-using-plyr.html의 샘플 코드를 사용하려고합니다. 그래서, 먼저, 자신의 예제 데이터를 복사 할 수 있습니다 :플라이어 그룹에 대한 샘플 공분산 행렬 계산

mydata=data.frame(X1=rnorm(30), X2=rnorm(30,5,2), 
SNP1=c(rep("AA",10), rep("Aa",10), rep("aa",10)), 
SNP2=c(rep("BB",10), rep("Bb",10), rep("bb",10)))

를 내가이 예에서 SNP2을 무시하고 단지 SNP1의 값은 그룹 구성원을 표시 척하기 위하여려고하고있다. 그렇다면 SNP1의 각 그룹에 대한 요약 통계 ("AA", "Aa", "aa")를 원할 수 있습니다.

은 그럼 각 변수의 수단을 계산하려면, 그것은 감각 (약간 자신의 코드를 수정) 사용한다 :

> ddply(mydata, c("SNP1"), function(df) 
data.frame(meanX1=mean(df$X1), meanX2=mean(df$X2))) 
    SNP1  meanX1 meanX2 
1 aa 0.05178028 4.812302 
2 Aa 0.30586206 4.820739 
3 AA -0.26862500 4.856006

그러나 나는 각 그룹의 샘플 공분산 행렬을 원한다면? 이상적으로는 3D 배열을 원합니다. 여기서 각 그룹에 대한 공분산 행렬을 갖고, 세 번째 차원은 해당 그룹을 나타냅니다. 나는 이전 코드의 수정 된 버전을 시도해보고 잘못된 결과를 내고 있다고 확신하는 다음과 같은 결과를 얻었다.

> daply(mydata, c("SNP1"), function(df) cov(cbind(df$X1, df$X2))) 
, , = 1 


SNP1   1   2 
    aa 1.4961210 -0.9496134 
    Aa 0.8833190 -0.1640711 
    AA 0.9942357 -0.9955837 

, , = 2 


SNP1   1  2 
    aa -0.9496134 2.881515 
    Aa -0.1640711 2.466105 
    AA -0.9955837 4.938320

것은 나는 3 차원의 희미한()가 3이 될 것이라고 생각하지만, 대신, 2. 정말이 각 그룹에 대한 공분산 행렬의 슬라이스 버전입니다. plyr를 사용

  [,1]  [,2] 
[1,] 1.4961210 -0.9496134 
[2,] -0.9496134 2.8815146

이 다음 내가 목록에() 형태로 원하는 걸 제공 : 우리가 수동으로 AA의 샘플 공분산 행렬을 계산하면, 우리가 얻을

> dlply(mydata, c("SNP1"), function(df) cov(cbind(df$X1, df$X2))) 
$aa 
      [,1]  [,2] 
[1,] 1.4961210 -0.9496134 
[2,] -0.9496134 2.8815146 

$Aa 
      [,1]  [,2] 
[1,] 0.8833190 -0.1640711 
[2,] -0.1640711 2.4661046 

$AA 
      [,1]  [,2] 
[1,] 0.9942357 -0.9955837 
[2,] -0.9955837 4.9383196 

attr(,"split_type") 
[1] "data.frame" 
attr(,"split_labels") 
    SNP1 
1 aa 
2 Aa 
3 AA

을하지만 같은 내가 말했듯이, 나는 이것을 3D 배열로 정말 좋아할 것입니다. 내가 daply() 또는 제안에 잘못 갔을 때 어떤 생각을 했습니까? 물론, dlply()에서 3D 배열로 목록을 유형 변환 할 수는 있지만, 시뮬레이션에서이 과정을 여러 번 반복 할 것이므로이 작업을 수행하지 않을 것입니다.

참고로 각 그룹에 대한 샘플 공분산 행렬을 제공하는 한 가지 방법 (http://www.mail-archive.com/[email protected]/msg86328.html)이 있지만 출력 된 개체는 부풀어 오른다.

미리 감사드립니다.

출처

2010-04-28 ramhiser

daply은 배열의 변수 을으로 만듭니다.

a <- daply(mydata, c("SNP1"), function(df) cov(cbind(df$X1, df$X2))) 
l <- dlply(mydata, c("SNP1"), function(df) cov(cbind(df$X1, df$X2)))

이 a[1, , ] 및 l[[1]]는 동일한 출력에 대응하도록한다. wkmor1에서 제안하는대로 aperm을 사용하여 크기를 재정렬 할 수 있지만 초기 양식이 사용자의 요구 사항에 맞지 않는 이유에 대해 자세히 알고 싶습니다.

출처

2010-04-28 12:21:34 hadley

해들리, 나는 위에 표시된 출력에 혼란스러워했습니다. 나는 daply()를 올바르게 사용했다는 것을 몰랐다. 내가 올바르게 한 것은 [1,,] 대신 [,,]을 사용하는 것입니다. 나는 [,, 1]을 사용할 것으로 예상했다. – ramhiser

어떻게 시합 ...

aperm(daply(mydata, c("SNP1"), function(df) cov(cbind(df$X1, df$X2))),perm=c(2,3,1))

'aperm는't '행렬에 그대로 배열하는 것입니다. perm 인자는 희미한 빛의 변화를 지정합니다.

출처

2010-04-28 08:32:36 wkmor1

아, 알겠습니다. 음, 고맙습니다. 건배! – ramhiser

플라이어 그룹에 대한 샘플 공분산 행렬 계산

답변

관련 문제