【R】遺伝子IDおよび遺伝子シンボルの相互変換【biomaRt】

1. biomaRtについて

biomaRtはRのライブラリで、

変換元の遺伝子IDの名称(ensembl_gene_idなど)
変換後の遺伝子IDの名称(uniprot_gn_symbolなど)
変換したい対象のIDの一覧(ENSP00000000233など)

を指定することで、①のフォーマットで記載された③のIDを②のフォーマットに変換するテーブルを出力できます。

2. biomaRtのインストール

まずは以下のようにRのコンソールを開き、BiocManagerを介してbiomaRtをインストールします。

if (!require("BiocManager", quietly = TRUE))

    install.packages("BiocManager")

BiocManager::install("biomaRt")

3. データベースの選択

biomaRtでは最初に変換に利用するデータベースをuseEnsemblで選択します。

本記事の目的の範囲ならば、汎用性が高いgenesを選択しておけば問題ありません。

方法は以下の通りです。

library(biomaRt)

ensembl <- useEnsembl(biomart = "genes")

4. データセットの選択

次にデータセットの選択を行います。biomaRtにおける「データセット」とは「種」に相当するものです。つまり、「データセットの選択」とは「どの種を用いるか選択する」プロセスのことです。

前節で指定したデータベースで解析したい種がどのような名称で用いられているか分からないので、まずは解析したい種に対応するデータセット名をserachDatasetsで検索します。その際、目的の種の学名の一部をpatternの引数に指定すると検索が易しくなります。

例えば、マウスのデータセット名を検索する場合は以下のようになります。

searchDatasets(mart = ensembl, pattern = "mus")

                      dataset                                   description

18     bmusculus_gene_ensembl                Blue whale genes (mBalMus1.v2)

104 mmoschiferus_gene_ensembl Siberian musk deer genes (MosMos_v2_BIUU_UCD)

108    mmusculus_gene_ensembl                          Mouse genes (GRCm39)

158       psimus_gene_ensembl       Greater bamboo lemur genes (Prosim_1.0)

182     smaximus_gene_ensembl                   Turbot genes (ASM1334776v1)

209   umaritimus_gene_ensembl                 Polar bear genes (UrsMar_1.0)

               version

18         mBalMus1.v2

104 MosMos_v2_BIUU_UCD

108             GRCm39

158         Prosim_1.0

182       ASM1334776v1

209         UrsMar_1.0

ヒトのデータセット名を検索する場合は以下のようになります。

searchDatasets(mart = ensembl, pattern = "sapiens")

                 dataset              description    version

80 hsapiens_gene_ensembl Human genes (GRCh38.p14) GRCh38.p14

データセット名が判明したら、useDatasetで以下のようにデータセットを指定します。

ensembl <- useDataset(dataset = "mmusculus_gene_ensembl", mart = ensembl) # マウス

ensembl <- useDataset(dataset = "hsapiens_gene_ensembl", mart = ensembl) # ヒト

5. IDの名称を検索

続いて、対応付けたいIDの名称を調べます。

前節で選択したデータセットにおいて、変換前後の遺伝子IDがどのような名称で登録されているかsearchFiltersで調べます。

例えば、

Ensembl Gene ID
Ensembl Protein ID
NCBI refseq ID
UniProt ID
遺伝子シンボル

の対応表を作成をする場合、以下のように進めます。

5.1. Ensembl ID

1と2のEnsembl ID (Ensembl Gene ID (ENSMUSG00000000001など) / Protein ID (ENSP00000000233など))の名称は以下のように調べます。

# マウス

ensembl <- useDataset(dataset = "mmusculus_gene_ensembl", mart = ensembl)

searchFilters(mart = ensembl, pattern = "ensembl.*id")

                            name

43               ensembl_gene_id

44       ensembl_gene_id_version

45         ensembl_transcript_id

46 ensembl_transcript_id_version

47            ensembl_peptide_id

48    ensembl_peptide_id_version

49               ensembl_exon_id

                                                        description

43                      Gene stable ID(s) [e.g. ENSMUSG00000000001]

44       Gene stable ID(s) with version [e.g. ENSMUSG00000000001.5]

45                Transcript stable ID(s) [e.g. ENSMUST00000000001]

46 Transcript stable ID(s) with version [e.g. ENSMUST00000000001.5]

47                   Protein stable ID(s) [e.g. ENSMUSP00000000001]

48    Protein stable ID(s) with version [e.g. ENSMUSP00000000001.5]

49                             Exon ID(s) [e.g. ENSMUSE00000097910]

# ヒト

ensembl <- useDataset(dataset = "hsapiens_gene_ensembl", mart = ensembl)

searchFilters(mart = ensembl, pattern = "ensembl.*id")

                            name

53               ensembl_gene_id

54       ensembl_gene_id_version

55         ensembl_transcript_id

56 ensembl_transcript_id_version

57            ensembl_peptide_id

58    ensembl_peptide_id_version

59               ensembl_exon_id

                                                      description

53                       Gene stable ID(s) [e.g. ENSG00000000003]

54       Gene stable ID(s) with version [e.g. ENSG00000000003.16]

