有的难度是想象出来的，这样这样再这样不就解决啦！

用户11414625

发布于 2025-03-28 17:09:32

3700

代码可运行

文章被收录于专栏：生信星球520生信星球520

运行总次数：0

代码可运行

1.问题数据

学员遇到一个有问题的基因表达芯片数据：

可以看到这个数据有两个问题：数值大小都在0附近，里面还有NA

2.探索原因

在GSE页面上点进去一个样本去看，可以看到这些值的说明，是标准化之后的数据。

标准化是不可逆的，不能用于做差异分析。

3.解决问题

拿到这样的数据，要么放弃，要么处理原始数据。

我的学生选择的是处理原始数据，可是呢，这个bgx数据格式我们没教欸！

4.所以什么是bgx呢？

只要略略搜索就可以知道，它是一种探针注释存储的格式。

https://support.bioconductor.org/p/69959/

到这里其实就可以打住了，因为我们并不需要使用这个探针注释。虽然这个GPL确实没有对应的注释R包，但从GSE页面上点击GPL编号进去，可以看到是有探针注释的。

（截图的时候没看到，其实后面还有一列叫symbol的，还是用symbol吧）

所以我们其实根本不需要bgx格式。但你却因为不会处理bgx，而止步不前咯。这不就是想象出来的难度嘛。

读取一下看看也行吧，中国人讲究一个来都来了。

#BiocManager::install("illuminaio")
library(illuminaio)
bgx = readBGX("GPL6883_HumanRef-8_V3_0_R0_11282963_A.bgx")

其中probes点开就是:

没啥用处。和GPL页面的表格一样，还少了探针ID列呢。

5.真诚是永远的必杀技

我就给你直接分析完得了呗！省的后面又有啥不会。

原始数据有两个，除了bgx，另一个txt文件就是表达矩阵了，这个表格长得还挺贴心的

1.表达矩阵

exp = read.delim("GSE16561_RAW.txt",row.names = 1,check.names = F)
range(exp)

## [1]   111.8617 72479.9900

exp = log2(exp+1)
exp[1:4,1:4]

##              3100083_Stroke 3100191_Stroke 3100068_Stroke 3100060_Stroke
## ILMN_1809034       7.884799       7.483076       7.673048       7.674788
## ILMN_1660305       8.195663       7.840273       7.607785       7.997668
## ILMN_1762337       7.490086       7.403282       7.185134       7.302428
## ILMN_2055271       7.689259       7.683558       7.448607       7.903191

boxplot(exp)

还挺不齐的，给他拉齐一下：

exp = limma::normalizeBetweenArrays(exp)
boxplot(exp)

2.临床信息

geo_download可以下载series.matrix.gz文件，解析里面的表达矩阵和分组信息以及GPL编号。虽然这个数据表达矩阵有问题，但临床信息表格没问题，转录组数据的临床信息也可以这样子读取。

GPL编号我们不需要了，直接下载他的探针注释文件，提取对应的列即可。

library(tinyarray)

pd = geo_download("GSE16561")$pd

## Warning in geo_download("GSE16561"): NA or NAN values detected

可以看到，title列和txt里的表达矩阵列名是一致的。

identical(colnames(exp),pd$title)

## [1] TRUE

3.分组信息

同一个分组对应同一个关键词，levels 设置对照组在前。

Group = factor(pd$description,levels = c("Control","Stroke"))
table(Group)

## Group
## Control  Stroke 
##      24      39

4.探针注释

a = data.table::fread("GPL6883-11606.txt")
colnames(a)

##  [1] "ID"                    "Species"               "Source"               
##  [4] "Search_Key"            "Transcript"            "ILMN_Gene"            
##  [7] "Source_Reference_ID"   "RefSeq_ID"             "Entrez_Gene_ID"       
## [10] "GI"                    "Accession"             "Symbol"               
## [13] "Protein_Product"       "Array_Address_Id"      "Probe_Type"           
## [16] "Probe_Start"           "SEQUENCE"              "Chromosome"           
## [19] "Probe_Chr_Orientation" "Probe_Coordinates"     "Cytoband"             
## [22] "Definition"            "Ontology_Component"    "Ontology_Process"     
## [25] "Ontology_Function"     "Synonyms"              "GB_ACC"

ids = a[,c("ID","Symbol")]
colnames(ids) =c("probe_id","symbol")

5.差异分析

dcp = get_deg_all(exp,Group,ids,entriz = F,logFC_cutoff = 0.585,cluster_cols = F)
dcp$plots

本文参与腾讯云自媒体同步曝光计划，分享自微信公众号。

原始发表：2025-03-27，如有侵权请联系 cloudcommunity@tencent.com 删除

表格

本文分享自生信星球微信公众号，前往查看

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

本文参与腾讯云自媒体同步曝光计划，欢迎热爱写作的你一起参与！

登录后参与评论

0 条评论

热度