TCGA | GEO | 文獻(xiàn)閱讀 | 數(shù)據(jù)庫(kù) | 理論知識(shí)

R語(yǔ)言 | Bioconductor

關(guān)于GEO數(shù)據(jù)庫(kù)的基本認(rèn)識(shí)和在線使用，可閱讀文章：GEO數(shù)據(jù)庫(kù)使用教程及在線數(shù)據(jù)分析工具。本文介紹GEOquery這個(gè)包，官方教程地址：

http://www./packages/release/bioc/vignettes/GEOquery/inst/doc/GEOquery.html

本文相當(dāng)于官方教程的譯文吧，再加上我個(gè)人的理解，可能翻譯不是那么準(zhǔn)確。請(qǐng)多多指教??！

1 GEO概述

NCBI基因表達(dá)綜合庫(kù)(GEO)可以用作各種高通量實(shí)驗(yàn)數(shù)據(jù)的公共存儲(chǔ)庫(kù)。這些數(shù)據(jù)包括測(cè)量mRNA，基因組DNA和蛋白質(zhì)豐度的基于單通道和雙通道微陣列的實(shí)驗(yàn)，以及非陣列技術(shù)，例如基因表達(dá)的串行分析(SAGE)，質(zhì)譜蛋白質(zhì)組學(xué)數(shù)據(jù)和高通量測(cè)序數(shù)據(jù)。

在GEO的最基本組織層次上，有四種基本實(shí)體類型。前三個(gè)(樣本，平臺(tái)和系列)由用戶提供；第四個(gè)是數(shù)據(jù)集，由GEO員工根據(jù)用戶提交的數(shù)據(jù)進(jìn)行編譯和整理。有關(guān)更多信息，請(qǐng)參見GEO主頁(yè)(https://www.ncbi.nlm./geo/)。

1.1平臺(tái)

平臺(tái)記錄描述了陣列上的元素列表(例如cDNA，寡核苷酸探針集，ORF，抗體)，或在該實(shí)驗(yàn)中可以檢測(cè)和定量的元素列表(例如SAGE標(biāo)簽，肽)。每個(gè)平臺(tái)記錄都分配有一個(gè)唯一且穩(wěn)定的GEO登錄號(hào)(GPLxxx)。平臺(tái)可以引用由多個(gè)提交者提交的許多樣本。

1.2樣品

樣品記錄描述了處理單個(gè)樣品的條件，進(jìn)行的操作以及從中得出的每個(gè)元素的豐度測(cè)量。每個(gè)樣本記錄都分配有一個(gè)唯一且穩(wěn)定的GEO登錄號(hào)(GSMxxx)。樣本實(shí)體只能引用一個(gè)平臺(tái)，并且可以包含在多個(gè)系列中。

1.3系列

系列記錄定義了一組相關(guān)的樣本，這些樣本被視為組的一部分，這些樣本的關(guān)聯(lián)方式以及它們是否有序和如何排序。系列作為一個(gè)整體提供了實(shí)驗(yàn)的重點(diǎn)和描述。系列記錄還可以包含描述提取的數(shù)據(jù)、總結(jié)結(jié)論或分析的表格。每個(gè)系列記錄都分配有一個(gè)唯一且穩(wěn)定的GEO登錄號(hào)(GSExxx)。系列記錄有兩種格式，由GEOquery獨(dú)立處理。較小的新GSEMatrix文件解析起來(lái)非常快。GEOquery使用一個(gè)簡(jiǎn)單的標(biāo)志來(lái)選擇使用GSEMatrix文件(請(qǐng)參見下文)。

1.4數(shù)據(jù)集

GEO數(shù)據(jù)集(GDSxxx)是GEO樣本數(shù)據(jù)的精選集。GDS記錄代表了生物學(xué)和統(tǒng)計(jì)上可比的GEO樣本的集合，并構(gòu)成了GEO數(shù)據(jù)顯示和分析工具套件的基礎(chǔ)。GDS中的樣本引用相同的平臺(tái)，也就是說(shuō)，它們共享一組通用的探針元素。假定以等效方式計(jì)算GDS中每個(gè)樣本的值測(cè)量值，也就是說(shuō)，諸如背景處理和規(guī)范化之類的考慮在整個(gè)數(shù)據(jù)集中是一致的。通過(guò)GDS子集提供反映實(shí)驗(yàn)設(shè)計(jì)的信息。

2.開始使用GEOquery

從GEO獲取數(shù)據(jù)確實(shí)非常容易。只需一個(gè)命令getGEO。這個(gè)函數(shù)解釋它的輸入以確定如何從GEO獲取數(shù)據(jù)，然后將數(shù)據(jù)解析成有用的R數(shù)據(jù)結(jié)構(gòu)。用法非常簡(jiǎn)單。

如果沒有安裝GEOquery包的話，先安裝。

if (!requireNamespace("BiocManager", quietly = TRUE))    install.packages("BiocManager")BiocManager::install("GEOquery")

加載這個(gè)包。

library(GEOquery)

現(xiàn)在，我們可以自由地訪問(wèn)任何GEO文件。下面代碼是使用GEOquery包打包的文件，而不是從網(wǎng)上下載。通常，我們會(huì)使用GEO登錄號(hào)下載數(shù)據(jù)，如代碼注釋中所述。

# If you have network access, the more typical way to do this# would be to use this:# gds gds "extdata/GDS507.soft.gz",package=

也就是說(shuō)，我們要從網(wǎng)絡(luò)上下載數(shù)據(jù)的話，使用下面代碼就行了。

gds "GDS507")