55                 Transcript stable ID(s) [e.g. ENST00000000233]

56 Transcript stable ID(s) with version [e.g. ENST00000000233.10]

57                    Protein stable ID(s) [e.g. ENSP00000000233]

58     Protein stable ID(s) with version [e.g. ENSP00000000233.5]

59                              Exon ID(s) [e.g. ENSE00000000001]

descriptionにIDの例が掲載されているので、変換したい対象のIDとパターンがマッチするか目視で確認して名称を同定します。

マウスを例にすると、ENSMUSG00000000001のようなIDと最も近いフォーマットは、descriptionの43行目のensembl_gene_idであるため、Ensembl Gene ID の名称はensembl_gene_idであると分かります。また、ENSP00000000233のようなIDの名称も同様に比較することでEnsembl Protein IDの名称はensembl_peptide_idとなることが分かります。

5.2. NCBI refseq ID

NCBI refseq IDについては以下のように検索できます。

# マウス

ensembl <- useDataset(dataset = "mmusculus_gene_ensembl", mart = ensembl)

searchFilters(mart = ensembl, pattern = "refseq")

                            name

27              with_refseq_mrna

28    with_refseq_mrna_predicted

29             with_refseq_ncrna

30   with_refseq_ncrna_predicted

31           with_refseq_peptide

32 with_refseq_peptide_predicted

76                   refseq_mrna

77         refseq_mrna_predicted

78                  refseq_ncrna

79        refseq_ncrna_predicted

80                refseq_peptide

81      refseq_peptide_predicted

                                          description

27                             With RefSeq mRNA ID(s)

28                   With RefSeq mRNA predicted ID(s)

29                            With RefSeq ncRNA ID(s)

30                  With RefSeq ncRNA predicted ID(s)

31                          With RefSeq peptide ID(s)

32                With RefSeq peptide predicted ID(s)

76              RefSeq mRNA ID(s) [e.g. NM_001001130]

77    RefSeq mRNA predicted ID(s) [e.g. XM_001002281]

78                RefSeq ncRNA ID(s) [e.g. NR_000002]

79      RefSeq ncRNA predicted ID(s) [e.g. XR_001538]

80           RefSeq peptide ID(s) [e.g. NP_001001130]

81 RefSeq peptide predicted ID(s) [e.g. XP_001004117]

# ヒト

ensembl <- useDataset(dataset = "hsapiens_gene_ensembl", mart = ensembl)

searchFilters(mart = ensembl, pattern = "refseq")

                             name

36               with_refseq_mrna

37     with_refseq_mrna_predicted

38              with_refseq_ncrna

39    with_refseq_ncrna_predicted

40            with_refseq_peptide

41  with_refseq_peptide_predicted

95                    refseq_mrna

96          refseq_mrna_predicted

97                   refseq_ncrna

98         refseq_ncrna_predicted

99                 refseq_peptide

100      refseq_peptide_predicted

                                           description

36                              With RefSeq mRNA ID(s)

37                    With RefSeq mRNA predicted ID(s)

38                             With RefSeq ncRNA ID(s)

39                   With RefSeq ncRNA predicted ID(s)

40                           With RefSeq peptide ID(s)

41                 With RefSeq peptide predicted ID(s)

95                  RefSeq mRNA ID(s) [e.g. NM_000014]

96     RefSeq mRNA predicted ID(s) [e.g. XM_003403597]

97                 RefSeq ncRNA ID(s) [e.g. NR_000005]

98    RefSeq ncRNA predicted ID(s) [e.g. XR_001736914]

99               RefSeq peptide ID(s) [e.g. NP_000005]

100 RefSeq peptide predicted ID(s) [e.g. XP_003403645]

NCBI refseq IDの名称はrefseq_peptideであることが分かります。

5.3. UniProt ID

UniProt IDは以下のようにして検索します。

# ヒト

ensembl <- useDataset(dataset = "hsapiens_gene_ensembl", mart = ensembl)

searchFilters(mart = ensembl, pattern = "uniprot")

                     name                                  description

48        with_uniprot_gn               With UniProtKB Gene Name ID(s)

49   with_uniprot_isoform                 With UniProtKB isoform ID(s)

50  with_uniprotswissprot              With UniProtKB/Swiss-Prot ID(s)

51   with_uniprotsptrembl                  With UniProtKB/TrEMBL ID(s)

107         uniprot_gn_id  UniProtKB Gene Name ID(s) [e.g. A0A023T6R1]

108     uniprot_gn_symbol   UniProtKB Gene Name symbol(s) [e.g. 5HT1A]

109       uniprot_isoform  UniProtKB isoform ID(s) [e.g. A0A096LP49-1]

110      uniprotswissprot UniProtKB/Swiss-Prot ID(s) [e.g. A0A024R1R8]

111       uniprotsptrembl     UniProtKB/TrEMBL ID(s) [e.g. A0A023T6R1]

UniProt IDの名称はuniprot_gn_idであることが分かります。

5.4. 遺伝子シンボル

遺伝子シンボル表記については以下のように検索します。

# マウス

# マウス・ヒト共通

searchFilters(mart = ensembl, pattern = "symbol")

                name                                        description

