在正式開始操作前做了一些資料的收集，正好搜到了18級學姐的CSDN文章，從她的作業(yè)流程中參考了很多↓

用Bionconductor的affy包處理.cel文件_affy包怎么用

還記得大二學統(tǒng)計的時候，段哥苦口婆心勸我們多看help，這次作業(yè)Help還help了我挺多的

整個作業(yè)的流程我用到的三個包

library(affy)
library(GEOquery)
library(dplyr)

1. 從Gene Expression Omnibus (GEO) 查詢GSE65496，并下載GSE65496_family.soft和原始數(shù)據(jù)(.cel格式)

老老實實查詢的話

紅框是要求的文件的信息，藍色是重要信息

2. 仔細研讀soft格式文件

純文本，\t是分隔符

打開慌神了，啥玩意(是一些其他信息)

好的，到這里才是這次實驗需要用到的信息

3. 試著用Bioconductor的affy包下函數(shù)處理.cel文件 (可以網(wǎng)頁搜索，或者上Bioconductor官網(wǎng)查詢)，得到探針表達值

untar("GSE65496_RAW.tar")
raw <- ReadAffy(celfile.path = "GSE65496_RAW")

關(guān)于ReadAffy()這個函數(shù)，看affy包內(nèi)容的時候其實沒看到，而是read.affybatch()。

看得出來我們就是要一個讀入cel文件的函數(shù)

ReadAffy是read.affybatch的一個包裝器，允許用戶使用部件讀取phenoData、MIAME信息和CEL文件。人們還可以定義文件來讀取PhenoData和MIAME信息。

不用把套娃的cel.tar一個一個解開了，一步到位

4. 下載該數(shù)據(jù)對應(yīng)的平臺數(shù)據(jù)，通過文本處理，將探針和gene symbol對應(yīng)，得到gene symbol的表達值

質(zhì)控

①灰度圖

很簡單的一步，但是耗時挺長。如果是par()成4行2列，可能會報figure margin too large的錯，網(wǎng)上查找了一下相關(guān)的解決方法 ↓

R繪圖： figure margins too large錯誤_figure margins

jpeg("GSE65496_grayscale.jpg", width=200*2, height=200*4)
par(mfcol=c(4,2))
for(i in 1:8){
  image(raw[,i])
}
dev.off()

GSE65496_grayscale.jpg

②RNA降解信息

閱讀包說明文檔

也是一個自帶函數(shù)

par(mfcol=c(1,1))
RNAdeg <- AffyRNAdeg(raw)
plotAffyRNAdeg(RNAdeg,col=topo.colors(8),lwd=2)
legend("topleft",sampleNames(raw),col=topo.colors(8),lwd=2,inset = 0.05,cex=0.3)

GSE65496_degregation.jpg

結(jié)果

背景矯正和歸一化(mas5或rma方法)

兩種方法就是函數(shù)名 mas5() rma()
原理上的區(qū)別的話↓

MAS5算法是由Affymetrix公司設(shè)計開發(fā)的，其對于芯片背景噪音的校正采用了比較直接的策略，即通過計算PM[i]-MM[ii]或PM/MM實現(xiàn)信噪分離。理論上講，MM探針只有非特異雜交信號，不會有特異性雜交，因此MM探針的信號值應(yīng)始終小于其對應(yīng)的PM探針信號值。但是，在實際分析過程中，經(jīng)常發(fā)現(xiàn)甚至多達三分之一的MM探針信號值高于對應(yīng)PM探針信號值。
RMA算法并不直接利用MM的信號去除背景噪音，而是基于20組探針的信號分布采用隨機模型來評估表達值，對于低噪號的實驗有更好的適用性。

我用下來RMA算法是明顯的更快速的，不過兩個結(jié)果都可以放上來看一下(兩種函數(shù)的結(jié)果都是list，不可以直接write.table()，我試過了。用一下expr()再輸出)

msa5

可以看到我苦苦探索的過程??

#write.csv(data.mas5,"GSE65496_mas5.csv")
#csv長度不夠用
#write.table(data.mas5,"GSE65496_mas5.txt",sep='\t')
#輸出結(jié)果不盡如人意，看了看老師的輸出結(jié)果，可能之后要稍微修改一下
#好吧 本來像個list 有專門的提取的
expr.mas5 <- exprs(data.mas5)
write.table(expr.mas5,"GSE65496_mas5.txt",sep='\t')

mas5結(jié)果：綜合考慮pm和mm信號，expression沒有取對數(shù)

rma

data.rma <- rma(raw)
#write.table(data.rma,"GSE65496_rma.txt",sep='\t')
expr.rma <- exprs(data.rma)
write.table(expr.rma,"GSE65496_rma.txt",sep='\t')

rma結(jié)果：RMA只使用pm信號，expression已經(jīng)對數(shù)處理過

探針號和gene symbol對應(yīng)，得到gene表達值

如果前面硬剛不用GEOqueory包還可以，這一步還是要用啊(這里的GPL570就是前面下載的SOFT文件，看我一開始geo頁面的藍框，看到那個platform了沒有？)

GPL570 <- getGEO("GPL570")
names <- Table(GPL570)[c(1,11)]

這里table()函數(shù)是不好使的，但是table輸了一半的時候看到有個Table()函數(shù)，來自GEOqueory。

剩下就是一些非常簡單但是繞來繞去的提取、連接工作了

mas5

expr.mas5 <-as.data.frame(expr.mas5)
expr.mas5$ID <- rownames(expr.mas5)
processed.mas5 <- inner_join(names,expr.mas5,by="ID")
processed.mas5 <- processed.mas5[,-1]
write.table(processed.mas5,"GSE65496_processed_mas5.txt",sep="\t")

和前面沒啥區(qū)別，就是探針名字變基因名了

rma

expr.rma <-as.data.frame(expr.rma)
expr.rma$ID <- rownames(expr.rma)
processed.rma <- inner_join(names,expr.rma,by="ID")
processed.rma <- processed.rma[,-1]
write.table(processed.rma,"GSE65496_processed_rma.txt",sep="\t")

和前面沒啥區(qū)別，就是探針名字變基因名了

總結(jié)

沒啥好總結(jié)的

其他參考過的

解決R語言高清圖片輸出
R語言中的顏色
R語言中Legend 函數(shù)的參數(shù)用法
R語言_legend()函數(shù)用法
玩轉(zhuǎn)數(shù)據(jù)可視化之R語言ggplot2：
R語言函數(shù)-image
R圖形參數(shù)-par() 函數(shù)詳解
cel格式的表達譜芯片數(shù)據(jù)如何讀??？
R-studio中如何利用Bioconductor中affy包處理基因探針數(shù)據(jù)
數(shù)據(jù)挖掘?qū)ｎ}Affymetrix表達譜芯片數(shù)據(jù)預處理
簡單七步教你處理芯片原始數(shù)據(jù)