同樣的，我們也可以通過(guò)樣品的登錄號(hào)獲取。

# If you have network access, the more typical way to do this# would be to use this:# gds gsm "extdata/GSM11805.txt.gz",package=

3.GEOquery數(shù)據(jù)結(jié)構(gòu)

GEOquery數(shù)據(jù)結(jié)構(gòu)實(shí)際上有兩種形式。第一種，包括GDS、GPL和GSM，這三種數(shù)據(jù)結(jié)構(gòu)比較類似，getGEO對(duì)他們請(qǐng)求是也類似。第四個(gè)GEOquery數(shù)據(jù)結(jié)構(gòu)是GSE，GSE是由GSM和GPL對(duì)象組合而成的復(fù)合數(shù)據(jù)類型。

3.1 GDS、GSM和GPL類

這些類中的每一個(gè)都由元數(shù)據(jù)標(biāo)頭(幾乎從SOFT格式標(biāo)頭中逐字獲取)和GEODataTable組成。GEODataTable有兩個(gè)簡(jiǎn)單的部分，一個(gè)Columns部分，它描述Table部分的列標(biāo)題。show每個(gè)類都有一個(gè)方法。例如，使用上面的gsm：

> head(Meta(gsm))$channel_count[1] "1"$comment[1] "Raw data provided as supplementary file"$contact_address[1] "715 Albany Street, E613B"$contact_city[1] "Boston"$contact_country[1] "USA"$contact_department[1] "Genetics and Genomics"

> Table(gsm)[1:5,]          ID_REF  VALUE ABS_CALL1 AFFX-BioB-5_at  953.9        P2 AFFX-BioB-M_at 2982.8        P3 AFFX-BioB-3_at 1657.9        P4 AFFX-BioC-5_at 2652.7        P5 AFFX-BioC-3_at 2019.5        P

> Columns(gsm)    Column1   ID_REF2    VALUE3 ABS_CALL                                                                 Description1                                                                           2                         MAS 5.0 Statistical Algorithm (mean scaled to 500)3 MAS 5.0 Absent, Marginal, Present call  with Alpha1 = 0.05, Alpha2 = 0.065

GPL類的行為與GSM類完全相同。但是，GDS類有更多與列方法關(guān)聯(lián)的信息。

