数据挖掘：从表达谱芯片原始数据（CEL）到探针注释

原创

生信探索

发布于 2023-02-26 11:15:17

1.9K1

CEL文件：探针的信号值和定位信息，是Affymetrix公司的芯片原始数据。

CEL files contain information on the probe set's intensity values, and a probe set represents a gene. Information about probes gets extracted from the image data by Affymetrix, an image analysis software.

library(affy)
library(limma)
library(stringr)
library(AnnoProbe)
library(magrittr)

0. 下载原始数据

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE66229

mkdir data
cd data
tar -xf GSE66229.TAR
rm GSE66229.TAR

1. 读取原始数据

rawdata <- affy::ReadAffy(celfile.path = "data")

2. rma标准化

rawdata %<>% affy::rma()
exprs <- Biobase::exprs(rawdata)
range(exprs, na.rm = TRUE) # 1.889917 14.620563 不超过50不需要log2转化
# 列重命名
colnames(exprs) <- stringr::str_split(string=colnames(exprs),pattern = "_", simplify = T)[, 1]

3. limma标准化

boxplot(exprs, outline = FALSE, notch = FALSE, las = 2)
#分位数校正
exprs %<>% limma::normalizeBetweenArrays()
boxplot(exprs, outline = FALSE, notch = FALSE, las = 2)
range(exprs, na.rm = TRUE) # 2.09520 14.30741
exprs %<>% as_tibble(rownames = "probe_id")

4. 获取探针注释文件

# 得到探针对应的基因名字
probe2Symbol <- AnnoProbe::idmap("GPL570")
# 看看前几行
head(probe2Symbol,3)
#         probe_id symbol
# 193731   1053_at   RFC2
# 193732    117_at  HSPA6
# 193733    121_at   PAX8

5.ID转换

transid <- function(probe2Symbol, exprs, method = "median") {
    probe2Symbol$probe_id %<>% as.character()
    exprs$probe_id %<>% as.character()
    exprs %>%
        dplyr::inner_join(probe2Symbol, by = "probe_id") %>%
        dplyr::select(-probe_id) %>% 
        dplyr::select(symbol, everything()) %>%
        dplyr::mutate(ref = apply(across(where(is.numeric)), 1, method)) %>%
        dplyr::arrange(desc(ref)) %>%
        dplyr::select(-ref) %>%
        dplyr::distinct(symbol, .keep_all = TRUE)
}

expression <- transid(probe2Symbol, exprs, method = "median")

原创声明：本文系作者授权腾讯云开发者社区发表，未经许可，不得转载。

如有侵权，请联系 cloudcommunity@tencent.com 删除。

analysis