一、GEO数据库数据下载和初步分析

1、R包安装

所使用R 的版本为 4.4.1

备注：低版本的R可能安装不上“GSEABase","GSVA","clusterProfiler”

1. # GEO数据分析，安装相关包
2. # GEO数据分析所用到的包如下
4. #BiocManager
5. #GEOquery
6. #KEGG.db
7. #GSEABase","GSVA","clusterProfiler
8. #limma","impute
9. #genefu","org.Hs.eg.db","hgu133plus2.db
10. #设置当前工作目录
11. #setwd(r"[D:/Bio_info/GEO Data Explore/CEO_EXPLORE_SWW/GEO]")
13. rm(list = ls())#删除工作空间中所有的对象
15. #设置R包下载的镜像
16. #options("repos" = c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/")) #对应清华源
17. #options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/") #对应中科大源
19. #查看目前所使用的镜像
20. #options()$repos
21. #options()$BioC_mirror
23. rm(list = ls())
24. options()$repos
25. options()$Bioc_mirror # 查看当前R所使用的镜像
26. options("repos" =  c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/")) #对应清华源
27. #options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/") #对应中科大源
28. options(BioC_mirror="https://mirrors.westlake.edu.cn/bioconductor")
29. options()$repos
30. options()$Bioc_mirror # 查看当前R所使用的镜像
32. #安装所用的包
33. if (!requireNamespace("BiocManager", quietly = TRUE))
34.   install.packages("BiocManager")
35. BiocManager::install("KEGG.db",ask = F,update = F)
36. BiocManager::install(c("GSEABase","GSVA","clusterProfiler" ),ask = F,update = F)
37. BiocManager::install(c("GEOquery","limma","impute" ),ask = F,update = F)
38. BiocManager::install(c("genefu","org.Hs.eg.db","hgu133plus2.db" ),ask = F,update = F)
40. install.packages('WGCNA')
41. install.packages(c("FactoMineR", "factoextra"))
42. install.packages(c("ggplot2", "pheatmap","ggpubr"))
44. #加载相应的包
45. library(GSEABase)
46. library(GSVA)
47. library(clusterProfiler)
48. library(genefu)
49. library(ggplot2)
50. library(ggpubr)
51. library(hgu133plus2.db)
52. library(limma)
53. library(org.Hs.eg.db)
54. library(pheatmap)

复制代码

2、数据下载

2. rm(list = ls())#删除工作空间中所有的对象
4. options(stringsASFactors =F) # 在调用as.data.frame 时，将StringAsFactors 设置为F,避免将charactor 转换成factor类型
5. #注意查看下载文件的大小，检查数据
7. GEO_file1 = 'GSE42872_eSet.Rdata'
9. #需要注意的是，如果你没有设置这个选项，GEOquery会默认使用基于curl和wget的自动选择下载方法。
10. #设置这个选项为'auto'，会使得GEOquery在可能的情况下优先使用基于curl的下载方法，因为它通常更快，且对防火墙的要求较低。
11. ## 这个包需要注意两个配置，一般来说自动化的配置是足够的。
12. #Setting.options('download.file.method.GEOquery' = 'auto')
13. #setting.options('GEOquery.inmemory.gpl' = FALSE)
15. #加载GEOquery包
16. library(GEOquery)
18. #获取数据
19. if( !file.exists(GEO_file1) ){
20.   #getGEO函数‌是GEOquery包中的一个重要函数，用于从GEO数据库下载数据。该函数的参数包括：
21.   #GE：决定下载的数据种类，可以是GDS/GSE/GSM/GPL，具体取决于提供的GEO参数。
22.   #filenam：如果已经下载好了文件，可以直接读取。
23.   #destdi：下载文件存放的地址，默认为工作目录。
24.   #GSElimit：用于限制下载的GSE数据集的数量。
25.   #GSEMatri：若为TRUE，则下载Matrix（GSE数据集的表达矩阵）文件；若为FALSE，则下载SOFT文件。默认为TRUE。
26.   #AnnotGP：若为TRUE，则下载GPL注释文件；若为FALSE，则不下载。默认为FALSE。
27.   #getGP：若为TRUE，则下载GPL注释文件；若为FALSE，则不下载。默认为TRUE。
28.   #parseCharacteristic：未在搜索结果中明确说明，可能是与解析数据特性相关的参数。
29.   file1 <- getGEO("GSE42872",
30.                   destdir = ".",
31.                   AnnotGPL = F, #注释文件
32.                   getGPL = F)  # 平台文件
33.   save(file1,file = GEO_file1) ## 保存文件到本地
34. }

复制代码

3、数据分析

1. #因为这个GEO数据集只有一个GPL平台，所以下载到的是一个含有一个元素的list
2. length(file1)
3. class(file1)
4. #使用exprs函数：exprs函数可以从矩阵、数组或数据框等多种数据结构中提取基因表达值。
5. #例如，如果expression_set是一个包含基因表达数据的数据结构，那么通过exprs(expression_set)可以提取出基因表达值，并返回一个向量、矩阵或数组形式的输出
6. b <- file1[[1]]

复制代码

1. #提取的数据如下
2. Geodata = exprs(file1[[1]])

复制代码

1. class(Geodata) #[1] "matrix" "array"
2. # GEO数据集的每一列代表一个样本。‌
3. # 在GEO数据集中，每一列并不是直接对应于某个特定的生物样本属性，而是代表一个具体的样本。
4. # 这些数据集通常包含大量的基因表达数据，其中每一行代表一个基因，而每一列则代表一个样本。
5. # 这样的数据结构使得研究人员能够分析特定基因在不同样本（如不同组织、不同条件下的样本）中的表达情况
6. # ，从而研究基因表达与特定生物学过程或疾病状态之间的关系。
7. dim(Geodata)  
8. #  [1] 33297     6  
9. # 33297个基因，6个样本

复制代码

4、下载GPL注释文件

1. # GPL6244  GPL平台注释
2. if(F){
3.   library(GEOquery)
4.   #Download GPL file, put it in the current directory, and load it:
5.   gpl <- getGEO('GPL6244', destdir=".")
6.   colnames(Table(gpl))  
7.   # [1] "ID"              "GB_LIST"         "SPOT_ID"         "seqname"         "RANGE_GB"        "RANGE_STRAND"    "RANGE_START"   
8.   # [8] "RANGE_STOP"      "total_probes"    "gene_assignment" "mrna_assignment" "category"
9.   # 正常应该有gene-Sympol 这一列，但是现在没有，所以需要切割字符串，获得geneid
10.   class(head(Table(gpl)))
11.   #[1] "data.frame"
12.   
13.   # 字符串切割获取id
14.   geneid_list <- strsplit(Table(gpl)$gene_assignment,"//")
15.   #class(geneid_list) #list
16.   gene_Symbol <- lapply(geneid_list,\(x) ifelse(length(x) < 2,x[[1]], x[[2]]))
17.   
18.   #获取探针名
19.   probe_id <- Table(gpl)$ID
20.   
21.   #class(gene_Symbol) list
22.   #class(probe_id) Integer
23.   
24.   #合并两个数据，保存在数据框中
25.   probetogene <- data.frame(
26.     `probe_id` <-probe_id,
27.     `gene_Symbol`<- unlist(gene_Symbol)
28.   )
29.   save(probetogene,file='probetogene.Rdata')
30. }

复制代码

1. class((Table(gpl))
2. #[1] "data.frame"
3. Table(gpl)

复制代码

 

1. # 字符串切割获取id
2. geneid_list <- strsplit(Table(gpl)$gene_assignment,"//")
3. #class(geneid_list) #list
5. # 获取基因名
6. gene_Symbol <- lapply(geneid_list,\(x) ifelse(length(x) < 2,x[[1]], x[[2]]))
8. Table(gpl)$gene_assignment

复制代码

注释：本文章是本人学习 B站 生信技能树-jimmy 课程 整理所得，如有错误，请见谅；如有侵权，请接洽本人，本人会立刻删除，谢谢。

免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！更多信息从访问主页：qidao123.com:ToB企服之家，中国第一个企服评测及商务社交产业平台。