[GEO]下载表达矩阵并检查是否完整

2019-05-09 本文已影响0人小洁忘了怎么分身

今天在备课GEO流程，想想下个周的今天我已经讲完了好开心。
本文参考生信技能树的《从GEO数据库下载得到表达矩阵一文就够》。
只想用jimmy大神的口音吐槽一句:“碗束实在是太差了，科研条件太恶略了”。
这不，下载个GEO数据，没有两三次是搞不定的，严重限制了我的进步。
回到正题。

如何获取表达矩阵

大神写过函数downGSE，下载数据、导出表达矩阵和分组信息的csv一气呵成。

downGSE <- function(studyID, destdir = ".") {

    library(GEOquery)
    eSet <- getGEO(studyID, destdir = destdir, getGPL = F)

    exprSet = exprs(eSet[[1]])
    pdata = pData(eSet[[1]])

    write.csv(exprSet, paste0(studyID, "_exprSet.csv"))
    write.csv(pdata, paste0(studyID, "_metadata.csv"))
    return(eSet)

}
downGSE("GSE32575")

getGEO下载的是一个eSet.Rdata和一个matrix.txt压缩文件。当碗束不给力的时候可以选择只下载matrix。
在网页上找到他就是：

一个两个还好，可以一个个复制，当做的多了就会觉得麻烦，于是大神再次出手，用get_GSE_links函数可以获取到这个连接。

get_GSE_links <- 
  function(studyID, down = F, destdir = "./") {
    ## studyID destdir ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE1nnn/GSE1009/matrix/GSE1009_series_matrix.txt.gz
    ## ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE1nnn/GSE1009/suppl/GSE1009_RAW.tar
    ## http://www.ncbi.nlm.nih.gov/geo/browse/?view=samples&mode=csv&series=1009
    
    supp_link = paste0("ftp://ftp.ncbi.nlm.nih.gov/geo/series/", substr(studyID, 1, nchar(studyID) - 3), "nnn/", 
                       studyID, "/suppl/", studyID, "_RAW.tar")
    meta_link = paste0("http://www.ncbi.nlm.nih.gov/geo/browse/?view=samples&mode=csv&series=", substr(studyID, 4, nchar(studyID)))
    matrix_link = paste0("ftp://ftp.ncbi.nlm.nih.gov/geo/series/", substr(studyID, 1, nchar(studyID) - 3), "nnn/", 
                         studyID, "/matrix/", studyID, "_series_matrix.txt.gz")
    
    print(paste0("The URL for raw data is : ", supp_link))
    print(paste0("The URL for metadata is : ", meta_link))
    print(paste0("The URL for matrix   is : ", matrix_link))
    
  }
get_GSE_links("GSE32575")

然后用download.file下载，用readr读取，完美。这个更有可能下载的完整一些。

download.file('ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE32nnn/GSE32575/matrix/GSE32575_series_matrix.txt.gz')
expr2 <- read.csv('GSE32575_series_matrix.txt',
                  sep = '\t',
                  row.names = 1,
                  comment.char = "!")