How to Get a Microarray Data From NCBI GEO
GEOquery を使って公共データベースから発現データを入手する
例としてこのデータを取り出します。
準備
1 2 | |
データをダウンロードします。ネットワークからダウンロードするので時間がかかります。この例では数秒から数分ぐらいだと思います。
library("GEOquery")## Loading required package: Biobase## Loading required package: BiocGenerics## Attaching package: 'BiocGenerics'## The following object(s) are masked from 'package:stats':
##
## xtabs## The following object(s) are masked from 'package:base':
##
## anyDuplicated, cbind, colnames, duplicated, eval, Filter, Find, get,
## intersect, lapply, Map, mapply, mget, order, paste, pmax, pmax.int, pmin,
## pmin.int, Position, rbind, Reduce, rep.int, rownames, sapply, setdiff,
## table, tapply, union, unique## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with 'browseVignettes()'. To
## cite Bioconductor, see 'citation("Biobase")', and for packages
## 'citation("pkgname")'.## Setting options('download.file.method.GEOquery'='auto')g <- getGEO("GSE20349")## Found 1 file(s)## GSE20349_series_matrix.txt.gz## File stored at:## /var/folders/cz/cny0ysmx205dnj0y2_34k8cc0000gn/T//RtmpqdGru3/GPL8492.softis(g)## [1] "list" "vector" "AssayData"g## $GSE20349_series_matrix.txt.gz
## ExpressionSet (storageMode: lockedEnvironment)
## assayData: 16331 features, 4 samples
## element names: exprs
## protocolData: none
## phenoData
## sampleNames: GSM509823 GSM509824 GSM509825 GSM509826
## varLabels: title geo_accession ... data_row_count (35 total)
## varMetadata: labelDescription
## featureData
## featureNames: 100009600_at 100012_at ... 99982_at (16331 total)
## fvarLabels: ID Gene_ID ... SPOT_ID (5 total)
## fvarMetadata: Column Description labelDescription
## experimentData: use 'experimentData(object)'
## Annotation: GPL8492is(g[[1]])## [1] "ExpressionSet" "eSet" "VersionedBiobase"
## [4] "Versioned"g[[1]]## ExpressionSet (storageMode: lockedEnvironment)
## assayData: 16331 features, 4 samples
## element names: exprs
## protocolData: none
## phenoData
## sampleNames: GSM509823 GSM509824 GSM509825 GSM509826
## varLabels: title geo_accession ... data_row_count (35 total)
## varMetadata: labelDescription
## featureData
## featureNames: 100009600_at 100012_at ... 99982_at (16331 total)
## fvarLabels: ID Gene_ID ... SPOT_ID (5 total)
## fvarMetadata: Column Description labelDescription
## experimentData: use 'experimentData(object)'
## Annotation: GPL8492データは、ExpressionSet が要素にはいった、List のデータ構造で保存されます。ExpressionSet は、Bioconductor で発現データを持つ際に良く使われるデータ構造です。この構造にしておけば、あらゆる発現解析の関数が利用できるようになります。なぜ List になっているかと言うと、検索で複数の結果が帰ってくる場合を想定しているためです。
これから発現量のテーブルを取り出してみます。
my.exprs <- exprs(g[[1]])
is(my.exprs)## [1] "matrix" "array" "structure" "vector"head(my.exprs)## GSM509823 GSM509824 GSM509825 GSM509826
## 100009600_at 5.688 5.568 5.874 5.510
## 100012_at 2.503 2.348 2.488 2.467
## 100017_at 8.028 8.233 6.987 7.004
## 100019_at 7.292 7.312 7.296 7.107
## 100034251_at 4.647 4.693 5.149 5.159
## 100036521_at 8.340 8.271 8.505 8.505注意
中に入っている数値についてはデータベースエントリをしっかり読んで把握しましょう。特に以下の点に注意しましょう。
- マイクロアレイのプラットフォームは?
- 列とサンプルの関係は?
- サンプル間の正規化をしてあるか?
- 競合ハイブリダイゼーション型のマイクロアレイの場合は、比率のデータになっている場合があるが、どちらが分母か?
実行環境
sessionInfo()## R version 2.15.1 (2012-06-22)
## Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
##
## locale:
## [1] ja_JP.UTF-8/ja_JP.UTF-8/ja_JP.UTF-8/C/ja_JP.UTF-8/ja_JP.UTF-8
##
## attached base packages:
## [1] datasets utils stats graphics grDevices methods base
##
## other attached packages:
## [1] GEOquery_2.23.5 Biobase_2.16.0 BiocGenerics_0.2.0
## [4] knitr_0.7.6 stringr_0.6.1 RColorBrewer_1.0-5
## [7] MASS_7.3-20 plyr_1.7.1 BiocInstaller_1.4.7
##
## loaded via a namespace (and not attached):
## [1] digest_0.5.2 evaluate_0.4.2 formatR_0.6 RCurl_1.91-1
## [5] tools_2.15.1 XML_3.9-4