GEO數(shù)據(jù)庫之芯片數(shù)據(jù)下載

風雨都停了 2018-06-17

展開全文

GEO數(shù)據(jù)庫是NCBI開發(fā)的基因表達數(shù)據(jù)庫，主要接收通過高通量測序、基因芯片等方法獲得的基因表達數(shù)據(jù)——這就方便大家利用他人數(shù)據(jù)發(fā)文章了。

而進行GEO數(shù)據(jù)挖掘的第一步就是進行數(shù)據(jù)下載，但是進入網(wǎng)站點點來下載，查詢搜索工作就不少，下載下來的數(shù)據(jù)還不一定能看懂。有沒有什么方法可以解決這個問題呢？當然有——R包GEOquery！下面就針對芯片數(shù)據(jù)，教大家用GEOquery包完成下載工作。

GEO數(shù)據(jù)

在下載之前要先了解GEO數(shù)據(jù)庫具體存放的四類數(shù)據(jù)：GSE、GDS、GSM、和GPL。

一個GSE Accession對應的是整個研究項目的系列的數(shù)據(jù)，可能涉及不同平臺；

一個GDS Accession對應的一個同一平臺的數(shù)據(jù)集；

一個GSM Accession對應單一樣品的數(shù)據(jù)信息，它只能是單一平臺的數(shù)據(jù)，往往，GSE 和GDS中會包含多個GSM數(shù)據(jù)；

一個GPL Accession，則對應一個platform信息。

R包安裝與加載

GEOquery

## try http:// if https:// URLs are not supported
source('https:///biocLite.R')
biocLite('GEOquery')

Biobase

## try http:// if https:// URLs are not supported
source('https:///biocLite.R')
biocLite('Biobase')

正確加載

library('Biobase')
library('GEOquery')
setwd('F:/GEO') ############有需要可以設置路徑

利用GSE Accession

通過閱讀文獻查找感興趣的GSE Accession，下載對應的表達數(shù)據(jù)和平臺信息等，可以利用GEOquery中的getGEO()函數(shù)下載series_matrix.txt。例如GSE70213：

> gse = getGEO('GSE70213', GSEMatrix =TRUE, destdir = '.', getGPL = T, AnnotGPL = T)
###destdir設置當前目錄，getGPL 和AnnotGPL都設置TRUE，可以下載和獲得平臺的注釋文件

gse為列表數(shù)據(jù)，對應的GSM是單平臺，則length為1，之后分別利用Biobase包中的exprs()、pData()和fData()獲得表達數(shù)據(jù)、樣品處理分組等信息、芯片平臺的設計注釋信息，還可以利用annotation()函數(shù)了解對應GPL Accession，譬如exprs()函數(shù)：

> exprSet=exprs(gse[[1]])
> head(exprSet,2)
GSM1720833 GSM1720834 GSM1720835 GSM1720836 GSM1720837 GSM1720838 GSM1720839 GSM1720840 GSM1720841 GSM1720842
10338001 2041.40800 2200.86100 2323.7600 3216.26300 2362.77500 2195.31800 2013.35900 2146.25800 1785.9460 2067.04100
10338002 63.78059 65.08438 58.3082 75.86145 66.95605 43.81526 49.11361 51.29279 48.9604 42.14286
GSM1720843 GSM1720844 GSM1720845 GSM1720846 GSM1720847 GSM1720848 GSM1720849 GSM1720850 GSM1720851 GSM1720852
10338001 1769.1150 1720.77400 1847.42900 2214.69800 2279.51500 2530.45600 2303.26400 2358.83400 1701.40000 1970.92400
10338002 42.5472 43.48373 64.34628 59.75188 57.48852 60.26423 54.81179 53.70885 57.86877 57.02808
GSM1720853 GSM1720854 GSM1720855 GSM1720856
10338001 1822.78600 2014.26000 1737.84200 2001.73400
10338002 59.26121 55.27306 54.36722 49.43959

