2016-12-07 3 views
2

내 데이터 집합은 3 개의 시간 지점에 기록 된 유전자 발현 값으로 구성됩니다. tukey 보정과 함께 anova 테스트를 적용하여 시간 지점에서 유전자의 차별적 인 표현을 찾습니다. 각각의 유전자 나 같은 비교 원한다 VS 3 유전자 평가시기 3 대 2 유전자 평가시기 2 대 유전자 평가시기 1 1데이터의 부분 집합 내의 ANOVA

내 데이터는 다음 포맷이다

> head(rf) 
       gene  expn timepoint rep 
2 EG620009 // EG620009 8.428851  x0hr 0 
3     LYPLA1 10.386500  x0hr 0 
21 EG620009 // EG620009 8.582346  x0hr 1 
31    LYPLA1 10.379710  x0hr 1 
22 EG620009 // EG620009 8.566248  x0hr 2 
32    LYPLA1 10.399080  x0hr 2 
    > tail(rf) 
       gene  expn timepoint rep 
23 EG620009 // EG620009 8.561409  x24hr 0 
33    LYPLA1 10.233400  x24hr 0 
24 EG620009 // EG620009 8.750639  x24hr 1 
34    LYPLA1 10.023780  x24hr 1 
25 EG620009 // EG620009 8.560267  x24hr 2 
35    LYPLA1 10.025980  x24hr 2 

내가 할 것 인 경우 :

TukeyHSD(aov(rf$expn ~ rf$timepoint * rf$gene)) 

이 나에게 유전자 즉 모든에서 모든 평가시기 사이의 비교를 줄 것이다. 유전자와 같은 비교를 포함하여 1 시점 대 유전자 b 시점 2

나는 유전자에 의해 부분 집합으로 데이터에 aov 함수를 적용하는 방법을 연구하려고 노력 해왔다. 나는 p 값을 출력으로 제공하고 아래와 같이 by 함수를 사용하여 각 유전자에 개별적으로 적용하려고하는 함수를 정의했습니다.

> gene.aov = function(x) {TukeyHSD(aov(expn ~ timepoint, data = x))} 
> aov.pval = function(y) {y$timepoint[,4]} 
> gene.pval = function(z) {aov.pval(gene.aov(z))} 
> pvals = by(rf$expn,list(rf$gene),gene.pval) 
> Error in eval(predvars, data, env) : 
    numeric 'envir' arg not of length one 

이것이 작동하지 않는 이유는 무엇입니까? 아니면 완전히 다른 방식으로이 문제에 접근해야합니까? 감사합니다. by 그것이 data.frame 또는 매트릭스, 당신이 전달하는 일에 첫 번째 인수의 기대 때문에

답변

1

가 작동하지 않는 것은 numeric 벡터이다 rf$exp입니다. 당신은 이것을 할 수 있었고 잘 동작 할 것입니다 (저는 가독성을 위해 여러 함수를 포기했습니다).

by(rf, rf$gene, function(x) {TukeyHSD(aov(expn ~ timepoint, data = x))}, simplify = F) 

rf$gene: EG620009 
    Tukey multiple comparisons of means 
    95% family-wise confidence level 

Fit: aov(formula = expn ~ timepoint, data = x) 

$timepoint 
       diff  lwr  upr  p adj 
x24hr-x0hr 0.09829 -0.123391 0.319971 0.2857424 

--------------------------------------------------------------------------- 
rf$gene: LYPLA1 
    Tukey multiple comparisons of means 
    95% family-wise confidence level 

Fit: aov(formula = expn ~ timepoint, data = x) 

$timepoint 
       diff  lwr  upr  p adj 
x24hr-x0hr -0.2940433 -0.4876756 -0.100411 0.0135193 
+0

대단히 감사합니다. – Sarah