> Columns(gds)[,1:3]     sample disease.state individual1  GSM11815           RCC        0352  GSM11832           RCC        0233  GSM12069           RCC        0014  GSM12083           RCC        0055  GSM12101           RCC        0116  GSM12106           RCC        0327  GSM12274           RCC          28  GSM12299           RCC          39  GSM12412           RCC          410 GSM11810        normal        03511 GSM11827        normal        02312 GSM12078        normal        00113 GSM12099        normal        00514 GSM12269        normal          115 GSM12287        normal          216 GSM12301        normal          317 GSM12448        normal          4

3.2 GSE類

GSE實(shí)例是GEO實(shí)例中最容易混淆的。一個(gè)GSE條目可以表示在任意數(shù)量的平臺(tái)上運(yùn)行的任意數(shù)量的樣本。與其他類一樣，GSE類有一個(gè)元數(shù)據(jù)部分。但是，它沒有GEODataTable。相反，它包含兩個(gè)列表，可以使用GPLList和GSMList方法訪問(wèn)，這兩個(gè)列表分別是GPL和GSM對(duì)象的列表。

# Again, with good network access, one would do:# gse "GSE781",GSEMatrix=FALSE)gse "extdata/GSE781_family.soft.gz",package=

> head(Meta(gse))$contact_address[1] "715 Albany Street, E613B"$contact_city[1] "Boston"$contact_country[1] "USA"$contact_department[1] "Genetics and Genomics"$contact_email[1] "mlenburg@bu.edu"$contact_fax[1] "617-414-1646"

Meta(gse)提取的信息就是網(wǎng)頁(yè)上的信息。你可以網(wǎng)頁(yè)上查看該數(shù)據(jù)集看看。

https://www.ncbi.nlm./geo/query/acc.cgi?acc=GSE781

我們可以利用GSMList函數(shù)提取GSM對(duì)象。

# GSE中包含的所有GSM對(duì)象的名稱names(GSMList(gse))

# 并獲取列表中的第一個(gè)GSM對(duì)象GSMList(gse)[[1]]

同樣，GPLList函數(shù)提取GPL對(duì)象。

names(GPLList(gse))

4. 轉(zhuǎn)化為 BioConductor的 ExpressionSets 和 limma MALists對(duì)象

GEO數(shù)據(jù)集(與其他一些GEO實(shí)例不同)非常類似于LIMMA數(shù)據(jù)結(jié)構(gòu)的MAList對(duì)象和BioBase數(shù)據(jù)結(jié)構(gòu)的ExpressionSet對(duì)象。因此，有兩個(gè)函數(shù)GDS2MA和GDS2eSet可以完成轉(zhuǎn)換任務(wù)。

4.1 獲取作為ExpressionSets的GSE系列矩陣文件

GEO系列是相關(guān)實(shí)驗(yàn)的集合。除了可以作為相當(dāng)大的軟格式文件提供之外，NCBI GEO還準(zhǔn)備了一個(gè)基于制表符分隔文本的更簡(jiǎn)單的格式文件。getGEO函數(shù)可以處理這種格式，并且可以相當(dāng)快地解析非常大的GSE。此解析返回的數(shù)據(jù)結(jié)構(gòu)是ExpressionSet列表。作為示例，我們下載并解析GSE2553。

gse2553 'GSE2553',GSEMatrix=TRUE)

> show(gse2553)$GSE2553_series_matrix.txt.gzExpressionSet (storageMode: lockedEnvironment)assayData: 12600 features, 181 samples   element names: exprs protocolData: nonephenoData  sampleNames: GSM48681 GSM48682 ... GSM48861 (181 total)  varLabels: title geo_accession ... data_row_count (30 total)  varMetadata: labelDescriptionfeatureData  featureNames: 1 2 ... 12600 (12600 total)  fvarLabels: ID PenAt ... Chimeric_Cluster_IDs (13 total)  fvarMetadata: Column Description labelDescriptionexperimentData: use 'experimentData(object)'  pubMedIds: 16230383 Annotation: GPL1977 

