SAM算法理解

bioconductor系列教程之一分析基因芯片下(數(shù)據(jù)分析)

SAM算法

SAM算法過程基于數(shù)據(jù)的噪音分析。一般的,信噪比會隨著信號的降低而減小。然而當人們針對每一個基因進行分析時,發(fā)現(xiàn)即使他們的信號在同一水平,每一個基因的表達噪音都不一樣。所以無法籠統(tǒng)地給出一個具體地線性方程來分析全部基因,必須一個一個基因單獨分析。假設有對照實驗I和U兩組,I和U實驗各有平行實驗w組,對于某個指定的基因i,定義相對差異值 (relative difference(d(i)))為基因i在實驗組I和U當中的平均值之差除以I和U的總體標準誤差均值。因為當標準誤差均值過低時,會產(chǎn)生一個較高的d(i)值,所以總體的標準誤差均值增加了一個較正項。具體公式如下:

其中上劃線表示平均值 ,?m和?n表示實驗組I和U的基因i表達的樣本標準差,n1和n2表示實驗組I和U當中基因i參與統(tǒng)計的探針的個數(shù)。S0為較正項,它可以由計算機自動給出,也可以由人手動指定。S0的取值標準是:將所有的d(i)值排序后100等分,每一等分先求得中位數(shù)絕對離差,然后計算這些中位數(shù)絕對離差值的變異系數(shù)。取s0分別等于100等分當中的s(i)的最小值,依次計算。最終s0的值就是在這些計算中,使變異系數(shù)最小的s0。

我們來比較一下它和t-test的異同:

[圖片上傳失敗...(image-92c592-1597045517277)]

比較之下,我們就會發(fā)現(xiàn)兩者十分的相同。左邊是兩組具有不同探針數(shù)和不同標準差的非平行實驗t-test計算公式,右邊就是具有不同探針數(shù)和不同標準差的非平行實驗的SAM公式,我們可以看到,它們的最終差別,只在那個s0附加項上。由此我們可以知道,實際上SAM就是t-test的一個具體應用。

在計算得到所有的d(i)值之后,重排樣品計算d(i)期待值dE(i)。然后用期待值dE(i)對實測值d(i)做圖。如果期待值和實測值相同,那么作圖點就應該落在斜率為1的直線上,如果不同,作圖點就落在其外。差別越大,距離該直線就越遠。給定一個D值(默認的話,可以設置為1.2),當期待值和實測值差別比D大時,就認為有顯著差異。這就是SAM算法。

下面我們就來使用R來實踐一下SAM算法。

<pre style="box-sizing: border-box; font-family: "Roboto Slab", Georgia, Times, serif; font-size: 1em; font-weight: 300; font-style: inherit; margin: -2px 0px 27px; padding: 0px; vertical-align: baseline; border: 0px; outline: 0px; line-height: 27px; overflow: auto; max-width: 100%; background: transparent; color: rgb(102, 102, 102);">##安裝SAM庫

install.packages("samr")

--- Please select a CRAN mirror for use in this session ---

also installing the dependency ‘impute’

trying URL 'http://software.rc.fas.harvard.edu/mirrors/R/bin/windows/contrib/2.11/impute_1.24.0.zip'

Content type 'application/zip' length 1205991 bytes (1.2 Mb)

opened URL

downloaded 1.2 Mb

trying URL 'http://software.rc.fas.harvard.edu/mirrors/R/bin/windows/contrib/2.11/samr_1.28.zip'

Content type 'application/zip' length 80418 bytes (78 Kb)

opened URL

downloaded 78 Kb

package 'impute' successfully unpacked and MD5 sums checked

package 'samr' successfully unpacked and MD5 sums checked

The downloaded packages are in

    C:\Documents and Settings\jianhong ou\Local Settings\Temp\RtmpHrqW1s\downloaded_packages

library(samr)

Loading required package: impute

表達定量

eset <- rma(Data)

Background correcting

Normalizing

Calculating Expression

eset.e <- exprs(eset)

</pre>

SAM分析。我們首先來了解一下samr函數(shù)。

?samr

我們注意到其參數(shù)resp.type有以下幾種類型:"quantitative" 定量反應,其中data參數(shù)當中的y必須是數(shù)值型數(shù)組; "Two class unpaired" 比較兩種不同條件下的結(jié)果,其中data參數(shù)當中的y必須是1或者2數(shù)組,即要么值屬于其中一種條件,要么屬于另外一種條件; "Survival" for censored survival outcome; "Multiclass" : 兩組以上實驗條件,其中data參數(shù)當中的y是自然數(shù)數(shù)組,指定實驗條件編號; "One class" 單組實驗,y是1的數(shù)組; "Two class paired" 比較兩種不同條件下的結(jié)果,并且實驗是成對進行的,其中data參數(shù)當中的y是-1,1,-2,2……-k,k這樣成對出現(xiàn)的數(shù)組; "Two class unpaired timecourse", "One class time course", "Two class.paired timecourse" or "Pattern discovery"與時間相關(guān)的比較,y值必須是kTimet這種形式,其中k是組編號,t是時間,Time就是字符Time,起始的時間點要為kTimetStart格式,結(jié)束的點要為kTimetEnd格式。

Table 1 Resp.type取值范圍

| Resp.type | 取值 |
| Quantitative | 實數(shù),例:27.4 或者-45.34 |
| Two class (unpaired) | 整數(shù)1,2 |
| Multiclass | 整數(shù)1,2,3,…… |
| Paired | 整數(shù)-1,1,-2,2等等成對出現(xiàn),通常認為負號代表未處理的對照組,正號代表實驗組;-1總是與1配對,-2與2配對直至-k與k |
| Survival data | (Time, status)這種成對的數(shù)字,比如(50,1),(120,0)。第一個數(shù)字代表時間,第二個數(shù)字代表狀態(tài)(1=死亡,0=存活) |
| One class | 整數(shù)1 |

在示例的數(shù)據(jù)中,實驗設計為estrogen存在與不存在,及作用10及48小時的比較。因為時間點過少,不足以形成timecourse,所以我們用多組來比較吧。

<pre style="box-sizing: border-box; font-family: "Roboto Slab", Georgia, Times, serif; font-size: 1em; font-weight: 300; font-style: inherit; margin: -2px 0px 27px; padding: 0px; vertical-align: baseline; border: 0px; outline: 0px; line-height: 27px; overflow: auto; max-width: 100%; background: transparent; color: rgb(102, 102, 102);">> gp=rep(1:4,each=2)

gnames <- geneNames(Data)

eset.list <- list(x=eset.e,y=gp,geneid=rownames(eset),

  • genenames=gnames,logged2=TRUE)

samr.results <- samr(eset.list,resp.type="Multiclass")

這里的delta的值會決定多少個基因會被篩選出來,具體可以查看samr.plot

delta=3

delta.table <- samr.compute.delta.table(samr.results)

siggenes.table <- samr.compute.siggenes.table(samr.results,delta,eset.list,delta.table)

siggenes.table$genes.up[,c(2,3,7:10)]

保存結(jié)果:

write.table(siggenes.table$genes.up,file="data",sep="\t")

</pre>

SAM算法因為其提供有 Microsoft Excel插件,所以使用非常廣泛,相反,其在R當中的應用卻顯得相對較少。有研究表明,SAM對FDR的控制并不是很好,這是它最主要的不足。

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點,簡書系信息發(fā)布平臺,僅提供信息存儲服務。

友情鏈接更多精彩內(nèi)容