虽然现在是高通量测序的时代,但是GEO、ArrayExpress等数据库储存并公开大量的基因表达芯片数据,还是会有大量的需求去处理芯片数据,并且建模或验证自己所研究基因的表达情况,芯片数据的处理也可能是大部分刚学生信的道友入门R语言数据处理的第一次实战,因此准备更新100个基因表达芯片或转录组高通量数据的处理。
获取处理好的数据集、代码、安装包等文章中出现的内容资源,可以添加小编微信:18983376561
数据信息检索可以看到GSE48452是基因表达芯片数据,因此可以使用GEOquery包下
图片
使用GEOquery包下载数据remotes::install_github('ScienceAdvances/using')using::using(tidyverse, GEOquery, magrittr, data.table, AnnoProbe, clusterProfiler, org.Hs.eg.db, org.Mm.eg.db)
注:using作用是一次性加载多个R包,不用写双引号,并且不在屏幕上打印包的加载信息
因为文件太大,在R内下载失败,可通过图片中的方法下载文件,GEOquery::getGEO直接读取本地的文件。
geo_accession <- "GSE48452"eSet <- getGEO(filename=stringr::str_glue('{geo_accession}_series_matrix.txt.gz'), AnnotGPL = F, getGPL = F)gpl <- eSet@annotation
图片
图片
处理表型数据这部分是很关键的,可以筛选一下分组表型信息,只保留自己需要的样本,作为后续分析的样本(根据自己的研究目的筛选符合要求的样本)
pdata <- pData(eSet)pdata %<>% dplyr::mutate( Sample = geo_accession, Sex=str_to_title(`Sex:ch1`), Fibrosis=`fibrosis:ch1`, BMI=`bmi:ch1`, Age=`age:ch1`, Group = case_when( `group:ch1`=='Control'~'Control', `group:ch1`=='Nash'~'NASH', TRUE~NA) ) %>% drop_na(Group) %>% dplyr::select(Sample,Group,Sex,BMI,Fibrosis,Age,everything())处理表达谱数据
数据大小不大于50不需要取log
exprs_mtx[is.na(exprs_mtx)] <- 0exprs_mtx <- exprs(eSet)if(max(exprs_mtx, na.rm = TRUE)<50 | min(exprs_mtx, na.rm = TRUE)<0){ message("基因表达最大值小于50或者最小值小于0不需要log转化")}else { message("基因表达最大值大于50需要log转化") exprs_mtx <- log2(exprs_mtx+1)}probe_exprs <- as.data.table(exprs_mtx, keep.rownames = "ProbeID")探针与基因Symbol对应关系
从AnnoProbe包中获取探针与GeneID对应关系
AnnoProbe::checkGPL(gpl)probe2symbol <- AnnoProbe::idmap(gpl = gpl, type = "bioc", mirror = "tencent", destdir = outdir) %>% dplyr::rename(ProbeID = probe_id, Feature = symbol) # pipe", "bioc", "soft"ID转换
把表达矩阵中的探针名转换为基因名;transid是我写的一个R函数,有需要可以联系我
fdata <- transid(probe2symbol, probe_exprs)保存数据
common_samples <- base::intersect(colnames(fdata),pdata$Sample)fdata %<>% select(all_of(c("Feature",common_samples)))fwrite(fdata, file = stringr::str_glue("{geo_accession}_{gpl}_fdata.csv.gz"))pdata %<>% dplyr::filter(Sample %in% common_samples)fwrite(pdata, file = stringr::str_glue("{geo_accession}_{gpl}_pdata.csv"))本站仅提供存储服务,所有内容均由用户发布,如发现有害或侵权内容,请点击举报。