Popularne formaty danych biologicznych

Magda Mielczarek
Magda Mielczarek 29 wrz, 10 minut czytania

Istnieje wiele biologicznych formatów zapisu danych, które zapewniają czytelność zarówno dla człowieka jak i oprogramowania do analizy bioinformatycznej. Aby zapewnić tę czytelność, dane muszą być zgodne ze ścisłymi zasadami formatowania, co pozwala na zrozumienie zawartości danych przez każdą osobę na świecie, jak również zautomatyzowaną obróbkę standardowym oprogramowaniem [1]. W tym artykule opisano 6 wybranych formatów danych, tj.: FASTA, FASTQ, SAM/BAM, VCF, BED oraz GFF3, które są szeroko stosowane w aplikacjach bioinformatycznych.

FASTA

Format FASTA jest jednym z najprostszych i standardowo używanych formatów do przechowywania sekwencji DNA, RNA lub białek. Sekwencje te mogą być różnej długości, od kilku par zasad do całych chromosomów. Plik rozpoczyna się od nagłówka pierwszej sekwencji biologicznej. Nagłówek zawsze zapoczątkowany jest znakiem większości (“>”), a tuż po nim zgaduje się numer identyfikacyjny, po którym z kolei może wystąpić opis sekwencji. Dodatkowy opis nie jest wymagany i jest ignorowany przez programy do analizy danych. Następne linie stanowią już właściwą sekwencję nukleotydów (lub białek) zapisanych za pomocą jednoliterowych symboli. Każda linia ma ograniczoną długość i zawiera od 60 do 80 znaków [2,3]. Typowym rozszerzeniem plików w formacie FASTA są „.fasta” lub „.fa”. Przykład sekwencji DNA człowieka (ID: M12531.1) [4] został przedstawiony na Rysunku 1. Prezentowany plik posiada tylko jedną sekwencję (fragment genomu mitochondrialnego), niemniej format FASTA może przechowywać ich więcej (format Multi-FASTA). W takim przypadku każda z sekwencji będzie poprzedzona nagłówkiem (Rysunek 2) [5]. Pliki FASTA mają relatywnie małe rozmiary. Przedstawione poniżej pliki zajmują odpowiednio 704 i 18 KB. W przypadku większej ilości danych, np. całego genomu człowieka rozmiar pliku osiąga ok. 3 GB.

Rysunek 1. Format FASTA. W nagłówku, po znaku większości występuje unikalny identyfikator sekwencji (M12531.1), po czym następuje jej opis, który jest opcjonalny. Kolejne linie stanowią sekwencję nukleotydową.

Rysunek 2. Format Multi-FASTA. Każda sekwencja nukleotydów jest poprzedzona nagłówkiem.
W każdym nagłówku występuje unikalny identyfikator sekwencji (M12531.1, MN413237.1).

FASTQ

