はじめに
RNA-seq解析では本当にたくさんのファイル形式を扱うことになります！
どのファイルフォーマットがどんなデータを保存していて、どんな形式で記載されていて、この列は何を意味しているのかを把握するのは本当に大変で、私も忘れるたびにいちいち調べていました。。。
この記事ではRNA-seq解析で使用するメジャーなファイル形式をまとめて紹介しています！
「あのファイル形式の中身どんなのだっけ？」とか「この列どんな意味だっけ？」ていう時に読んでいただけるととても嬉しいです！
一番最初に各ファイル形式の簡単な概要をまとめてご紹介して、後ほど各項目で詳細にご紹介します！
FASTQ: 生のRNA-seqリードを保存するための最も一般的なファイル形式です。テキストベースのフォーマットで、各リードの塩基配列とクオリティスコアを別々の行に保存する。FASTQファイルの拡張子は通常”.fastq “または”.fq “。
SAM（Sequence Alignment/Map）: 参照ゲノムに対するRNA-seqリードのアラインメントを格納するテキストベースのフォーマット。リード名、参照名、アライメントの開始位置と終了位置、アライメントスコアなどの情報が含まれる。
BAM (Binary Alignment/Map) : アライメントされたRNA-seqリードを保存するためのバイナリフォーマット。同じ情報をテキストベースで保存するSAM（Sequence Alignment/Map）フォーマットを圧縮したものです。BAMファイルは一般的にSAMファイルより小さく、読み書きが速いため、アラインメントされたRNA-seqデータの保存によく使用される。
SRA (Sequence Read Archive): BAMをベースにした圧縮ファイルで、NCBI（National Center for Biotechnology Information）が高スループットシーケンシングデータ（次世代シーケンシングプラットフォームによって生成されるデータなど）を保存するために使用するファイルフォーマットです。 生データや処理済みデータ、リード配列、品質スコア、サンプル情報などのメタデータを格納することができる。圧縮されているためファイルサイズを小さく、効率的なデータ転送や保存が可能。
GTF/GFF: GTF (Gene Transfer Format) , GFF (General Feature Format) は、遺伝子、エクソン、イントロンなどの参照ゲノム上のフィーチャーに関する情報を格納するフォーマット。参照ゲノムのアノテーション情報を格納するために使用されることが多い。
BED (Browser Extensible Data) : 参照ゲノム上の遺伝子やエクソンなどのフィーチャーの位置を格納するテキストベースの簡易フォーマットです。RNA-seq解析で同定された発現量の異なる遺伝子の位置を格納するためによく利用される。
FASTQ
FASTQ形式は生のRNA-seqリードを保存するための最も一般的なファイル形式です。テキストベースのフォーマットで、各リードの塩基配列とクオリティスコアを別々の行に保存します。FASTQファイルの拡張子は通常”.fastq “または”.fq “です。
塩基配列やアミノ酸配列を保存するFASTA形式にSequence Qualityを記述できるように拡張したものです。ちなみにFASTA は、”FAST-Aye”（ファストエー）と発音するらしいので、FASTQの読み方は”FAST-Qu”（ファストキュー）でしょうか。
FASTQファイルは容量を見て貰えばわかるように数億行から成るファイルなので、エクセルなどで見ようとするとフリーズの原因になります。terminalのheadコマンドなどで確認することをお勧めします。
FASTQファイルは4行で1セットであり、各行の記載内容は以下です
1行目:リードの識別子。「@」記号で始まる。
2行目:実際に読み取った塩基配列。
3行目: 1行目と同じリード識別子。「+」記号で始まる。
4行目:読み取った塩基のクオリティスコア。
シーケンスラン番号、レーン番号、リード番号などの情報が含まれます。
2行目:実際に読み取った塩基配列。
大文字でも小文字でもどちらで表現しても大丈夫です。
3行目: 1行目と同じリード識別子。「+」記号で始まる。
はじめの文字が「@」から「+」に変わっただけで内容は1行目と全く同じです。
この行はFASTQファイルによっては含まれていないものもあります。この行がある場合、配列と品質情報を分離するために使用され、それ以上の解析には使用されません。
4行目:読み取った塩基のクオリティスコア。
このスコアはPhred quality scoreシステムを用いてエンコードされていて、各塩基の塩基判定が不正確である確率を反映したスコアを一文字一文字が表現しています。スコアはASCII文字としてエンコードされ、スコアが高いほど高いASCII valueで表されます。
ちなみにクオリティスコアは最低の0が「!」で表現され、最大の40が「I」で表現されます。これはASCIIコード表を見ていただいた方が面白いと思いますので、こちらのサイトなどをご参照ください。（https://ja.wikipedia.org/wiki/ASCII）
クオリティスコア最大の40（「I」）だとその塩基が正しい確率（推定ベースコール精度）が99.99%で、クオリティスコア30（「?」）だと99.9％のようです。
FASTQ形式の簡単な例
@SRR1234.1
AGCTAGCTAGCTAGCTAGCTAGCTAGCTAGCTAGCT
+SRR1234.1
IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII実際のファイル例 (SRR1039508_1.fastq)
@SRR1039508.1 HWI-ST177:290:C0TECACXX:1:1101:1225:2130 length=63
CATTGCTGATACCAANNNNNNNNGCATTCCTCAAGGTCTTCCTCCTTCCCTTACGGAATTACA
+SRR1039508.1 HWI-ST177:290:C0TECACXX:1:1101:1225:2130 length=63
HJJJJJJJJJJJJJJ########00?GHIJJJJJJJIJJJJJJJJJJJJJJJJJHHHFFFFFD
@SRR1039508.2 HWI-ST177:290:C0TECACXX:1:1101:1311:2131 length=63
CCCTGGACTGCTTCTTGAAAAGTGCCATCCAAACTCTATCTTTGGGGAGAGTATGATAGAGAT
+SRR1039508.2 HWI-ST177:290:C0TECACXX:1:1101:1311:2131 length=63
HJJJJJJJJJJJJJJJJJIIJIGHIJJJJJJJJJJJJJJJJJJJJJJGHHIDHIJJHHHHHHF
Peter J. A. Cock, Christopher J. Fields, Naohisa Goto, Michael L. Heuer, Peter M. Rice, The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants, Nucleic Acids Research, Volume 38, Issue 6, 1 April 2010, Pages 1767–1771, https://doi.org/10.1093/nar/gkp1137
SAM（Sequence Alignment/Map）
SAM: SAM（Sequence Alignment/Map）は、参照ゲノムに対するRNA-seqリードのアラインメントを格納するテキストベースのフォーマットです。リード名、参照名、アライメントの開始位置と終了位置、アライメントスコアなどの情報が含まれています。
SAM (Sequence Alignment/Map) は、ハイスループットなシーケンシングリードを参照ゲノムにアライメントして保存するためのファイルフォーマットです。SAMファイルは、アライメント情報をタブ区切りで格納したテキストファイルで、テキスト編集ツールを用いて簡単に閲覧・操作することができます。
SAMファイルは、1行に1つのアライメントを含み、各行は11の必須フィールドといくつかのオプションフィールドで構成されています。
各フィールドの意味は以下です！
- QNAME (Query template NAME): クエリー・テンプレートネームと呼ばれ、リードの識別子を含む
- FLAG: アライメントのプロパティを記述するビット単位のフラグ
- RNAME (Reference sequence NAME): 参照配列名、リードがアライメントされる参照配列の名前
- POS (1-based leftmost mapping POSition): 1ベースの左端のマッピング位置、アライメントを開始する参照配列上の位置
- MAPQ (MAPping Quality): マッピングの品質、アライメントの品質をPhred-scaled probabilityで表す
- CIGAR: CIGAR文字列、アライメントのコンパクトな記述、一致、挿入、削除された塩基の数を表す
- RNEXT (Reference name of NEXT read): テンプレート内の次のリードの参照名、テンプレート内の次のセグメントの参照配列の名前
- PNEXT (Position of NEXT read): テンプレート内の次のリードの位置、テンプレート内の次のセグメントの参照配列上の位置
- TLEN (Template LENgth): 観察されたテンプレートの長さ、テンプレートの観察された長さ
- SEQ (SEQuence): リードの塩基配列
- QUAL: ASCII of Phred-scaled base quality、リードの塩基の品質スコア
それでは実際にSAMファイルの例を見ていきましょう！こちらのページ(https://genome.sph.umich.edu/wiki/SAM#Example_SAM)で用意してくれているSAMファイルの例を参考にしています。
Header Linesの例
@HD VN:1.0 SO:coordinate
@SQ SN:1 LN:249250621 AS:NCBI37 UR:file:/data/local/ref/GATK/human_g1k_v37.fasta M5:1b22b98cdeb4a9304cb5d48026a85128
@SQ SN:2 LN:243199373 AS:NCBI37 UR:file:/data/local/ref/GATK/human_g1k_v37.fasta M5:a0d9851da00400dec1098a9255ac712e
@SQ SN:3 LN:198022430 AS:NCBI37 UR:file:/data/local/ref/GATK/human_g1k_v37.fasta M5:fdfd811849cc2fadebc929bb925902e5
@RG ID:UM0098:1 PL:ILLUMINA PU:HWUSI-EAS1707-615LHAAXX-L001 LB:80 DT:2010-05-05T20:00:00-0400 SM:SD37743 CN:UMCORE
@RG ID:UM0098:2 PL:ILLUMINA PU:HWUSI-EAS1707-615LHAAXX-L002 LB:80 DT:2010-05-05T20:00:00-0400 SM:SD37743 CN:UMCORE
@PG ID:bwa VN:0.5.4
@PG ID:GATK TableRecalibration VN:1.0.3471 CL:Covariates=[ReadGroupCovariate, QualityScoreCovariate, CycleCovariate, DinucCovariate, TileCovariate], default_read_group=null, default_platform=null, force_read_group=null, force_platform=null, solid_recal_mode=SET_Q_ZERO, window_size_nqs=5, homopolymer_nback=7, exception_if_no_tile=false, ignore_nocall_colorspace=false, pQ=5, maxQ=40, smoothing=1Alignmentsの例
1:497:R:-272+13M17D24M 113 1 497 37 37M 15 100338662 0 CGGGTCTGACCTGAGGAGAACTGTGCTCCGCCTTCAG 0;==-==9;>>>>>=>>>>>>>>>>>=>>>>>>>>>> XT:A:U NM:i:0 SM:i:37 AM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:37
19:20389:F:275+18M2D19M 99 1 17644 0 37M = 17919 314 TATGACTGCTAATAATACCTACACATGTTAGAACCAT >>>>>>>>>>>>>>>>>>>><<>>><<>>4::>>:<9 RG:Z:UM0098:1 XT:A:R NM:i:0 SM:i:0 AM:i:0 X0:i:4 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:37
19:20389:F:275+18M2D19M 147 1 17919 0 18M2D19M = 17644 -314 GTAGTACCAACTGTAAGTCCTTATCTTCATACTTTGT ;44999;499<8<8<<<8<<><<<<><7<;<<<>><< XT:A:R NM:i:2 SM:i:0 AM:i:0 X0:i:4 X1:i:0 XM:i:0 XO:i:1 XG:i:2 MD:Z:18^CA19
9:21597+10M2I25M:R:-209 83 1 21678 0 8M2I27M = 21469 -244 CACCACATCACATATACCAAGCCTGGCTGTGTCTTCT <;9<<5><<<<><<<>><<><>><9>><>>>9>>><> XT:A:R NM:i:2 SM:i:0 AM:i:0 X0:i:5 X1:i:0 XM:i:0 XO:i:1 XG:i:2 MD:Z:35
この例では3つのアライメントが含まれていて、各アライメントは各フィールドに以下のような値を保持しています！
| Field | Alignment 1 | Alignment 2 | Alignment 3 | Alignment 4 |
| QNAME | 1:497:R:-272+13M17D24M | 19:20389:F:275+18M2D19M | 19:20389:F:275+18M2D19M | 9:21597+10M2I25M:R:-209 |
| FLAG | 113 | 99 | 147 | 83 |
| RNAME | 1 | 1 | 1 | 1 |
| POS | 497 | 17644 | 17919 | 21678 |
| MAPQ | 37 | 0 | 0 | 0 |
| CIGAR | 37M | 37M | 18M2D19M | 8M2I27M |
| RNEXT | 15 | = | = | = |
| PNEXT | 100338662 | 17919 | 17644 | 21469 |
| TLEN | 0 | 314 | ||
| SEQ | CGGGTCTGACCTGAGGAGAACTGTGCTCCGCCTTCAG | TATGACTGCTAATAATACCTACACATGTTAGAACCAT | GTAGTACCAACTGTAAGTCCTTATCTTCATACTTTGT | CACCACATCACATATACCAAGCCTGGCTGTGTCTTCT |
| QUAL | 0;==-==9;>>>>>=>>>>>>>>>>>=>>>>>>>>>> | >>>>>>>>>>>>>>>>>>>><<>>><<>>4::>>:<9 | ;44999;499<8<8<<<8<<><<<<><7<;<<<>><< | <;9<<5><<<<><<<>><<><>><9>><>>>9>>><> |
| TAGs | XT:A:U NM:i:0 SM:i:37 AM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:37 | RG:Z:UM0098:1 XT:A:R NM:i:0 SM:i:0 AM:i:0 X0:i:4 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:37 | XT:A:R NM:i:2 SM:i:0 AM:i:0 X0:i:4 X1:i:0 XM:i:0 XO:i:1 XG:i:2 MD:Z:18^CA19 | |
Heng Li, Bob Handsaker, Alec Wysoker, Tim Fennell, Jue Ruan, Nils Homer, Gabor Marth, Goncalo Abecasis, Richard Durbin, 1000 Genome Project Data Processing Subgroup, The Sequence Alignment/Map format and SAMtools, Bioinformatics, Volume 25, Issue 16, August 2009, Pages 2078–2079, https://doi.org/10.1093/bioinformatics/btp352
BAM (Binary Alignment/Map)
BAM: BAM (Binary Alignment/Map) は、アライメントされたRNA-seqリードを保存するためのバイナリフォーマットです。同じ情報をテキストベースで保存するSAM（Sequence Alignment/Map）フォーマットを圧縮したものなので、内容についてはこの記事の上のSAMの内容をご参照ください。BAMファイルは一般的にSAMファイルより小さく、読み書きが速いため、アラインメントされたRNA-seqデータの保存によく使用されます。
内容はSAMと同じ！
SAM（Sequence Alignment/Map）フォーマットを圧縮したもの
SAMファイルよりファイルサイズが小さく、読み書きが速い！
Heng Li, A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data, Bioinformatics, Volume 27, Issue 21, November 2011, Pages 2987–2993, https://doi.org/10.1093/bioinformatics/btr509
SRA (Sequence Read Archive)
SRAは大量のハイスループットシーケンスデータを保存・配布するために使われるファイル形式で、BAM形式をベースにしており、バイナリで圧縮されており、一般にダウンロード可能です！大規模なシーケンスデータセットの保存と共有に利用されていて、下流の解析のために他の一般的なファイル形式に簡単に変換することができます。
以下がSRAファイル形式のまとめです！
- NCBIが提供する、シーケンスリードアーカイブ（SRA）データベースのフォーマット
- BAM形式をベースにしており、ファイルサイズの縮小、データ転送や格納の効率化のために圧縮されている
- 読み取り配列、品質スコア、サンプル情報などのメタデータなどの、生データや処理済みデータを格納
- SRAファイルは、FASTQやBAMなどの他のフォーマットに変換可能
- NCBIのSRAデータベースにアップロードされ、誰にでもアクセス可能
- SRAファイルを開くには、特別なソフトウェア（SRAツールキット）が必要
SRAファイルはNational Center for Biotechnology Information (NCBI) データベースの一部であり、一般にダウンロードが可能です。NCBIが作成したSRA Toolkitのfastq-dumpコマンドで簡単にダウンロードできます！
SRA形式は、NCBIのBAM（Binary Alignment/Map） 形式に基づいていて、生のシーケンスリードに加えて、リード識別子、品質スコア、ベースコールなどの追加情報が含まれています！
SRA ファイルはバイナリファイルであり、効率的な保存と検索のために圧縮とインデックス化が行われています。
NCBIのウェブサイトからダウンロードし、SRA Toolkitなどの専用ツールを用いてFASTQ、BAM、SAMなどの一般的なファイル形式に変換することができます！
また、SRA ファイルは階層構造で構成されています！SRAファイルの階層構造は、「実験」、「ラン」、そして「リード」の3つのレベルに分かれています！各実験は複数のランを持ち、各ランは1組の生シーケンスリードに対応します。これにより、大量のシーケンスデータを効率的に管理することができます！
https://www.ncbi.nlm.nih.gov/sra/docs/
https://academic.oup.com/nar/article/39/suppl_1/D19/2505848
https://academic.oup.com/nar/article/50/D1/D387/6438001
GTF (Gene Transfer Format) と GFF (General Feature Format)
GTF (Gene Transfer Format) と GFF (General Feature Format) は、遺伝子、エクソン、イントロンなどの参照ゲノム上の特徴（Feature）に関する情報を格納するフォーマットです。参照ゲノムのアノテーション情報を格納するために使用されることが多いです！
また、GTF と GFF ファイルは似た情報を含んでいますが、フォーマット間には重要な違いがあります！GTFは特にモデル生物の参照ゲノムに使用するために設計されたフォーマットである一方、GFFはあらゆる生物に使用できるより一般的なフォーマットです！
両方のフォーマットは9つのタブ区切りカラムから構成されていて、各列は情報を含んでいます！
- 配列名。染色体の名前のことが多い。一番染色体ならchr1や1のように表記される
- アノテーションのソース
- 特徴の種類 (gene, exon, transcriptなど )
- 特徴の開始位置
- 特徴の終了位置
- スコア（e-valueや bitscore など）アノテーションの信頼性を示す。”.”はこのアノテーションにはスコアが割り当てられていないことを意味する
- ストランド(Strand)（+または-）特徴が配列上でどちらの方向に存在するかを示す。”+” は順方向、”-” は逆方向
- フレーム（コーディング領域用）タンパク質配列を構成するDNA配列の三つのフレーム（0, 1, 2）のうち、特徴がどのフレームに属しているかを示す。”.”はフレームが割り当てられていないことを意味する
- 追加の属性（遺伝子 ID またはトランスクリプト IDなど）
それでは単純化した例と実際のファイル例を見ていきましょう！
GTF/GFF形式の簡単な例
chr1 gene exon 1000 2000 . + . gene_id "gene1"; transcript_id "transcript1"; exon_number "1";この例ではchromoseme1から、「gene」というソースから得られたアノテーションです。
exonが1000〜2000の領域に存在していることを示します。
スコアは”.”なので割り当てられておらず、ストランドは「＋」なので順方向に存在し、フレームは”.”なので割り当てられていないです。
9列目の追加の情報には「gene_id "gene1"; transcript_id "transcript1"; exon_number "1"」と記載されていますが、これはgene_idが”gene1″という遺伝子で、transcript_idが”transcript1″で exon_numberが “1”ということを示しています！
9列目の追加情報は遺伝子のIDに関する情報を詰め込むのでかなり長い文字列になることが多いみたいです。
実際のファイル例 (Homo_sapiens.GRCh38.108.chromosome.1.gtf)
1 ensembl_havana gene 1471765 1497848 . + . gene_id "ENSG00000160072"; gene_version "20"; gene_name "ATAD3B"; gene_source "ensembl_havana"; gene_biotype "protein_coding";
1 ensembl_havana transcript 1471765 1497848 . + . gene_id "ENSG00000160072"; gene_version "20"; transcript_id "ENST00000673477"; transcript_version "1"; gene_name "ATAD3B"; gene_source "ensembl_havana"; gene_biotype "protein_coding"; transcript_name "ATAD3B-206"; transcript_source "ensembl_havana"; transcript_biotype "protein_coding"; tag "CCDS"; ccds_id "CCDS30"; tag "basic"; tag "Ensembl_canonical"; tag "MANE_Select";
1 ensembl_havana exon 1471765 1472089 . + . gene_id "ENSG00000160072"; gene_version "20"; transcript_id "ENST00000673477"; transcript_version "1"; exon_number "1"; gene_name "ATAD3B"; gene_source "ensembl_havana"; gene_biotype "protein_coding"; transcript_name "ATAD3B-206"; transcript_source "ensembl_havana"; transcript_biotype "protein_coding"; tag "CCDS"; ccds_id "CCDS30"; exon_id "ENSE00003889014"; exon_version "1"; tag "basic"; tag "Ensembl_canonical"; tag "MANE_Select";
1 ensembl_havana CDS 1471885 1472089 . + 0 gene_id "ENSG00000160072"; gene_version "20"; transcript_id "ENST00000673477"; transcript_version "1"; exon_number "1"; gene_name "ATAD3B"; gene_source "ensembl_havana"; gene_biotype "protein_coding"; transcript_name "ATAD3B-206"; transcript_source "ensembl_havana"; transcript_biotype "protein_coding"; tag "CCDS"; ccds_id "CCDS30"; protein_id "ENSP00000500094"; protein_version "1"; tag "basic"; tag "Ensembl_canonical"; tag "MANE_Select";
1 ensembl_havana start_codon 1471885 1471887 . + 0 gene_id "ENSG00000160072"; gene_version "20"; transcript_id "ENST00000673477"; transcript_version "1"; exon_number "1"; gene_name "ATAD3B"; gene_source "ensembl_havana"; gene_biotype "protein_coding"; transcript_name "ATAD3B-206"; transcript_source "ensembl_havana"; transcript_biotype "protein_coding"; tag "CCDS"; ccds_id "CCDS30"; tag "basic"; tag "Ensembl_canonical"; tag "MANE_Select";
https://asia.ensembl.org/info/website/upload/gff.html
http://mblab.wustl.edu/GTF22.html
https://www.ncbi.nlm.nih.gov/genbank/genomes_gff/
BED (Browser Extensible Data)（追記予定）
BED BED (Browser Extensible Data) は、参照ゲノム上の遺伝子やエクソンなどのフィーチャーの位置を格納するテキストベースの簡易フォーマットです。RNA-seq解析で同定された発現量の異なる遺伝子の位置を格納するためによく利用されます！
参考
Peter J. A. Cock, Christopher J. Fields, Naohisa Goto, Michael L. Heuer, Peter M. Rice, The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants, Nucleic Acids Research, Volume 38, Issue 6, 1 April 2010, Pages 1767–1771, https://doi.org/10.1093/nar/gkp1137
SAM参考
Heng Li, Bob Handsaker, Alec Wysoker, Tim Fennell, Jue Ruan, Nils Homer, Gabor Marth, Goncalo Abecasis, Richard Durbin, 1000 Genome Project Data Processing Subgroup, The Sequence Alignment/Map format and SAMtools, Bioinformatics, Volume 25, Issue 16, August 2009, Pages 2078–2079, https://doi.org/10.1093/bioinformatics/btp352
BAM参考
Heng Li, A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data, Bioinformatics, Volume 27, Issue 21, November 2011, Pages 2987–2993, https://doi.org/10.1093/bioinformatics/btr509
SRA参考
https://www.ncbi.nlm.nih.gov/sra/docs/
https://academic.oup.com/nar/article/39/suppl_1/D19/2505848
https://academic.oup.com/nar/article/50/D1/D387/6438001
GTF/GFF参考
https://asia.ensembl.org/info/website/upload/gff.html
http://mblab.wustl.edu/GTF22.html
https://www.ncbi.nlm.nih.gov/genbank/genomes_gff/