두 표본 집단의 차이가 의미가 있나 알아볼 때 보통
t-검정을 많이 씁니다.
t-검정에서는 값의 분포가 정규분포라는 가정이 있어서, 자료의
분포를 모르거나 정규분포로 바꾸기 매우 난감한 경우에는 t-검정을
할 수 없어서 비모수 검정을 사용하는데요, 이쪽으로 가장 인기있는
검정법은 아무래도 Mann-Whitney U test입니다.
어느날 MW-U로 재미나게~♪ 통계를 하다가~ 아니!! 결과가 엄청 편향되어 나오는 것입니다!
비모수라더니!! 마이크로RNA
마이크로어레이로 나온
결과를 분석하고 있었는데, 마이크로RNA가 염색체에서 떼로 몰려 있어서
한 놈이 올라오면 같이 우루루 올라오는 특성이 있다보니 군집이 엄청 큰
놈들에 의해 전체가 흔들려서, 아 무슨 좋은 방법이 없을까! 하다가
군집을 이룬 자료들에 쓸 수 있게 변형한 것
(Rosner and Grove, 1999,
Datta and Satten, 2005,
Haataja et al., 2009)
들을 발견했습니다. +_+
오…. 다 괜찮아 보이지만, 결국은 코드가 공개돼 있는
Datta와 Satten의 것으 로 해서 일단 결과를 뽑았습니다. 매우 만족스럽군요! ^_^
그런데 문제는 순수 R로 구현되어 있다보니
속도가 엄청나게 느려서, 자료가 별로 크지도 않은데 거의 2시간씩 돌려야
돼서 시험삼아 돌려볼 때도 마음을 크게 먹고 돌려야 해서, 바로바로
결과를 확인하는 재미가 없었습니다.
그래서 쭉 벼르고 있다가, 주말에 여자친구가 기말고사 공부해야 해서, 같이 공부하는
척 하느라 짬이 난 김에
C로 새로 코딩했습니다 (다운로드).
R이나 SciPy같은데서 쉽게 쓸 수 있게 외부 의존성을 없애려고, Z-통계치에서
p-value구하는 부분은 빼서 의존성을 줄였습니다. 그냥 libm만 있으면 됩니다.
(파이썬에서는 scipy.stats.norm.pdf를 써서 Z에서 p-value를 구할 수 있습니다.)
대략 돌려봤더니 2시간 걸리던게 30초로 줄었군요! 이히히 ^__^*
간단하게 컴파일한 다음에 파이썬에서 쓰려면 이렇게 쓸 수 있습니다~
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
|
<span class=“kn”>import</span> <span class=“nn”>ctypes</span>
<span class=“n”>cranksum</span> <span class=“o”>=</span> <span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>cdll</span><span class=“o”>.</span><span class=“n”>LoadLibrary</span><span class=“p”>(</span><span class=“s”>'./cranksum.so'</span><span class=”p”>)</span>
<span class=“n”>crs</span> <span class=“o”>=</span> <span class=“n”>cranksum</span><span class=“o”>.</span><span class=“n”>clustered_rank_sum</span>
<span class=“n”>crs</span><span class=“o”>.</span><span class=“n”>restype</span> <span class=“o”>=</span> <span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>POINTER</span><span class=“p”>(</span><span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>c_double</span> <span class=“o”>*</span> <span class=“mf”>4</span><span class=“p”>)</span>
<span class=“n”>freecrs</span> <span class=“o”>=</span> <span class=“n”>cranksum</span><span class=“o”>.</span><span class=“n”>free_clustered_rank_sum_result</span>
<span class=“k”>def</span> <span class=“nf”>clustered_rank_sum</span><span class=“p”>(</span><span class=“n”>X</span><span class=“p”>,</span> <span class=“n”>grp</span><span class=“p”>,</span> <span class=“n”>cluster</span><span class=“p”>):</span>
<span class=“n”>N</span> <span class=“o”>=</span> <span class=“nb”>len</span><span class=“p”>(</span><span class=“n”>X</span><span class=“p”>)</span>
<span class=“n”>ret</span> <span class=“o”>=</span> <span class=“n”>crs</span><span class=“p”>((</span><span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>c_double</span> <span class=“o”>*</span> <span class=“n”>N</span><span class=“p”>)(</span><span class=“o”>*</span><span class=“n”>X</span><span class=“p”>),</span> <span class=“p”>(</span><span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>c_int</span> <span class=“o”>*</span> <span class=“n”>N</span><span class=“p”>)(</span><span class=“o”>*</span><span class=“n”>grp</span><span class=“p”>),</span>
<span class=“p”>(</span><span class=“n”>ctypes</span><span class=“o”>.</span><span class=“n”>c_int</span> <span class=“o”>*</span> <span class=“n”>N</span><span class=“p”>)(</span><span class=“o”>*</span><span class=“n”>cluster</span><span class=“p”>),</span> <span class=“n”>N</span><span class=“p”>)</span><span class=“o”>.</span><span class=“n”>contents</span>
<span class=“n”>r</span> <span class=“o”>=</span> <span class=“p”>{</span><span class=“s”>'S'</span><span class=”p”>:</span> <span class=”n”>ret</span><span class=”p”>[</span><span class=”mf”>0</span><span class=”p”>],</span> <span class=”s”>'E.S'</span><span class=”p”>:</span> <span class=”n”>ret</span><span class=”p”>[</span><span class=”mf”>1</span><span class=”p”>],</span> <span class=”s”>'Var.S'</span><span class=”p”>:</span> <span class=”n”>ret</span><span class=”p”>[</span><span class=”mf”>2</span><span class=”p”>],</span> <span class=”s”>'z.stat'</span><span class=”p”>:</span> <span class=”n”>ret</span><span class=”p”>[</span><span class=”mf”>3</span><span class=”p”>]}</span>
<span class=“n”>freecrs</span><span class=“p”>(</span><span class=“n”>ret</span><span class=“p”>)</span>
<span class=“k”>return</span> <span class=“n”>r</span>
<span class=“n”>X</span> <span class=“o”>=</span> <span class=“p”>[</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>4</span><span class=“p”>,</span><span class=“mf”>2</span><span class=“p”>,</span><span class=“mf”>4</span><span class=“p”>,</span><span class=“mf”>6</span><span class=“p”>,</span><span class=“mf”>7</span><span class=“p”>,</span><span class=“mf”>4</span><span class=“p”>,</span><span class=“mf”>7</span><span class=“p”>,</span><span class=“mf”>8</span><span class=“p”>]</span>
<span class=“n”>grp</span> <span class=“o”>=</span> <span class=“p”>[</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>]</span>
<span class=“n”>cluster</span> <span class=“o”>=</span> <span class=“p”>[</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>1</span><span class=“p”>,</span><span class=“mf”>2</span><span class=“p”>,</span><span class=“mf”>2</span><span class=“p”>,</span><span class=“mf”>2</span><span class=“p”>,</span><span class=“mf”>2</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>,</span><span class=“mf”>0</span><span class=“p”>]</span>
<span class=“n”>print</span> <span class=“n”>clustered_rank_sum</span><span class=“p”>(</span><span class=“n”>X</span><span class=“p”>,</span> <span class=“n”>grp</span><span class=“p”>,</span> <span class=“n”>cluster</span><span class=“p”>)</span>
|
으음.. R에 붙이는 것은… 아직 잘 몰라서… (나중에..;)
저도 miRNA chip data를 분석하고 있는데, 한가지 궁금한게 있어서요. 위 내용을 보면 일반적인 mRNA expression data는 normalization을 한 후 보통 t-test와 같은 통계기법을 이용하는데, miRNA의 경우는 다른가요? 만약 다르다면 어떤 점을 고려해서 분석해야 하나요? 저는 여태까지 일반적인 microarray data와 같다고 생각했는데 걱정이 되네요. =_=;;
오. 저는 사실 miRNA expression을 중간에 여러 단계로 다시 가공해서 쓰다보니 분산이 완전 제멋대로 나와서 t-test는 도저히 쓸 수 없는 상황이더라구요~ 크흐; 아마도 보통 mRNA array에서 하는 시나리오라면 t-test를 해도 괜찮지 않을까요~ ^_^
참, 그러고보니까 miR-466같은 것은 워낙 클러스터가 커서, expression이 단체로 우루루 쫙 올라갔다가 우루루 쫙 내려오고 그러는 경향이 전체적으로 영향을 많이 줍니다. replicate가 많아서 그냥 probeset단위 t-test를 하실 때는 고려하지 않으셔도 될 것 같은데, gene set analysis같은 것 하실 때는 좀 고려를 해야하지 않을까 싶네요~
그렇군요. 우선 일반적인 방법으로 분석을 한 후 각 gene의 location에 따른 발현변화의 연관성을 조사해 봐야 될 것 같네요. 좋은 정보 감사합니다. ^ ^