Format FASTQ służy do przechowywania danych wygenerowanych za pomocą technologii sekwencjonowania wysokoprzepustowego. Każda sekwencja nukleotydowa (tzw. „odczyt”, z ang. „read”) jest opisana przez 4 rekordy. Pierwsza, typowo zaczyna się od znaku “@”, po którym znajdują się informacje opisujące konkretny odczyt (jego identyfikator, atrybuty techniczne sekwencjonowania takie jak np. konkretne miejsce sekwencjonowania na płytce sekwencjonującej). Druga linia to odczytana sekwencja nukleotydowa. Trzecia zawiera znak “+” i opcjonalnie powtórzoną linię pierwszą. Linia czwarta zawiera jednoliterowy opis jakości (tzw. phred quality score) dla poszczególnych nukleotydów, zatem liczba znaków opisujących jakość jest równa liczbie nukleotydów. Jakość jest kodowana przez szereg znaków kodu ASCII i w zależności od platformy sekwencjonowania przybiera różne wartości. Zawsze jednak wyższa wartość jakości koresponduje z niższym prawdopodobieństwem popełnienia błędu w odczycie nukleotydu [6]. Kontrolę jakości sekwencji na podstawie tych znaków można przeprowadzić za pomocą wielu narzędzi, jednym z nich jest niezwykle popularny program FASTQC [7]. Standardowe rozszerzenie plików w formacie FASTQ to “.fq” lub “.fastq”. Co więcej, ze względu na ogromny rozmiar danych sekwencjonowania, pliki FASTQ są często kompresowane za pomocą narzędzi gzip lub bzip2, tak że końcowe rozszerzenie plików to „.fq.gz”, „.fastq.gz”, „.fq.bz2”, lub „.fastq.bz2”. Przykład pojedynczego odczytu (DRR000003; mRNA człowieka) w formacie FASTQ przedstawiono na Rysunku 3. [8]. Rozmiary plików w formacie FASTQ mogą osiągać do kilkudziesięciu, czy nawet kilkuset GB na próbkę (np. jeden pełny genom człowieka). Przetwarzanie i analiza danych takich rozmiarów jest możliwa na serwerach o dużej mocy obliczeniowej i dużej przestrzeni dyskowej.

Rysunek 3. Format FASTQ. Każda sekwencja nukleotydów (odczyt) jest poprzedzona nagłówkiem.
W nagłówku występuje identyfikator odczytu jak również informacje techniczne odnośnie samego procesu sekwencjonowania. Trzecią linię stanowi separator „+”, a czwarta zawiera znaki ASCII kodujące jakość odczytanego przez sekwenator nukleotydu.

Sequence Alignment/Map Format (SAM)