> show(pData(phenoData(gse2553[[1]]))[1:5,c(1,6,8)])                                                                 titleGSM48681                      Patient sample ST18, DermatofibrosarcomaGSM48682                           Patient sample ST410, Ewing SarcomaGSM48683                            Patient sample ST130, Sarcoma, NOSGSM48684 Patient sample ST293, Malignant Peripheral Nerve Sheath TumorGSM48685                             Patient sample ST367, Liposarcoma         type                         source_name_ch1GSM48681  RNA                     DermatofibrosarcomaGSM48682  RNA                           Ewing SarcomaGSM48683  RNA                            Sarcoma, NOSGSM48684  RNA Malignant Peripheral Nerve Sheath TumorGSM48685  RNA                             Liposarcoma

4.2 將GDS轉(zhuǎn)換為ExpressionSet

eset do.log2=TRUE)

現(xiàn)在，eset是一個(gè)ExpressionSet包含相同的信息作為GEO數(shù)據(jù)集，包括樣品信息，我們可以在這里看到：

> head(eset)[,1:5]                GSM1071862 GSM1071863 GSM1071864 GSM1071865 GSM10718661053_3p_at        2.453620   2.045711   2.780113   2.417370   2.260668117_3p_at         4.745937   4.030957   5.365418   5.536504   3.2758811494_3p_f_at      5.756770   5.614358   7.382213   6.912064   7.0344211552275_3p_s_at   2.070568   2.657432   3.915170   3.250838   2.4281281552281_3p_at     6.523435   6.436356   7.873713   6.567717   6.9138351552296_3p_at     1.614668   3.412571   1.741802   2.345888   4.213399

4.3 轉(zhuǎn)換GDS為MAList

ExpressionSet通常未獲取任何注釋信息(GEO稱其為平臺(tái)信息)，但是，很容易獲得此信息。首先，我們需要知道此GDS使用的平臺(tái)。然后，再次通過(guò)getGEO將獲得我們所需的東西。

> Meta(gds)$platform[1] "GPL97"

gpl "extdata/GPL97.annot.gz",package="GEOquery"))

因此，gpl現(xiàn)在包含來(lái)自GEO的GPL5信息。與ExpressionSetlimma 不同，limma MAList確實(shí)存儲(chǔ)了基因注釋信息，因此我們可以利用GDS2MA來(lái)新創(chuàng)建含gpl的GPL類。

MA 

> class(MA)[1] "MAList"attr(,"package")[1] "limma"

現(xiàn)在，MA屬于MAList類，不僅包含數(shù)據(jù)，還包含與GDS507相關(guān)的樣本信息和基因信息。

ano "genes"]]

> ano[1:5,1:4]           ID                                         Gene title     Gene symbol          Gene ID1 200000_s_at                       pre-mRNA processing factor 8           PRPF8            105942   200001_at                            calpain small subunit 1          CAPNS1              8263   200002_at                              ribosomal protein L35           RPL35            112244 200003_s_at              microRNA 6805///ribosomal protein L28 MIR6805///RPL28 102465483///61585   200004_at eukaryotic translation initiation factor 4 gamma 2          EIF4G2             1982

4.4 GSE轉(zhuǎn)化為 ExpressionSet

gsmplatforms function(x) {Meta(x)$platform_id})

> head(gsmplatforms)$GSM11805[1] "GPL96"$GSM11810[1] "GPL97"$GSM11814[1] "GPL96"$GSM11815[1] "GPL97"$GSM11823[1] "GPL96"$GSM11827[1] "GPL97"

確實(shí)，有兩個(gè)GPL，即GPL96和GPL97，作為它們的平臺(tái)(我們可以通過(guò)查看GPLList來(lái)確定它們gse)。我們可以過(guò)濾原始的GSMList，使其僅包含那些具有GPL96平臺(tái)的GSM，并將此列表用于進(jìn)一步處理。

