1、R包安装

所使用R 的版本为 4.4.1

备注:低版本的R可能安装不上“GSEABase","GSVA","clusterProfiler”
# GEO数据分析,安装相关包
# GEO数据分析所用到的包如下

#BiocManager 
#GEOquery
#KEGG.db
#GSEABase","GSVA","clusterProfiler
#limma","impute
#genefu","org.Hs.eg.db","hgu133plus2.db
#设置当前工作目录
#setwd(r"[D:/Bio_info/GEO Data Explore/CEO_EXPLORE_SWW/GEO]")

rm(list = ls())#删除工作空间中所有的对象

#设置R包下载的镜像
#options("repos" = c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/")) #对应清华源
#options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/") #对应中科大源

#查看目前所使用的镜像
#options()$repos
#options()$BioC_mirror

rm(list = ls())
options()$repos
options()$Bioc_mirror # 查看当前R所使用的镜像
options("repos" =  c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/")) #对应清华源
#options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/") #对应中科大源
options(BioC_mirror="https://mirrors.westlake.edu.cn/bioconductor")
options()$repos
options()$Bioc_mirror # 查看当前R所使用的镜像

#安装所用的包
if (!requireNamespace("BiocManager", quietly = TRUE))
  install.packages("BiocManager")
BiocManager::install("KEGG.db",ask = F,update = F)
BiocManager::install(c("GSEABase","GSVA","clusterProfiler" ),ask = F,update = F)
BiocManager::install(c("GEOquery","limma","impute" ),ask = F,update = F)
BiocManager::install(c("genefu","org.Hs.eg.db","hgu133plus2.db" ),ask = F,update = F)

install.packages('WGCNA')
install.packages(c("FactoMineR", "factoextra"))
install.packages(c("ggplot2", "pheatmap","ggpubr"))

#加载相应的包
library(GSEABase)
library(GSVA)
library(clusterProfiler)
library(genefu)
library(ggplot2)
library(ggpubr)
library(hgu133plus2.db)
library(limma)
library(org.Hs.eg.db)
library(pheatmap)

2、数据下载


rm(list = ls())#删除工作空间中所有的对象

options(stringsASFactors =F) # 在调用as.data.frame 时,将StringAsFactors 设置为F,避免将charactor 转换成factor类型
#注意查看下载文件的大小,检查数据 

GEO_file1 = 'GSE42872_eSet.Rdata'

#需要注意的是,如果你没有设置这个选项,GEOquery会默认使用基于curl和wget的自动选择下载方法。
#设置这个选项为'auto',会使得GEOquery在可能的情况下优先使用基于curl的下载方法,因为它通常更快,且对防火墙的要求较低。
## 这个包需要注意两个配置,一般来说自动化的配置是足够的。
#Setting.options('download.file.method.GEOquery' = 'auto')
#setting.options('GEOquery.inmemory.gpl' = FALSE)

#加载GEOquery包
library(GEOquery)

#获取数据
if( !file.exists(GEO_file1) ){
  #getGEO函数‌是GEOquery包中的一个重要函数,用于从GEO数据库下载数据。该函数的参数包括:
  #GE:决定下载的数据种类,可以是GDS/GSE/GSM/GPL,具体取决于提供的GEO参数。
  #filenam:如果已经下载好了文件,可以直接读取。
  #destdi:下载文件存放的地址,默认为工作目录。
  #GSElimit:用于限制下载的GSE数据集的数量。
  #GSEMatri:若为TRUE,则下载Matrix(GSE数据集的表达矩阵)文件;若为FALSE,则下载SOFT文件。默认为TRUE。
  #AnnotGP:若为TRUE,则下载GPL注释文件;若为FALSE,则不下载。默认为FALSE。
  #getGP:若为TRUE,则下载GPL注释文件;若为FALSE,则不下载。默认为TRUE。
  #parseCharacteristic:未在搜索结果中明确说明,可能是与解析数据特性相关的参数。
  file1 <- getGEO("GSE42872",
                  destdir = ".",
                  AnnotGPL = F, #注释文件
                  getGPL = F)  # 平台文件
  save(file1,file = GEO_file1) ## 保存文件到本地
}

3、数据分析

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

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

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

4、下载GPL注释文件

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

 

 # 字符串切割获取id
geneid_list <- strsplit(Table(gpl)$gene_assignment,"//")
#class(geneid_list) #list

# 获取基因名
gene_Symbol <- lapply(geneid_list,\(x) ifelse(length(x) < 2,x[[1]], x[[2]]))

Table(gpl)$gene_assignment

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

Logo

腾讯云面向开发者汇聚海量精品云计算使用和开发经验,营造开放的云计算技术生态圈。

更多推荐