目标二：数据预处理，探针ID转换，探针去重
1、数据预处理
自动log2化，解释
ex <- exprSet
qx <- as.numeric(quantile(ex, c(0., 0.25, 0.5, 0.75, 0.99, 1.0), na.rm=T))
LogC <- (qx[5] > 100) ||
(qx[6]-qx[1] > 50 && qx[2] > 0) ||
(qx[2] > 0 && qx[2] < 1 && qx[4] > 1 && qx[4] < 2)
## 开始判断
if (LogC) {
ex[which(ex <= 0)] <- NaN
## 取log2
exprSet <- log2(ex)
print("log2 transform finished")
}else{
print("log2 transform not needed")
}
以上代码在GEO官网可查。
library(limma)
boxplot(exprSet,outline=FALSE, notch=T, las=2) ##出箱线图
### 该函数默认使用quntile 矫正差异
exprSet=normalizeBetweenArrays(exprSet)
boxplot(exprSet,outline=FALSE, notch=T, las=2)
## 这步把矩阵转换为数据框很重要
class(exprSet) ##注释：此时数据的格式是矩阵（Matrix）
exprSet <- as.data.frame(exprSet)
2、探针基因名转换
platformMap.txt中有常见的平台各R注释包的对应关系，是大佬整理好的。直接读取即可。
platformMap <- data.table::fread("/home/data/vip11t50/ZH/GSE42872(2)/platformMap.txt",data.table = F)
## 平台的名称如何知道?（GEO官网也可以查到）
index <- "GPL6244"
## 数据储存在bioc_package这一列中
paste0(platformMap$bioc_package[grep(index,platformMap$gpl)],".db")
## 安装R包，可以直接安装，这里用了判断
if(!requireNamespace("hugene10sttranscriptcluster.db")){
options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/")
BiocManager::install("hugene10sttranscriptcluster.db",update = F,ask = F)
}
##备用方法：在必应上搜索R包名称（如hugene10sttranscriptcluster）加“.db”，可找到下载方法。
## 加载R包
library(hugene10sttranscriptcluster.db)
## 获取探针和基因的对应关系：这是探针注释的关键步骤
probe2symbol_df <- toTable(get("hugene10sttranscriptclusterSYMBOL"))
## 探针有多少个？
length(unique(probe2symbol_df$probe_id)) #结果：19870个。
## 这么多行中，基因名称有重复的么？
length(unique(probe2symbol_df$symbol)) #结果：18859行。
3、探针转换以及去重，获得最终的表达矩阵
tips:
可以用#拆分体会，学会分步探索排查。
符号%>%，是管道操作，其作用是将前一步的结果直接传参给下一步的函数，从而省略了中间的赋值步骤，可以大量减少内存中的对象，节省内存。
library(dplyr)
library(tibble)
exprSet <- exprSet %>%
## 行名转列名,因为只有变成数据框的列,才可以用inner_join
rownames_to_column("probe_id") %>%
## 合并探针的信息
inner_join(probe2symbol_df,by="probe_id") %>%
## 去掉多余信息
dplyr::select(-probe_id) %>% #备注：若跑select时报Error则可加dplyr::试试
## 重新排列
dplyr::select(symbol,everything()) %>%
## rowMeans求出行的平均数(这边的.代表上面传入的数据)
## .[,-1]表示去掉出入数据的第一列，然后求行的平均值
mutate(rowMean =rowMeans(.[,-1])) %>% #备注：.=exprSet
## 把表达量的平均值按从大到小排序
arrange(desc(rowMean)) %>%
## 去重，symbol留下第一个
distinct(symbol,.keep_all = T) %>%
## 反向选择去除rowMean这一列
dplyr::select(-rowMean) %>%
## 列名转行名
column_to_rownames("symbol")
保存数据
save(exprSet,file = "/home/data/vip11t50/ZH/GSE42872(2)/exprSet_rmdup.Rdata"）
可以用load("exprSet_rmdup.Rdata")打开，或左上角File-Open File
btw，把R保存为csv文件可用write.csv(exprSet,file="excel1.csv")