gsmlist = Filter(function(gsm) {Meta(gsm)$platform_id=='GPL96'},GSMList(gse))

> length(gsmlist)[1] 17

所以，現(xiàn)在我們想知道哪一列代表我們想要提取的數(shù)據(jù)。查看單個(gè)GSM的表的前幾行可能會(huì)給我們一個(gè)概念(順便說(shuō)一句，GEO使用了一個(gè)約定，即包含每個(gè)數(shù)組的單個(gè)度量的列稱為值列，如果我們不知道其他哪些列最相關(guān)，可以使用它)。

> Table(gsmlist[[1]])[1:5,]          ID_REF  VALUE ABS_CALL1 AFFX-BioB-5_at  953.9        P2 AFFX-BioB-M_at 2982.8        P3 AFFX-BioB-3_at 1657.9        P4 AFFX-BioC-5_at 2652.7        P5 AFFX-BioC-3_at 2019.5        P

> Columns(gsmlist[[1]])[1:5,]       Column                                                                Description1      ID_REF                                                                           2       VALUE                         MAS 5.0 Statistical Algorithm (mean scaled to 500)3    ABS_CALL MAS 5.0 Absent, Marginal, Present call  with Alpha1 = 0.05, Alpha2 = 0.065NA       NA.1     

我們將使用VALUE列。然后，我們想做一個(gè)矩陣，這些值如下:

# get the probeset orderingprobesets 1]])$ID# make the data matrix from the VALUE columns from each GSM# being careful to match the order of the probesets in the platform# with those in the GSMsdata.matrix do.call({tab mymatch return(tab$VALUE[mymatch])}))data.matrix 2,data.matrix 

> data.matrix[1:5,1:5]      GSM11805  GSM11814  GSM11823  GSM11830  GSM12067[1,] 10.926963 11.105254 11.275019 11.438636 11.424376[2,]  5.749534  7.908092  7.093814  7.514122  7.901470[3,]  7.066089  7.750205  7.244126  7.962896  7.337176[4,] 12.660353 12.479755 12.215897 11.458355 11.397568[5,]  6.195741  6.061776  6.565293  6.583459  6.877744

require(Biobase)# go through the necessary steps to make a compliant ExpressionSetrownames(data.matrix) colnames(data.matrix) pdata rownames(pdata) pheno as(pdata,eset2 new(eset2

注意，我們要做的match是確保值和平臺(tái)信息的順序相同。最后，制作ExpressionSet對(duì)象。

5.從GEO獲取原始數(shù)據(jù)

NCBI GEO接受(但并非總是需要)原始數(shù)據(jù)，例如.CEL文件，.CDF文件，圖像等。有時(shí)，快速訪問(wèn)此類數(shù)據(jù)很有用。單個(gè)函數(shù)getGEOSuppFiles可以將GEO加入作為參數(shù)，并將下載與該加入相關(guān)的所有原始數(shù)據(jù)。默認(rèn)情況下，該函數(shù)將在當(dāng)前工作目錄中創(chuàng)建一個(gè)目錄，以存儲(chǔ)所選GEO入藏的原始數(shù)據(jù)。結(jié)合使用簡(jiǎn)單的sapply語(yǔ)句或其他循環(huán)結(jié)構(gòu)，getGEOSuppFiles可以以一種非常簡(jiǎn)單的方式快速輕松地獲取原始數(shù)據(jù)，而無(wú)需了解GEO原始數(shù)據(jù)URL的細(xì)節(jié)。

a 'GSM1137', fetch_files = FALSE)

上面內(nèi)容來(lái)自官方文檔的教程，加上我個(gè)人的理解，后續(xù)我們?cè)?strong>GEO專輯實(shí)戰(zhàn)中進(jìn)行應(yīng)用。