注釋信息的獲取可以進行探針和基因的對應，方便后續(xù)分析。通過exprs()、pData()和fData()獲得的數(shù)據(jù)都可以利用write.table等進行文件保存。

利用GDS Accession

GDS數(shù)據(jù)同樣可以利用getGEO()函數(shù)下載soft文件。例如GDS5881：

> gds = getGEO('GDS5881', GSEMatrix =TRUE, destdir = '.', getGPL = T, AnnotGPL = T)
###destdir設置當前目錄

gds可以利用GEOquery包中的Table()獲取表達數(shù)據(jù)，并利用Meta()獲得描述信息，其中Meta(gds)$platform可以獲得GPL Accession。

> exprSet=Table(gds)
> head(exprSet,1)
ID_REF IDENTIFIER GSM1720845 GSM1720846 GSM1720847 GSM1720848 GSM1720849 GSM1720850 GSM1720851 GSM1720852 GSM1720853
1 10344614 Gm2889 48.4971 47.252 39.3331 49.9048 36.8313 41.9501 37.5569 38.1924 46.0668
GSM1720854 GSM1720855 GSM1720856
1 34.689 38.5762 32.2618

> Meta(gset)$platform
[1] 'GPL6246'

針對getGEO返回的gds——GDS數(shù)據(jù)，可以利用GEOquery包中GDS2Set()和GDS2MA()轉(zhuǎn)變?yōu)闉镋xpressionSets 和limma MALists。

> gds2eSet=GDS2eSet(gds)
> MA=GDS2MA(gds)

再針對返回的gds2eSet,利用exprs()、pData()和fData()同樣可以獲得表達數(shù)據(jù)、樣品處理分組信息、芯片平臺的設計注釋信息。返回的MA中涉及大量的描述信息，其中MA$tragets也是樣品處理信息。

利用GSM Accession

利用GSM Accession下載的是單樣本的表達數(shù)據(jù)，例如GSM1720833：

> gsm = getGEO('GSM1720833', GSEMatrix =TRUE, destdir = '.', getGPL = T, AnnotGPL = T)

針對gsm，同樣是利用GEOquery包中的Table()獲取表達數(shù)據(jù)，并利用Meta()獲得描述信息，而獲取對應的GSE Accession 和GPL Accsesion利用Meta(gsm)$series_id和Meta(gsm)$platform_id。

利用GPL Accession

針對芯片平臺，利用GPL Accession下載得到的數(shù)據(jù)是芯片的設計和注釋信息，可以獲得探針組和基因的對應關系，利用Table()函數(shù)可以顯示annotation，例如GPL6246：

> gpl = getGEO('GPL6246', GSEMatrix =TRUE, destdir = '.', getGPL = T, AnnotGPL = T)
> ann=Table(gpl)
> head(ann,2)
       ID          Gene title Gene symbol   Gene ID UniGene title UniGene symbol UniGene ID
1 10344614 predicted gene 2889      Gm2889 100040658                                        
2 10344616                                                                                  
                                                                                                                                                                                                                                                                                                                                                                                                                                                            Nucleotide Title
1 Mus musculus blastocyst blastocyst cDNA, RIKEN full-length enriched library, clone:I1C0009C06 product:hypothetical DeoxyUTP pyrophosphatase/Aspartyl protease, retroviral-type family profile/Retrovirus capsid, C-terminal/Peptidase aspartic/Peptidase aspartic, active site containing protein, full insert sequence///Mus musculus blastocyst blastocyst cDNA, RIKEN full-length enriched library, clone:I1C0042P10 product:hypothetical protein, full insert sequence
2                                                                                                                                                                                                                                                                                                                                                                                                                                                                           
                   GI   GenBank Accession Platform_CLONEID Platform_ORF      Platform_SPOTID Chromosome location
1 74211482///74217103 AK145513///AK145782                               chr1:3054233-3054733                  18
2                                                                       chr1:3102016-3102125                    
  Chromosome annotation GO:Function GO:Process GO:Component GO:Function ID GO:Process ID GO:Component ID
1         Chromosome 18                                                                                 
2