结肠腺癌细胞系过表达apoM的芯片数据分析

生信技能树

发布于 2022-07-26 10:27:19

6701

发布于 2022-07-26 10:27:19

比如GSE162325这个数据集，它比较新：https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE162325，所以如果你使用我的AnnoProbe包里面的 geoChina("GSE162325") 函数会失败，因为我最近没有空去同步这些新的表达量芯片数据集。

而且平台好奇怪：GPL23126 [Clariom_D_Human] Affymetrix Human Clariom D Assay [transcript (gene) version] ，不过是很标准的两个分组：

GSM4949747 NC-1
GSM4949748 NC-2
GSM4949749 NC-3
GSM4949750 Overexpression-1
GSM4949751 Overexpression-2
GSM4949752 Overexpression-3

因为比较新，所以没有被 AnnoProbe 收录，我是自己去其页面下载了 GSE162325_series_matrix.txt.gz 文件，网页地址是：

https://ftp.ncbi.nlm.nih.gov/geo/series/GSE162nnn/GSE162325/matrix/GSE162325_series_matrix.txt.gz

然后使用 GEOquery 加载它，如下所示的代码：

rm(list = ls())  ## 魔幻操作，一键清空~
options(stringsAsFactors = F)

getOption('timeout')
options(timeout=10000)

library(AnnoProbe)
library(GEOquery) 
if(F){
  gset <- geoChina("GSE162325")
  gset
}

gset <- getGEO("GSE162325",
               filename = 'GSE162325_series_matrix.txt.gz',
               AnnotGPL = FALSE, getGPL = FALSE)
               
# trying URL 'https://ftp.ncbi.nlm.nih.gov/geo/series/GSE162nnn/GSE162325/matrix/GSE162325_series_matrix.txt.gz'
#  cannot open URL 'https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?targ=self&acc=GPL23126&form=text&view=full'

可以看到，确实比较麻烦，而且表达量芯片对应的探针信息也需要自己注释：https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GPL23126

可以看到：

探针信息也需要自己注释

在 https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?targ=self&acc=GPL23126&form=text&view=full 可以下载一个 GPL23126.txt的文本文件，大约是1G，内容有点多，所以需要简单的处理一下，我这里使用了shell命令：

cut -f 1,8 GPL23126.txt |grep -v "^!"|grep -v "^#" >  GPL23126_ids.txt

然后又是R里面的解析即可，如下所示的代码：


library(data.table)
library(stringr)
# b=fread('GPL23126.txt',data.table = F)
# cut -f 1,8 GPL23126.txt |grep -v "^!"|grep -v "^#" >  GPL23126_ids.txt
b=fread('GPL23126_ids.txt',data.table = F)
head(b)
ids=data.frame(
  ID=b$ID,
  symbol=str_split(b$gene_assignment,' // ',simplify = T)[,2]
)
head(ids)
ids=ids[ids$symbol != '',]
dat=dat[rownames(dat) %in% ids$ID,]
ids=ids[match(rownames(dat),ids$ID),]
head(ids) 

colnames(ids)=c('probe_id','symbol')  
ids$probe_id=as.character(ids$probe_id)
rownames(dat)=ids$probe_id
dat[1:4,1:4] 

ids=ids[ids$probe_id %in%  rownames(dat),]
dat[1:4,1:4]   
dat=dat[ids$probe_id,]

后续也是常规分析，如下所示：