求多重校正的perl腳本

yang7473453

FDR錯誤控制法是Benjamini于1995年提出一種方法,通過控制FDR(FalseDiscovery Rate)來決定P值的域值.假設(shè)你挑選了R個差異表達的基因，其中有S個是真正有差異表達的，另外有V個其實是沒有差異表達的，是假陽性的。實踐中希望錯誤比例Q＝V/R平均而言不能超過某個預(yù)先設(shè)定的值（比如0.05），在統(tǒng)計學(xué)上，這也就等價于控制FDR不能超過5％.
對所有候選基因的p值進行從小到大排序，則若想控制fdr不能超過q，則只需找到最大的正整數(shù)i，使得p(i)<=(i*q)/m.然后，挑選對應(yīng)p(1),p(2),...,p(i)的基因做為差異表達基因，這樣就能從統(tǒng)計學(xué)上保證fdr不超過q。因此，F(xiàn)DR的計算公式如下：
q-value(i)=p(i)*length(p)/rank(p)


表示看完這些我還是不會寫perl，不知道有沒有大神會寫，幫寫一個吧，急求�。。。�！

yinyuemi

回復(fù) 1# yang7473453

大概是這個意思

1. #!/usr/bin/env perl
   
2. my @test_data=(0.0003,0.0001,0.02);
   
3. my %sorted_data = map{$_=>++$i}sort @test_data;
   
4. @p_adjust_data = map{$_ * ($#test_data+1) / $sorted_data{$_}}@test_data;
   
5. print join "\t", "@test_data"';
   
6. print "\n";
   
7. print join "\t", "@p_adjust_data";
   
8. 
   

復(fù)制代碼

rubyish

yinyuemi 發(fā)表于 2014-02-28 17:28
回復(fù) 1# yang7473453

大概是這個意思

這個有意思~學(xué)習(xí)了！

q1208c

表示看完樓主寫的, 我覺得我需要回去重新學(xué)一下語文了.

完全沒理解是什么意思.

pitonas

(⊙０⊙) ~ 比較茫然, 看完以后也確實什么都沒理解


假設(shè)你挑選了R個差異表達的基因，其中有S個是真正有差異表達的，另外有V個其實是沒有差異表達的

大牛, 幫忙指點一下, 舉個栗子:

舉個栗子,什么是 有差異表達 ?
什么是 沒有差異表達 ?

回復(fù) 2# yinyuemi


   

yang7473453

非常感謝，很有幫助回復(fù) 2# yinyuemi


   

yinyuemi

回復(fù) 5# pitonas


    這是統(tǒng)計學(xué)的知識(http://en.wikipedia.org/wiki/False_discovery_rate), 一種控制假陽性率(FDR)的方法, 全稱應(yīng)該是benjamini hochberg.
樓主應(yīng)該是搞生物信息的, 兩種條件下生物的基因的表達量可能相同,也可能不一樣, 前者就是 表達無差異, 后者即是 表達差異