GTF (gene feature format) ファイルは遺伝子アノテーションを記述しているファイルである。ファイルの 1 行は 1 つの feature を表す。ここでいう feature は、exon、intron、tRNA などのゲノム上の領域を指す。
各行には 9 項目のデータが記述されている。各項目は「タブ」によって区切られている。各項目は以下の順で記述される。9 項目をすべて記述することになっている。不明な場合はピリオド(.)を記入する。
- seqname - 染色体名
- source- データの種類またはプログラム名など
- feature - intron、exon、CDS、RNA、miRNA など
- start - feature の開始位置
- end - feature の終了位置
- score - ?
- strand - 鎖の方向( + または - )
- frame - コドンの開始位置( 0、1 または 2 )
- attribute - 属性など。Ensembl ID、遺伝子汎用名、エキソンの番号などが記載されている
実際の GTP ファイルをテキストエディタで開くと、以下のようになっている。
GL000213.1 protein_coding CDS 138767 139287 . - 0 gene_id "ENSG00000237375"; transcript_id "ENST00000327822"; exon_number "1"; gene_name "BX072566.1"; gene_biotype "protein_coding"; transcript_name "BX072566.1-201"; protein_id "ENSP00000329990"; GL000213.1 protein_coding start_codon 139285 139287 . - 0 gene_id "ENSG00000237375"; transcript_id "ENST00000327822"; exon_number "1"; gene_name "BX072566.1"; gene_biotype "protein_coding"; transcript_name "BX072566.1-201"; GL000213.1 protein_coding exon 134276 134390 . - . gene_id "ENSG00000237375"; transcript_id "ENST00000327822"; exon_number "2"; gene_name "BX072566.1"; gene_biotype "protein_coding"; transcript_name "BX072566.1-201"; GL000213.1 protein_coding CDS 134276 134390 . - 1 gene_id "ENSG00000237375"; transcript_id "ENST00000327822"; exon_number "2"; gene_name "BX072566.1"; gene_biotype "protein_coding"; transcript_name "BX072566.1-201"; protein_id "ENSP00000329990"; GL000213.1 protein_coding exon 133943 134116 . - . gene_id "ENSG00000237375"; transcript_id "ENST00000327822"; exon_number "3"; gene_name "BX072566.1"; gene_biotype "protein_coding"; transcript_name "BX072566.1-201";
例えば 1 行目を見ると、これ ID が GL000213.1 染色体上のアノテーションであることがわかる。そして、この染色体上の 138767 番目の塩基から 139287 番目の塩基までがタンパク質をコーディングする領域(protein_coding)で、実際に翻訳される部分(CDS)である。また、この部分は遺伝子 ENSG00000237375 の一部である。
awk を利用した GTF ファイルの処理
GTF ファイルの記述規則はタブ区切りで非常に明快である。Perl や Python などを利用することで簡単に処理できる。ここでは、awk とよばれるスクリプトを利用して、GTP ファイルから必要な情報を取り出す方法を示す。利用する GTF ファイルは Ensembl FTP サイトからダウンロードできる。ダウンロードしたファイルを展開して利用する。
Ensembl ID と遺伝子名の対応リストの作成
awk を利用して行全体に対してマッチングを行い、後方参照により遺伝子 ID と名前を出力する。
gawk '{print gensub(/^.+gene_id "([^"]+)".+gene_name "([^"]+)".+$/,"\\1\t\\2", "")}' Homo_sapiens.GRCh37.68.gtf
awk スクリプトの実行結果。
ENSG00000265283 MIR3118-5 ENSG00000237375 BX072566.1 ENSG00000237375 BX072566.1 ENSG00000237375 BX072566.1
重複する部分があるため、awk の連想配列を利用して重複を削除する。
gawk '
{print gensub(/^.+gene_id "([^"]+)".+gene_name "([^"]+)".+$/,"\\1\t\\2", "")}
' Homo_sapiens.GRCh37.68.gtf |
awk '
BEGIN{FS="\t"}
{if (a[$0]++ == 0) print $0}
'
awk の実行結果。
ENSG00000265283 MIR3118-5 ENSG00000237375 BX072566.1 ENSG00000238432 U6 ENSG00000262510 C21orf74
エキソンの情報を抽出
GTF ファイルから遺伝子名、エキソンの番号と位置を抽出する例。
gawk '
BEGIN{FS="\t";OFS="\t"}
{print gensub(/^.+exon_number "([^"]+)".+gene_name "([^"]+)".+$/,"\\2\t\\1", ""), $4, $5}
' Homo_sapiens.GRCh37.68.gtf
awk の実行結果。遺伝子名、エキソンの番号、開始位置、終了位置の順で出力される。
MIR3118-5 1 104742 104817 BX072566.1 1 138767 139339 BX072566.1 1 138767 139287 BX072566.1 1 139285 139287 BX072566.1 2 134276 134390 BX072566.1 2 134276 134390
RNA のアノンテーション情報を抽出
ヒトゲノムの GTF ファイルから RNA に関するアノンテーション情報だけを抽出する例。
awk 'BEGIN{FS="\t"}$2 ~ /RNA/{print $0}' Homo_sapiens.GRCh37.68.gtf
さらに対象を絞って miRNA のアノンテーション情報だけを取り出す場合は次のようにする。
awk 'BEGIN{FS="\t"}$2 ~ /^miRNA$/{print $0}' Homo_sapiens.GRCh37.68.gtf
結果。
GL000213.1 miRNA exon 104742 104817 . + . gene_id "ENSG00000265283"; transcript_id "ENST00000578976"; exon_number "1"; gene_name "MIR3118-5"; gene_biotype "miRNA"; transcript_name "MIR3118-5-201"; HSCHR21_4_CTG1_1 miRNA exon 34836774 34836866 . + . gene_id "ENSG00000266367"; transcript_id "ENST00000585000"; exon_number "1"; gene_name "AC207981.1"; gene_biotype "miRNA"; transcript_name "AC207981.1-201"; HSCHR1_2_CTG31 miRNA exon 155221944 155222032 . - . gene_id "ENSG00000262403"; transcript_id "ENST00000572520"; exon_number "1"; gene_name "AL713999.2"; gene_biotype "miRNA"; transcript_name "AL713999.2-201"; HSCHR6_MHC_QBL miRNA exon 28743697 28743804 . + . gene_id "ENSG00000265877"; transcript_id "ENST00000580289"; exon_number "1"; gene_name "AL662890.13"; gene_biotype "miRNA"; transcript_name "AL662890.13-201";