LoginSignup
5
2

More than 5 years have passed since last update.

GEOから遺伝子発現データを取ってくる【R】

Last updated at Posted at 2018-06-20

GEOとは

https://www.ncbi.nlm.nih.gov/geo/
遺伝子発現データ(マイクロアレイ,次世代シークエンス)のレポジトリです
transcriptomeだけではなく, メチル化アレイのデータもあります

データの取得

Rを使ってダウンロードします

データの概要

GSE26174
ラットの網膜に電気刺激(transcorneal electrical stimulation, TES)を与えた時の効果を調べたもの

Willmann G, Schäferhoff K, Fischer MD, Arango-Gonzalez B et al. Gene expression profiling of the retina after transcorneal electrical stimulation in wild-type Brown Norway rats. Invest Ophthalmol Vis Sci 2011 Sep 29;52(10):7529-37. PMID: 21873684

処理群4サンプル, 対照群4サンプル 合計8サンプル

コード

ライブラリの読み込み

library(Biobase)
library(GEOquery)
library(data.table)

GEOのウェブサイトからデータを取得します
getGEO()の引数 GSEMatrix = F とするとサンプルごとの数値データ, GSEMatrix = T とするとデータセットの概要(GEOのaccessionのページではseries matrix file)を取得します

gsem <- getGEO("GSE26174",GSEMatrix = T)
gse <- getGEO("GSE26174",GSEMatrix = F) 

サンプル名とか実験条件の取得

gsemはlistオブジェクトなので, listの1番目を取り出します

x <- pData(phenoData(gsem[[1]]))

得られた表はこんな感じ

> names(x)
 [1] "title"                   "geo_accession"           "status"                  "submission_date"        
 [5] "last_update_date"        "type"                    "channel_count"           "source_name_ch1"        
 [9] "organism_ch1"            "characteristics_ch1"     "characteristics_ch1.1"   "characteristics_ch1.2"  
[13] "molecule_ch1"            "extract_protocol_ch1"    "label_ch1"               "label_protocol_ch1"     
[17] "taxid_ch1"               "hyb_protocol"            "scan_protocol"           "data_processing"        
[21] "platform_id"             "contact_name"            "contact_email"           "contact_institute"      
[25] "contact_address"         "contact_city"            "contact_zip/postal_code" "contact_country"        
[29] "supplementary_file"      "supplementary_file.1"    "data_row_count"          "age:ch1"                
[33] "gender:ch1"              "strain:ch1"             

> x[1:3]
                     title geo_accession                status
GSM642729 TES Control 1 4h     GSM642729 Public on Nov 21 2011
GSM642730 TES Control 2 4h     GSM642730 Public on Nov 21 2011
GSM642801 TES Control 3 4h     GSM642801 Public on Nov 21 2011
GSM642802 TES Control 4 4h     GSM642802 Public on Nov 21 2011
GSM642803         TES 1 4h     GSM642803 Public on Nov 21 2011
GSM642804         TES 2 4h     GSM642804 Public on Nov 21 2011
GSM642805         TES 3 4h     GSM642805 Public on Nov 21 2011
GSM642806         TES 4 4h     GSM642806 Public on Nov 21 2011

ここで必要になるのは各サンプルのGEOaccessionと実験条件です
このデータセットでは"title"列に実験条件が書かれていますが, それぞれのデータセットによって表記がいろいろなので, 必ずseries matrixを確認します

数値データの取得

gseは各サンプルのデータが入ったGSEインスタンスです
gse@gsemにマイクロアレイの数値データ(プロセシング済み)が入っています
Table()でプローブIDと対応する数値データを取り出し, サンプルごとにforループで読んでいきます

y <- gse@gsms

dat <- Table(y[[1]])
dat <- dat[,1:2]

for (i in 2:length(y)){
  sub_dat <- Table(y[[i]])
  dat <- cbind(dat,sub_dat[,2])
}
dat <- data.frame(dat,row.names = 1)
colnames(dat) <- x$geo_accession

サンプルごとの数値データをdata.frameとして得ることができました

> dat[1:3,]
         GSM642729 GSM642730 GSM642801 GSM642802 GSM642803 GSM642804 GSM642805 GSM642806
10700001  10.34770  10.22240  10.23050  10.37600  10.55910  10.29450  10.30840  10.16820
10700002   2.77305   2.85680   2.80159   2.64945   2.72091   2.77305   2.72123   2.77305
10700003   8.86587   8.91388   8.85737   8.87960   9.09890   8.91388   8.85572   8.82685

保存と次回の読み込み

Rdata形式にしておくと読み込みがはやいのでsave()で保存します

save(dat,file = "GSE26714.Rdata")

次回からはload()で読み込むとgetGEO()から戻らなくてよくなります(しかも速い)

load(file = "GSE26714.Rdata")
5
2
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
5
2