62       hgnc_symbol                      HGNC symbol(s) [e.g. OR5BS1P]

65        mgi_symbol                 MGI symbol(s) [e.g. 0610005C13Rik]

87 uniprot_gn_symbol UniProtKB Gene Name symbol(s) [e.g. 0610009B22Rik]

uniprot_gn_symbolが名称です。

6. 変換のテスト

一度に全ての対応表を作成する前に、まずは1件について対応表を作成してみます。

対応表の作成はgetBMで実施します。以下のようにattributesに前節で調べたID名称をベクトルで指定します。filtersには変換元のID名称を、valuesには変換対象のIDを指定します。

# ヒト

getBM(

  attributes=c("ensembl_peptide_id", "ensembl_gene_id", "refseq_peptide", "uniprot_gn_id", "uniprot_gn_symbol"),

  filters = "ensembl_peptide_id",

  values = c("ENSP00000000233"),

  mart= ensembl

 ensembl_peptide_id ensembl_gene_id refseq_peptide uniprot_gn_id

1    ENSP00000000233 ENSG00000004059      NP_001653        A4D0Z3

2    ENSP00000000233 ENSG00000004059      NP_001653        P84085

3    ENSP00000000233 ENSG00000004059      NP_001653        C9J1Z8

  uniprot_gn_symbol

1              ARF5

2              ARF5

3              ARF5

指定した遺伝子ID/シンボルの対応が得られました。

7. IDの変換

予め変化したい遺伝子リストを作成しておき、例えばnamelist.txtのようにして保存しておきます。

namelist.txtの内容は、以下のような1行目がヘッダで2行目以降に遺伝子IDリストが掲載されている形式にしておきます。

ensembl_id

ENSP00000000233

ENSP00000000234

ENSP00000000235

...

次に以下のようにしてnamelist.txtを読み込んで、遺伝子IDの部分だけ選択します。

d <- read.csv("namelist.txt")

v <- d[, 1]

続いて、以下のようにして変換を実行します。

res <- getBM(

 	attributes=c("ensembl_peptide_id", "ensembl_gene_id", "refseq_peptide", "uniprot_gn_id", "uniprot_gn_symbol"),

 	filters = "ensembl_peptide_id",

  	values = v,

  	mart= ensembl

最後に以下のようにして保存をします。

write.csv(res, "id_conversion/ensemblP_to_symbols.csv", quote=FALSE, row.names=FALSE)

8. ここまでのまとめ

上記を1つにまとめたスクリプトは以下です。

8.1. マウス

# install

if (!require("BiocManager", quietly = TRUE))

    install.packages("BiocManager")

BiocManager::install("biomaRt")

# select database and datasets

ensembl <- useEnsembl(biomart = "genes")

ensembl <- useDataset(dataset = "mmusculus_gene_ensembl", mart = ensembl)

# id conversion

d <- read.csv("namelist.txt")

v <- d[, 1]

res <- getBM(

 	attributes=c("ensembl_peptide_id", "ensembl_gene_id", "refseq_peptide", "uniprot_gn_id", "uniprot_gn_symbol"),

 	filters = "ensembl_peptide_id",

  	values = v,

  	mart= ensembl

# save

write.csv(res, "id_conversion/ensemblP_to_symbols.csv", quote=FALSE, row.names=FALSE)

8.2. ヒト

# install

if (!require("BiocManager", quietly = TRUE))

    install.packages("BiocManager")

BiocManager::install("biomaRt")

# select database and datasets

ensembl <- useEnsembl(biomart = "genes")

ensembl <- useDataset(dataset = "hsapiens_gene_ensembl", mart = ensembl)

# id conversion

d <- read.csv("namelist.txt")

v <- d[, 1]

res <- getBM(

 	attributes=c("ensembl_peptide_id", "ensembl_gene_id", "refseq_peptide", "uniprot_gn_id", "uniprot_gn_symbol"),

 	filters = "ensembl_peptide_id",

  	values = v,

  	mart= ensembl

# save

write.csv(res, "id_conversion/ensemblP_to_symbols.csv", quote=FALSE, row.names=FALSE)

9. Pythonから結果を読み込む

以下のようにしてRのbiomaRtで作成した遺伝子ID/シンボルの対応表をPythonから読み込んでみます。

import pandas as pd

df_eid_sym = pd.read_csv('id_conversion/ensemblP_to_symbols.csv', index_col=0)

df_eid_sym.head()

ensembl_peptide_id	uniprot_gn_id	uniprot_gn_symbol
ENSP00000354687	U5Z754	ND1
ENSP00000354687	P03886	MT-ND1
ENSP00000354499	U5YWV7	COX1
ENSP00000354499	P00395	MT-CO1
ENSP00000354665	U5Z977

以下のようにして辞書にすれば、Ensembl Protein IDから遺伝子シンボルを取得することができるようになります。

</figcaption> </figure>

eid_to_sym = df_eid_sym['uniprot_gn_symbol'].to_dict()

10. まとめ

本記事では、RのbiomaRtを用いて遺伝子ID/シンボルの変換テーブルを作成する手順を示しました。

以上になります。

目次