--------[p.46] -----------------------------------------
The GO annotation of genes
GO= Gene Ontology
http://geneontology.org/ 遺伝子オントロジーのプロジェクトサイト
遺伝子機能情報を生物種横断的に利用する為の、統一的な語彙オントロジをプロジェクトで規定していく。
UniProt(EBI), Mouse Genome Informatics(Jackson Lab.), Yeast Genome Database(Stanford Univ)などの老舗のデータベースサイトメンバがプロジェクトを推進している。
オントロジー＝有向グラフにより上位・下位概念を規定した語彙リスト
--------------------------------------------------------------------
＜GOの特徴1: 3区分 ＞
Gene Ontology(GO)は、次の3区分に分けられる。
1)Biological process=(BP)
2)Molecular function=(MF)
3)Cellular component=(CC)
---------------------------------------------------------------------
＜GOの特徴2: 注釈付けの「型」をエビデンスコードで表す＞
http://geneontology.org/page/guide-go-evidence-codes
Inferred from Electronic Annotation (IEA) ：機械注釈推定
Inferred from Experiment (EXP)： 実験による推定
Inferred from Sequence or structural Similarity (ISS)：配列の類似性による推定
------------------------------------------------------------------------
＜GOの特徴3: モデル生物の注釈データ＞
http://geneontology.org/page/download-annotations
===============================================
ここではorg.Hs.eg.dbパッケージを利用して遺伝子のリストへのGOの注釈付け（Annotationと呼ぶ）を行っていく。
＜org.Hs.eg.dbのヘルプ＞
http://bioconductor.org/packages/release/data/annotation/manuals/org.Hs.eg.db/man/org.Hs.eg.db.pdf
> library(org.Hs.eg.db)
要求されたパッケージ AnnotationDbi をロード中です
要求されたパッケージ stats4 をロード中です
要求されたパッケージ BiocGenerics をロード中です
要求されたパッケージ parallel をロード中です
略
要求されたパッケージ IRanges をロード中です
要求されたパッケージ DBI をロード中です
> myEIDs <- c("1","10","100")
遺伝子番号としてEntrez ID番号を設定する
> myEIDs
[1] "1" "10" "100"
> myGO <- org.Hs.egGO[[myEIDs[1]]]
番号1だけGOに変換、
①GO番号,②オントロジー,③エビデンスコードのセットを返す
> myGO 中身を確認
$`GO:0008150`
$`GO:0008150`$GOID
[1] "GO:0008150" ＜ ①GO番号
$`GO:0008150`$Evidence ＜②エビデンスコード
[1] "ND"
$`GO:0008150`$Ontology ＜③オントロジー
[1] "BP"
$`GO:0005576`
$`GO:0005576`$GOID
[1] "GO:0005576"
$`GO:0005576`$Evidence
[1] "IDA"
$`GO:0005576`$Ontology
[1] "CC"
$`GO:0005615`
$`GO:0005615`$GOID
[1] "GO:0005615"
$`GO:0005615`$Evidence
[1] "IDA"
$`GO:0005615`$Ontology
[1] "CC"
$`GO:0070062`
$`GO:0070062`$GOID
[1] "GO:0070062"
$`GO:0070062`$Evidence
[1] "IDA"
$`GO:0070062`$Ontology
[1] "CC"
$`GO:0072562`
$`GO:0072562`$GOID
[1] "GO:0072562"
$`GO:0072562`$Evidence
[1] "IDA"
$`GO:0072562`$Ontology
[1] "CC"
$`GO:0003674`
$`GO:0003674`$GOID
[1] "GO:0003674"
$`GO:0003674`$Evidence
[1] "ND"
$`GO:0003674`$Ontology
[1] "MF"
> myGO_ALL <- mget(myEIDs,org.Hs.egGO)
mget関数を使えば、一度にGO変換可能
> summary(myGO_ALL)
Length Class Mode
1 6 -none- list
10 5 -none- list
100 64 -none- list
summaryで何個のGO注釈を持つか個数を確認できる
> head(myGO_ALL)
$`1`
$`1`$`GO:0008150`
$`1`$`GO:0008150`$GOID
[1] "GO:0008150"
略
> GOgenes <- org.Hs.egGO2ALLEGS[["GO:0008150"]]
逆にGO番号からEntrez ID(遺伝子番号)に変換
> summary(GOgenes)
Length Class Mode
39588 character character
> head(GOgenes)
ND IDA IEA TAS TAS TAS
"1" "2" "2" "2" "9" "10"
39588個のデータがエビデンスコードと共に表示される
＜GO:0008150はヒット数多なので1909でテスト＞
> GOgenes2 <- org.Hs.egGO2ALLEGS[["GO:0001909"]]
> summary(GOgenes2)
Length Class Mode
90 character character
> GOgenes2
IMP IMP IEA IEA ISS
"164" "409" "563" "567" "1054"
IEA IMP IEA ISS IEA
"1075" "1130" "1511" "1512" "1991"
ISS IDA IMP IEA IDA
"2013" "2150" "3002" "3077" "3105"
IEA IDA IEA IEA IDA
"3105" "3106" "3106" "3107" "3133"
IEA IEA IEA IEA TAS
"3134" "3135" "3136" "3140" "3383"
IEA IDA IDA ISS ISS
"3575" "3592" "3593" "3593" "3594"
IMP IEA IMP IEA IBA
"3811" "3902" "3916" "4068" "4277"
TAS IEA IEA IMP IEA
"4843" "5027" "5052" "5272" "5532"
IEA ISS IDA IMP IDA
"5777" "5788" "5817" "5817" "5819"
IMP IEA IMP IEA IDA
"5819" "5873" "6318" "6347" "6375"
IEA IEA IMP IMP IMP
"6776" "6777" "6813" "6845" "8417"
IEA ISS TAS NAS IMP
"8676" "8772" "9437" "10383" "10666"
IDA ISS IEA IDA IDA
"10859" "11035" "11334" "22914" "23705"
ISS ISS IEA IDA NAS
"23705" "51561" "54210" "56253" "57823"
IDA IDA IDA IDA IEA
"59067" "79465" "80328" "80329" "127544"
IBA IDA IDA ISS IBA
"135250" "135250" "146850" "149233" "154064"
IEA NAS IMP IDA TAS
"201294" "203068" "259197" "353091" "387836"
IEA TAS IDA IMP IDA
"653361" "100431172" "100507436" "100507436" "100528032"
> GOterm <- c("GO:0008150","GO:0001909")
> GOgenes_ALL <- mget(GOterm, org.Hs.egGO2ALLEGS)
mgetを使えば、一括で処理可能
======[p.48]=========================================
GO.dbパッケージには、GOのTermの番号・親子関係が格納されている。
> biocLite(c("AnnotationDbi", "GO.db"), type="source") ＜注意！本と異なります＞
BioC_mirror: http://bioconductor.org
Using Bioconductor version 3.1 (BiocInstaller 1.18.4), R version 3.2.2.
Installing package(s) ‘AnnotationDbi’, ‘GO.db’
trying URL 'http://bioconductor.org/packages/3.1/bioc/src/contrib/AnnotationDbi_1.30.1.tar.gz'
Content type 'application/x-gzip' length 4266883 bytes (4.1 MB)
downloaded 4.1 MB
trying URL 'http://bioconductor.org/packages/3.1/data/annotation/src/contrib/GO.db_3.1.2.tar.gz'
Content type 'application/x-gzip' length 30192432 bytes (28.8 MB)
downloaded 28.8 MB
*
Creating a generic function for 'nchar' from package 'base' in package 'S4Vectors'
* DONE (GO.db)
Update all/some/none? [a/s/n]:
n
> library(GO.db)
＜以下本と異なります＞
> tmp<-GOID("GO:0001909")
> tmp
GO:0001909
"GO:0001909"
>
> tmp_term <- Term(tmp)
> tmp_term
GO:0001909
"leukocyte mediated cytotoxicity"
>
> tmp_ont <- Ontology(tmp)
> tmp_ont
GO:0001909
"BP"
> tmp_syn <- Synonym(tmp)
> tmp_syn
$`GO:0001909`
[1] "immune cell mediated cell death" "immune cell mediated cell killing"
[3] "immune cell mediated cytotoxicity" "leucocyte mediated cytotoxicity"
> tmp_sec <- Secondary(tmp)
> tmp_sec
$`GO:0001909`
[1] NA NA NA NA
> tmp_def <- Definition(tmp)
> tmp_def
GO:0001909
"The directed killing of a target cell by a leukocyte."
========[p.48]=======================
The GO enrichment of genes
Enrichment analysis(またはHypergeometric test)で、発現遺伝子のセットがどのGOにあてはまるかを統計的にテストします。
＜Enrichment解析の説明＞
まずゲノムワイドな遺伝子セットについてGOのTerm別で分布をカウントしておきます。知りたいクエリの遺伝子セット（例：Microarray解析でのUp-regulated genesセット）について、GO Termの分布をカウントします。GOTerm毎に、ゲノムワイドな遺伝子セット(background dataset)と比較して有意に多いかP値で表します。分布の検定には、超幾何分布検定(hypergeometric test)もしくはFisherの正確確率検定(Fisher's exact test)が使われます。