Ponieważ format FASTQ nie zawiera informacji odnośnie pochodzenia odczytów, w celu ustalenia kolejności ich występowania w genomie zachodzi potrzeba przyrównania (zmapowania) tych sekwencji do genomu referencyjnego (genom – wzór reprezentujący dany gatunek). Wynikiem tego procesu jest plik w formacie SAM. Format ten składa się z nagłówka (opcjonalny) oraz części właściwej. Linie nagłówka zaczynają się od znaku „@” i zawierają informacje uzupełniające do danych z części właściwej pliku, są to np. wersja genomu referencyjnego, nazwa programu użytego do przyrównania itp. Część właściwa ma 11 obowiązkowych pól, które opisują w jaki sposób odczyty przyrównały się do genomu (1 linia odpowiada jednemu odczytowi przyrównanemu w danej lokalizacji). Kolumny są rozdzielone tabulatorami i przechowują następujące informacje: (1) identyfikator przyrównanego odczytu (QNAME), (2) dokładną charakterystykę sposobu przyrównania odczytu zakodowaną w postaci flagi bitowej (np. orientacja nici DNA, dekoder flag jest dostępny pod adresem https://broadinstitute.github.io/picard/explain-flags.html; FLAG), (3) nazwę fragmentu genomu referencyjnego, do którego został przyrównany odczyt (np. numer chromosomu, RNAME), (4) pozycję pierwszej pary zasad przyrównanego odczytu (POS), (5) jakość przyrównania określająca prawdopodobieństwo pomyłki (phred quality score; MAPQ), (6) informację odnośnie liczby nukleotydów przyrównanych do danego miejsca uwzględniając także przerwy w tym przyrównaniu (np. 100M50D oznacza, że sekwencja 100 nukleotydów jak taka sama jak w genomie referencyjnym, podczas gdy 50 nukleotydów uległo delecji; CIGAR), (7) nazwę fragmentu genomu referencyjnego, do którego został przyrównany odczyt sparowany z odczytem, którego dotyczy dany wiersz (jeżeli dane są typu paired-end czyli fragment DNA jest sekwencjonowany w obu kierunkach; RNEXT) oraz (8) pozycję jego pierwszego nukleotydu (PNEXT), (9) dystans mierzony od pozycji pierwszego nukleotydu przyrównanego odczytu do ostatniej pozycji jego pary (jeżeli dane są typu paired-end; TLEN), (10) nukleotydy składające się na odczyt (SEQ) i ich (11) jakość (tzw. phred quality score; QUAL). Warto pamiętać, że nie wszystkie odczyty muszą przyrównać się do genomu referencyjnego jak również jeden odczyt może zostać przyrównany do kilku miejsc w genomie. W ostatnim przypadku jeden odczyt będzie opisany przez tyle linii do ilu miejsc został przyrównany. Pliki w formacie SAM są bogate w informacje i przez to osiągają ogromne rozmiary (kilkadziesiąt – kilkaset GB). Z tego powodu nie rekomenduje się zapisywania wyniku procesu przyrównania w formacie SAM, lecz do jego binarnego odpowiednika czyli formatu BAM, który zawiera dokładnie te same informacje, ale ma znacznie mniejszy rozmiar (ok. 27 % oryginalnego rozmiaru pliku SAM). Pliki BAM nie są przystosowane do odczytu przez tradycyjne edytory tekstu i obsługuje się je za pomocą narzędzi pakietu SAMtools [9,10].

Rysunek 4. Format SAM. Rysunek przedstawia rekord dla jednego przyrównanego odczytu (kolumna 1). Ponieważ 1 linia „nie mieści” się na ekranie komputera, nastąpiło zawijanie tekstu (nie zmieniło to jednak formatowania danych i wciąż zajmują 1 wiersz). Odczyt zlokalizowany jest na chromosomie 13 (3), a pozycja pierwszego nukleotydu przyrównanego odczytu to 34 692 686 (4). Wszystkie 151 nukleotydów jest takich samych jak w genomie referencyjnym (6) i przyrównały się one do nici antysensownej (2). Fragment danych nie zawiera nagłówka.

Variant Call Format (VCF)

Jak sama nazwa wskazuje, Variant Call Format służy do przechowywania polimorfizmów genetycznych (np. SNP, insercji czy delecji) oraz informacji o genotypach tj. informacji czy mamy do czynienia z hetero- czy homozygotą pod kątem danego polimorfizmu. Rzeczywista lista polimorfizmów jest poprzedzona liniami tworzącymi nagłówek rozpoczynający się od ciągu znaków „##”. Nagłówek zawiera opisy chromosomów i zawiera definicję zmiennych używanych do opisu polimorfizmów w polu „INFO”. Po nagłówku następuje pojedynczy wiersz zaczynający się od „#”, który definiuje układ kolejnych ośmiu kolumn opisujących polimorfizmy. Kolumny są rozdzielone tabulatorami i definiują następujące wartości: (1) identyfikator chromosomu lub innej struktury genomu (CHROM), (2) pozycja polimorfizmu wyrażona w parach zasad (POS), (3) jego identyfikator (ID), (4) nukleotyd(y) genomu referencyjnego w danej pozycji (REF), (5) nukleotyd reprezentujący allel/allele alternatywne występujący/e w genomie badanym (ALT), (6) jakość określająca prawdopodobieństwo pomyłki przy detekcji tego polimorfizmu (QUAL), (7) oznaczenie czy polimorfizm przeszedł określone standardy jakości (FILTER) oraz (8) dodatkowe informacje (INFO). Zawartość ostatniej kolumny jest zależna od oprogramowania, które służyło do detekcji polimorfizmów, może zawierać np. informacje o rodzaju polimorfizmu (np. SNP, CNV, indel, duplikacja), genotypie (pole GT) i jego jakości (pole GQ), głębokości pokrycia w danej pozycji genomu (pole DP) i wiele innych. Brakujące wartości w kolumnach są kodowane przez kropkę [11]. Przykład formatu VCF zawierającego znane polimorfizmy człowieka przedstawiono poniżej (Rysunek 5). Dane pochodzą z bazy danych Ensembl [12,13] i są dostępne pod adresem https://www.ensembl.org/info/data/ftp/index.html. Rozmiar pliku w formacie VCF zależy od tego ile polimorfizmów i dla ilu osób przechowywane są dane. Przykładowy plik, który zawiera ponad 34 miliony polimorfizmów SNP i indel dla 1 człowieka osiąga rozmiar 3.9 GB.

Rysunek 5. Fragment pliku w formacie VCF (z pominięciem nagłówka). Każda linia odpowiada jednemu polimorfizmowi. W podanym przykładzie polimorfizmy typu SNV (Single Nucleotide Variation) znajdują się na chromosomie 1 (kolumna 1) w pozycjach od 10 001 do 10 020 (kolumna 2). Na pozycji 10 020 znajduje się również indel. Kolumna 3 zawiera ID polimorfizmu (jeśli jest on znany w bazie danych), kolumna 4 to nukleotyd, który występuje w genomie referencyjnym, a kolumna 5 – w genomie badanym. W powyższym przykładzie kolumny 6 i 7 nie zawierają informacji odnośnie jakości, co jest reprezentowane przez znak kropki. Kolumna 8 wskazuje na bazę danych, z której pochodzą polimorfizmy (dbSNP) oraz na typ polimorfizmu (SNV, indel).

Browser Extensible Data (BED)

Format BED standardowo wykorzystuje się do przechowywania koordynat regionów genomu i powiązanych z nimi adnotacjami. Te regiony mogą np. stanowić różne funkcjonalne struktury genomu jak geny, transkrypty czy polimorfizmy genetyczne. Podobnie jak w VCF, w formacie BED na jeden region genomu przypada jedna linia jego opisu. Format ten nie musi domyślnie zawierać nagłówka i składa się przynajmniej z 3 obligatoryjnych kolumn, które stanowią odpowiednio (1) numer chromosomu (lub identyfikator innego fragmentu genomu np. takiego o nieznanej lokalizacji; chrom), miejsce gdzie (2) zaczyna się (chromStart) i (3) kończy (chromEnd) dany region na chromosomie. Kolumny od 4 do 12 są opcjonalne i zawierają: (4) etykietę (nazwę regionu; name), (5) score przyjmujący wartości od 0 do 1000 i którego definicja jest zależna od rodzaju danych i oprogramowania, za pomocą którego wygenerowano plik, (6) orientacja nici DNA (dodatnia “+”, ujemna “-” lub “.” jeśli nie ma wskazanej nici; strand). Kolumny od 7 do 9 odnoszą się do reprezentacji graficznej danego regionu w przeglądarkach genomowych i opisują: (7) pozycję początkową (thickStart), od której adnotacja jest wyświetlana w pogrubiony sposób (np. kodon początkowy genu), (8) pozycję końcową (np. pogrubiony kodon stop; thickEnd) oraz (9) model przestrzeni barw (itemRgb) w postaci R,G,B (np. 255,0,0) określający kolor wyświetlania adnotacji. Kolumny od 10 do 12 odnoszą się do podelementów regionu genomu, którego dotyczy dana linia i zawierają: (10) liczbę podelementów (np. liczba eksonów danego genu; blockCount), (11) długość tych podelementów (blockSizes) oraz (12) lokalizację początkową każdego z nich (blockStarts) [14, 15]. Rysunek 6 przedstawia fragment pliku w formacie BED pochodzący z bazy danych Ensembl [16]. Ze względu na prostszą strukturę, rozmiary pliku w formacie BED są znacznie mniejsze niż wspominane wcześniej FASTQ czy SAM. Przykładowo, plik zawierający koordynaty wszystkich genów zlokalizowanych w genomie człowieka ma 33 MB.

Rysunek 6. Format BED. Każda linia odpowiada jednemu regionowi genomu. Pierwsze 3 kolumny są obligatoryjne i odpowiadają lokalizacji danego fragmentu w genomie, kolumna 4 to nazwa regionu, 5 – score, 6 – informacja odnośnie nici DNA, na której dany region się znajduje. Kolumny od 7 do 9 odnoszą się do reprezentacji graficznej danego regionu w przeglądarkach genomowych.

General Feature Format (GFF3)

Każda linia nagłówka w formacie GFF3 charakteryzuje się występowanie dwóch znaków „#”. Informacje w nagłówku mogą dotyczyć np. wersji formatu czy informacji, do jakiego regionu genomu odnosi się cały plik. Pozostała część pliku zawiera 9 kolumn oddzielonych od siebie tabulatorami. Kolumny te odpowiadają za opis sekwencji nukleotydowych i białkowych. Informacje zawarte w kolumnach to: (1) identyfikator sekwencji, (2) źródło danych (np. baza danych czy procedura/algorytm ich uzyskania), (3) typ cząsteczki/fragmentu genomu (np. mRNA czy ekson), (4) początkowa oraz (5) końcowa lokalizacja danej cząsteczki w genomie, (6) score, którego definicja jest zależna od rodzaju danych i oprogramowania, za pomocą którego wygenerowano plik (typowo E-value przy określeniu podobieństwa sekwencji lub P-value do predykcji genów), (7) orientacja nici DNA (dodatnia “+”, ujemna “-” lub “.” jeśli nie ma wskazanej nici), (8) phase występujący w przypadku sekwencji kodujących, przyjmujący wartości od 0 do 2, gdzie 0 oznacza, że pierwszy nukleotyd jest początkiem kodonu, 1 że drugi, a 2 że trzeci, (9) atrybuty czyli dodatkowe informacje o sekwencji (niektóre z atrybutów są predefiniowane, np. ID) [17,18]. Na rysunku 7 pokazane są 3 wiersze dotyczące genu OR4F5 u człowieka. Rozmiar pliku dla pełnego genomu człowieka osiąga około 1,6 GB.

Rysunek 7. Format GFF3. Przedstawione 3 wiersze (podobnie jak przy formacie SAM nastąpiło zawijanie tekstu) dotyczące genu OR4F5. Każdy wiesz zaczyna się od „1”, co wskazuje na nazwę chromosomu. Opis ma strukturę hierarchiczną, pierwszy wiersz dotyczy najbardziej ogólnej jednostki czyli genu (gene_id=ENSG00000186092), a poniższe jego transkryptu (transcript_id=ENST00000641515) oraz eksonu (exon_id=ENSE00003812156) (dodatkowe atrybuty kolumny 9). Dane pochodzą z bazy danych Ensembl i były ręcznie adnotowane przez zespół Human and Vertebrate Analysis and Annotation (havana; 2) [19]. Każdy z regionów ma taką samą lokalizację początkową, czyli pozycję 65 419 pz (4), lecz inna końcową, tj. 71 585 pz w przypadku genu i transkryptu oraz 65 433 pz w przypadku eksonu (5).

Umieszczone w tym artykule formaty danych zostały opisane w nieprzypadkowej kolejności. W przypadku przetwarzania i analizy danych pochodzących z sekwencjonowania wysokoprzepustowego (następnej generacji), będziemy mieli do czynienia ze wszystkimi z tych formatów. Jeżeli za cel przyjmiemy detekcję i adnotację polimorfizmów w genomie człowieka będziemy musieli podjąć następujące kroki: (1) kontrola jakości i czyszczenie danych (FASTQ), (2) przyrównanie krótkich odczytów (FASTQ) do genomu referencyjnego (FASTA) czego efektem będzie plik w formacie BAM i (3) detekcja polimorfizmów genetycznych (VCF). W przypadku (4) adnotacji pomocne mogą okazać się pliki dodatkowe w formatach VCF, BED czy GFF3. Omówiony schemat został przedstawiony na rysunku nr 8.

Rysunek 8. Schemat przedstawiający kroki przetwarzania danych od surowych danych z sekwenatara do polimorfizmów genetycznych (białe pola). Formaty plików wejściowych (In) oraz wyjściowych (Out) zostały umieszczone w szarych polach.

Źródła:

[1] Mielczarek M, Czech B, et al. Extraordinary Command Line: Basic Data Editing Tools for Biologists Dealing with Sequence Data. The Open Bioinformatics Journal 2020. 13 137-145; doi: 10.2174/1875036202013010137.

[2] Lipman DJ, Pearson WR. Rapid and sensitive protein similarity searches. Science. 1985 Mar 22;227(4693):1435-41. doi: 10.1126/science.2983426.

[3] Xiong Jin. Podstawy Bioinformatyki. Wydawnictwo Uniwersytetu Warszawskiego. 2011.

[4] Walberg MW, Clayton DA. In vitro transcription of human mitochondrial DNA. Identification of specific light strand transcripts from the displacement loop region. J Biol Chem. 1983 Jan 25;258(2):1268-75. PMID: 6571694.

[5] Starikovskaya E, Shalaurova S, Dryomov S, et al. Mitochondrial DNA Variation of Leber’s Hereditary Optic Neuropathy in Western Siberia. Cells. 2019;8(12):1574. Published 2019 Dec 4. doi:10.3390/cells8121574.

[6] Cock PJ, Fields CJ et al. The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants. Nucleic Acids Res. 2010;38(6):1767-1771. doi:10.1093/nar/gkp1137.

[7] Andrews S. FastQC: A Quality Control Tool for High Throughput Sequence Data. 2010. Available online at: http://www.bioinformatics.babraham.ac.uk/projects/fastqc/.

[8] Okuda T, Wakaguri H, Suzuki Y, Sugano S. Monitoring endoplasmic reticulum stress responsive mRNAs by RNA sequencing. Gene. 2012 May 25;500(1):63-72. doi: 10.1016/j.gene.2012.03.048. Epub 2012 Mar 17. PMID: 22465528.

[9] Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R; 1000 Genome Project Data Processing Subgroup. The Sequence Alignment/Map format and SAMtools. Bioinformatics. 2009 Aug 15;25(16):2078-9. doi: 10.1093/bioinformatics/btp352. Epub 2009 Jun 8. PMID: 19505943; PMCID: PMC2723002.

[10] https://samtools.github.io/hts-specs/SAMv1.pdf

[11] Danecek P, Auton A et al. The variant call format and VCFtools. Bioinformatics. 2011;27(15):2156-2158. doi:10.1093/bioinformatics/btr330.

[12] Hunt SE, McLaren W et al. Ensembl variation resources. Database (Oxford). 2018 Jan 1;2018:bay119. doi: 10.1093/database/bay119.

[13] www.ensembl.org.

[14] https://genome.ucsc.edu/FAQ/FAQformat.html#format1.

[15] Aaron R. Quinlan et al. BEDTools: a flexible suite of utilities for comparing genomic features, Bioinformatics, Volume 26, Issue 6, 15 March 2010, Pages 841–842, https://doi.org/10.1093/bioinformatics/btq033.

[16] https://www.ensembl.org/info/website/upload/bed.html#required.

[17] https://github.com/The-Sequence-Ontology/Specifications/blob/master/gff3.md).

[18] https://www.ensembl.org/info/website/upload/gff3.html

[19] Hancock J. HAVANA (Human and Vertebrate Analysis and Annotation). Wiley Online Library, 2004. 10.1002/9780471650126.dob0940.


Fakty i Mity Genetyki tworzone są przez pasjonatów, specjalistów w swoich dziedzinach.
Ten artykuł czytasz za darmo, bez reklam, bez spamu. Doceń naszą pracę i postaw nam wirtualną kawę 🙂
Dziękujemy! – Wasza Redakcja FiMG

Postaw mi kawę na buycoffee.to


Podziel się: