Bölüm 1: Yönteme Genel Bakış¶
Yapay zeka destekli çeviri - daha fazla bilgi ve iyileştirme önerileri
Varyant çağırma, bir genom dizisindeki varyasyonları referans genoma göre belirlemeyi amaçlayan bir genomik analiz yöntemidir. Burada, tüm genom dizileme verilerinde kısa germline varyantları, yani SNP'ler ve indel'leri çağırmak için tasarlanmış araç ve yöntemleri kullanacağız.

Tam bir varyant çağırma pipeline'ı genellikle referansa haritalama (bazen genom hizalaması olarak da adlandırılır) ve varyant filtreleme ile önceliklendirme dahil olmak üzere pek çok adım içerir. Bu kursta işleri basit tutmak adına yalnızca varyant çağırma kısmına odaklanacağız.
Yöntemler¶
Tüm genom dizileme örneklerinde germline SNP'leri ve indel'leri belirlemek için varyant çağırmayı uygulamanın iki yolunu göstereceğiz. Önce, her örnekten bağımsız olarak varyant çağıran basit bir örnek başına yaklaşım ile başlayacağız. Ardından, birden fazla örneği birlikte analiz ederek daha doğru ve bilgilendirici sonuçlar üreten daha gelişmiş bir ortak çağırma yaklaşımını göstereceğiz.
Her iki yaklaşım için herhangi bir iş akışı kodu yazmaya başlamadan önce, komutları bazı test verileri üzerinde manuel olarak deneyeceğiz.
Veri Seti¶
Aşağıdaki veri ve ilgili kaynakları sağlıyoruz:
- Bir referans genomu: İnsan kromozom 20'sinin küçük bir bölgesinden (hg19/b37) oluşan ve yardımcı dosyalarını (indeks ve dizi sözlüğü) içeren referans genomu.
- Üç tüm genom dizileme örneği: Bir aile üçlüsüne (anne, baba ve oğul) karşılık gelen ve dosya boyutlarını küçük tutmak amacıyla kromozom 20 üzerinde küçük bir veri dilimine indirgenmiş örnekler. Bu, referans genoma önceden haritalanmış Illumina kısa okuma dizileme verisidir ve BAM formatında (Binary Alignment Map; SAM yani Sequence Alignment Map'in sıkıştırılmış bir sürümü) sağlanmaktadır.
- Genomik aralıklar listesi: Örneklerimizin varyant çağırmaya uygun veriye sahip olduğu genomdaki koordinatlar; BED formatında sağlanmaktadır.
Yazılım¶
Kullanılan iki ana araç şunlardır: Dizi hizalama dosyalarını işlemek için yaygın olarak kullanılan bir araç seti olan Samtools ve Broad Institute'ta geliştirilen varyant keşfi araçları seti olan GATK (Genome Analysis Toolkit).
Bu araçlar GitHub Codespaces ortamında kurulu değildir; bu nedenle Seqera Containers hizmeti aracılığıyla alınan konteynerler üzerinden kullanacağız (bkz. Hello Containers).
İpucu
pwd komutunu yazdığınızda gösterilen yolun son kısmının genomics olması için nf4-science/genomics dizininde olduğunuzdan emin olun.
1. Örnek Başına Varyant Çağırma¶
Örnek başına varyant çağırma, her örneği bağımsız olarak işler: varyant çağırıcı, bir seferde tek bir örneğin dizileme verilerini inceleyerek örneğin referanstan farklılaştığı konumları belirler.
Bu bölümde, örnek başına varyant çağırma yaklaşımını oluşturan iki komutu test ediyoruz: Samtools ile bir BAM dosyasını indeksleme ve GATK HaplotypeCaller ile varyant çağırma. Bu kursta Bölüm 2'de bir Nextflow iş akışına sarmalayacağımız komutlar bunlardır.
- Samtools kullanarak bir BAM girdi dosyası için indeks dosyası oluşturma
- VCF (Variant Call Format) formatında örnek başına varyant çağrıları üretmek için indekslenmiş BAM dosyası üzerinde GATK HaplotypeCaller'ı çalıştırma
Başlangıç olarak iki komutu yalnızca bir örnek üzerinde test ediyoruz.
1.1. Samtools ile BAM Girdi Dosyasını İndeksleme¶
İndeks dosyaları, biyoinformatik dosya formatlarının yaygın bir özelliğidir; GATK gibi araçların tüm dosyayı okumak zorunda kalmadan verinin bir alt kümesine erişmesini sağlayan, ana dosyanın yapısına ilişkin bilgiler içerirler. Bu dosyaların ne kadar büyük olabileceği göz önüne alındığında bu özellik büyük önem taşır.
BAM dosyaları çoğunlukla indeks olmadan sağlanır; bu nedenle pek çok analiz iş akışının ilk adımı samtools index kullanarak bir indeks oluşturmaktır.
Bir Samtools konteyneri indirip etkileşimli olarak başlatacak ve BAM dosyalarından biri üzerinde samtools index komutunu çalıştıracağız.
1.1.1. Samtools Konteynerini İndirme¶
Samtools konteyner imajını indirmek için docker pull komutunu çalıştırın:
Komut çıktısı
1.20--b5dfbd93de237464: Pulling from library/samtools
6360b3717211: Pull complete
2ec3f7ad9b3c: Pull complete
7716ca300600: Pull complete
4f4fb700ef54: Pull complete
8c61d418774c: Pull complete
03dae77ff45c: Pull complete
aab7f787139d: Pull complete
4f4fb700ef54: Pull complete
837d55536720: Pull complete
897362c12ca7: Pull complete
3893cbe24e91: Pull complete
d1b61e94977b: Pull complete
c72ff66fb90f: Pull complete
0e0388f29b6d: Pull complete
Digest: sha256:bbfc45b4f228975bde86cba95e303dd94ecf2fdacea5bfb2e2f34b0d7b141e41
Status: Downloaded newer image for community.wave.seqera.io/library/samtools:1.20--b5dfbd93de237464
community.wave.seqera.io/library/samtools:1.20--b5dfbd93de237464
Bu imajı daha önce indirmediyseniz tamamlanması bir dakika sürebilir. İşlem tamamlandığında konteyner imajının yerel bir kopyasına sahip olursunuz.
1.1.2. Samtools Konteynerini Etkileşimli Olarak Başlatma¶
Konteyneri etkileşimli olarak çalıştırmak için docker run komutunu -it bayraklarıyla kullanın.
-v ./data:/data seçeneği, araçların girdi dosyalarına erişebilmesi için yerel data dizinini konteynere bağlar.
İstemcinizin (base) root@a1b2c3d4e5f6:/tmp# gibi bir şeye dönüştüğünü fark edeceksiniz; bu, artık konteynerin içinde olduğunuzu gösterir.
/data/bam altındaki dizi verisi dosyalarını görebildiğinizi doğrulayın:
Artık ilk komutu denemeye hazırsınız.
1.1.3. İndeksleme Komutunu Çalıştırma¶
Samtools belgeleri, bir BAM dosyasını indekslemek için çalıştırılacak komut satırını bize verir.
Yalnızca girdi dosyasını sağlamamız yeterlidir; araç, girdi dosya adına .bai ekleyerek çıktı için otomatik olarak bir ad oluşturacaktır.
samtools index komutunu bir veri dosyası üzerinde çalıştırın:
Komut terminalde herhangi bir çıktı üretmez; ancak artık orijinal BAM girdi dosyasıyla aynı dizinde reads_mother.bam.bai adlı bir dosya görmeniz gerekir.
Dizin içeriği
Bu, ilk adımın testini tamamlar.
1.1.4. Samtools Konteynerinden Çıkma¶
Konteynerden çıkmak için exit yazın.
İstemciniz, konteyneri başlatmadan önceki haline dönmüş olmalıdır.
1.2. GATK HaplotypeCaller ile Varyant Çağırma¶
Az önce indekslediğimiz BAM dosyası üzerinde gatk HaplotypeCaller komutunu çalıştırmak istiyoruz.
1.2.1. GATK Konteynerini İndirme¶
Önce GATK konteyner imajını indirmek için docker pull komutunu çalıştıralım:
Komut çıktısı
Bazı katmanlar, daha önce indirdiğimiz Samtools konteyner imajıyla paylaşıldığından Already exists olarak görünür.
4.5.0.0--730ee8817e436867: Pulling from library/gatk4
6360b3717211: Already exists
2ec3f7ad9b3c: Already exists
7716ca300600: Already exists
4f4fb700ef54: Already exists
8c61d418774c: Already exists
03dae77ff45c: Already exists
aab7f787139d: Already exists
4f4fb700ef54: Already exists
837d55536720: Already exists
897362c12ca7: Already exists
3893cbe24e91: Already exists
d1b61e94977b: Already exists
e5c558f54708: Pull complete
087cce32d294: Pull complete
Digest: sha256:e33413b9100f834fcc62fd5bc9edc1e881e820aafa606e09301eac2303d8724b
Status: Downloaded newer image for community.wave.seqera.io/library/gatk4:4.5.0.0--730ee8817e436867
community.wave.seqera.io/library/gatk4:4.5.0.0--730ee8817e436867
İki konteyner imajı katmanlarının büyük bölümünü paylaştığından bu işlem ilk indirmeden daha hızlı olacaktır.
1.2.2. GATK Konteynerini Etkileşimli Olarak Başlatma¶
Samtools için yaptığımız gibi, veri dizini bağlı şekilde GATK konteynerini etkileşimli olarak başlatın.
İstemciniz, artık GATK konteynerinin içinde olduğunuzu gösterecek şekilde değişir.
1.2.3. Varyant Çağırma Komutunu Çalıştırma¶
GATK belgeleri, bir BAM dosyası üzerinde varyant çağırma işlemi gerçekleştirmek için çalıştırılacak komut satırını bize verir.
BAM girdi dosyasını (-I), referans genomunu (-R), çıktı dosyası için bir ad (-O) ve analiz edilecek genomik aralıkların listesini (-L) sağlamamız gerekir.
Ancak indeks dosyasının yolunu belirtmemize gerek yoktur; araç, yerleşik adlandırma ve birlikte konumlandırma kuralına dayanarak onu aynı dizinde otomatik olarak arayacaktır.
Aynı durum referans genomun yardımcı dosyaları (indeks ve dizi sözlüğü dosyaları, *.fai ve *.dict) için de geçerlidir.
gatk HaplotypeCaller \
-R /data/ref/ref.fasta \
-I /data/bam/reads_mother.bam \
-O /data/vcf/reads_mother.vcf \
-L /data/ref/intervals.bed
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar HaplotypeCaller -R /data/ref/ref.fasta -I /data/bam/reads_mother.bam -O reads_mother.vcf -L /data/ref/intervals.bed
00:27:50.687 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
00:27:50.854 INFO HaplotypeCaller - ------------------------------------------------------------
00:27:50.858 INFO HaplotypeCaller - The Genome Analysis Toolkit (GATK) v4.5.0.0
00:27:50.858 INFO HaplotypeCaller - For support and documentation go to https://software.broadinstitute.org/gatk/
00:27:50.858 INFO HaplotypeCaller - Executing as root@a1fe8ff42d07 on Linux v6.10.14-linuxkit amd64
00:27:50.858 INFO HaplotypeCaller - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
00:27:50.859 INFO HaplotypeCaller - Start Date/Time: February 8, 2026 at 12:27:50 AM GMT
00:27:50.859 INFO HaplotypeCaller - ------------------------------------------------------------
00:27:50.859 INFO HaplotypeCaller - ------------------------------------------------------------
00:27:50.861 INFO HaplotypeCaller - HTSJDK Version: 4.1.0
00:27:50.861 INFO HaplotypeCaller - Picard Version: 3.1.1
00:27:50.861 INFO HaplotypeCaller - Built for Spark Version: 3.5.0
00:27:50.862 INFO HaplotypeCaller - HTSJDK Defaults.COMPRESSION_LEVEL : 2
00:27:50.862 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
00:27:50.862 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
00:27:50.863 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
00:27:50.864 INFO HaplotypeCaller - Deflater: IntelDeflater
00:27:50.864 INFO HaplotypeCaller - Inflater: IntelInflater
00:27:50.864 INFO HaplotypeCaller - GCS max retries/reopens: 20
00:27:50.864 INFO HaplotypeCaller - Requester pays: disabled
00:27:50.865 INFO HaplotypeCaller - Initializing engine
00:27:50.991 INFO FeatureManager - Using codec BEDCodec to read file file:///data/ref/intervals.bed
00:27:51.016 INFO IntervalArgumentCollection - Processing 6369 bp from intervals
00:27:51.029 INFO HaplotypeCaller - Done initializing engine
00:27:51.040 INFO NativeLibraryLoader - Loading libgkl_utils.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_utils.so
00:27:51.042 INFO NativeLibraryLoader - Loading libgkl_smithwaterman.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_smithwaterman.so
00:27:51.042 INFO SmithWatermanAligner - Using AVX accelerated SmithWaterman implementation
00:27:51.046 INFO HaplotypeCallerEngine - Disabling physical phasing, which is supported only for reference-model confidence output
00:27:51.063 INFO NativeLibraryLoader - Loading libgkl_pairhmm_omp.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_pairhmm_omp.so
00:27:51.085 INFO IntelPairHmm - Flush-to-zero (FTZ) is enabled when running PairHMM
00:27:51.086 INFO IntelPairHmm - Available threads: 10
00:27:51.086 INFO IntelPairHmm - Requested threads: 4
00:27:51.086 INFO PairHMM - Using the OpenMP multi-threaded AVX-accelerated native PairHMM implementation
00:27:51.128 INFO ProgressMeter - Starting traversal
00:27:51.136 INFO ProgressMeter - Current Locus Elapsed Minutes Regions Processed Regions/Minute
00:27:51.882 WARN InbreedingCoeff - InbreedingCoeff will not be calculated at position 20_10037292_10066351:3480 and possibly subsequent; at least 10 samples must have called genotypes
00:27:52.969 INFO HaplotypeCaller - 7 read(s) filtered by: MappingQualityReadFilter
0 read(s) filtered by: MappingQualityAvailableReadFilter
0 read(s) filtered by: MappedReadFilter
0 read(s) filtered by: NotSecondaryAlignmentReadFilter
0 read(s) filtered by: NotDuplicateReadFilter
0 read(s) filtered by: PassesVendorQualityCheckReadFilter
0 read(s) filtered by: NonZeroReferenceLengthAlignmentReadFilter
0 read(s) filtered by: GoodCigarReadFilter
0 read(s) filtered by: WellformedReadFilter
7 total reads filtered out of 1867 reads processed
00:27:52.971 INFO ProgressMeter - 20_10037292_10066351:13499 0.0 35 1145.7
00:27:52.971 INFO ProgressMeter - Traversal complete. Processed 35 total regions in 0.0 minutes.
00:27:52.976 INFO VectorLoglessPairHMM - Time spent in setup for JNI call : 0.003346916
00:27:52.976 INFO PairHMM - Total compute time in PairHMM computeLogLikelihoods() : 0.045731709
00:27:52.977 INFO SmithWatermanAligner - Total compute time in native Smith-Waterman : 0.02 sec
00:27:52.981 INFO HaplotypeCaller - Shutting down engine
[February 8, 2026 at 12:27:52 AM GMT] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.04 minutes.
Runtime.totalMemory()=203423744
Günlük çıktısı oldukça ayrıntılıdır; bu nedenle yukarıdaki örnekte en ilgili satırları vurguladık.
Çıktı dosyaları olan reads_mother.vcf ve indeks dosyası reads_mother.vcf.idx, konteynerdeki çalışma dizininizde oluşturulur.
VCF dosyası varyant çağrılarını içerir; bunu birazdan göreceğiz. İndeks dosyası ise BAM indeks dosyasıyla aynı işlevi görür: araçların tüm dosyayı yüklemeden verinin alt kümelerini aramasına ve almasına olanak tanır.
VCF metin tabanlı bir format olduğundan ve bu dosya küçük bir test dosyası olduğundan, içeriğini açıp görüntülemek için cat reads_mother.vcf komutunu çalıştırabilirsiniz.
Dosyanın başına kadar geri dönerseniz, pek çok meta veri satırından oluşan bir başlık ve ardından her biri bir satırda yer alan varyant çağrılarının listesini bulacaksınız.
Dosya içeriği (kısaltılmış)
Yukarıdaki örnek çıktıda, ardından gelen tablo verilerinin sütun adlarını veren son başlık satırını vurguladık. Her veri satırı, örneğin dizileme verilerinde tespit edilen olası bir varyantı tanımlar. VCF formatını yorumlamaya ilişkin rehberlik için bu faydalı makaleye bakabilirsiniz.
1.2.4. Çıktı Dosyalarını Taşıma¶
Konteynerin içinde kalan her şey gelecekteki çalışmalar için erişilemez olacaktır.
BAM indeks dosyası doğrudan bağlı dosya sistemindeki /data/bam dizininde oluşturuldu; ancak VCF dosyası ve indeksi oluşturulmadı; bu nedenle bu iki dosyayı manuel olarak taşımamız gerekiyor.
Dizin içeriği
Bu işlem tamamlandığında tüm dosyalar normal dosya sisteminizde erişilebilir hale gelir.
1.2.5. GATK Konteynerinden Çıkma¶
Konteynerden çıkmak için exit yazın.
İstemciniz normale dönmüş olmalıdır. Bu, örnek başına varyant çağırma testini tamamlar.
Bir iş akışı olarak yazın!
Bu analizi bir Nextflow iş akışı olarak uygulamaya başlamak istiyorsanız doğrudan Bölüm 2'ye geçebilirsiniz. Bölüm 3'e geçmeden önce ikinci test turunu tamamlamak için geri dönmeniz yeterli olacaktır.
2. Bir Kohort Üzerinde Ortak Çağırma¶
Az önce kullandığımız varyant çağırma yaklaşımı, örnek başına varyant çağrıları üretir. Her örneğin varyantlarına ayrı ayrı bakmak için bu yeterlidir; ancak sınırlı bilgi sağlar. Varyant çağrılarının birden fazla örnek arasında nasıl farklılaştığına bakmak çoğunlukla daha ilgi çekicidir. GATK bu amaç için ortak varyant çağırma adı verilen alternatif bir yöntem sunar.
Ortak varyant çağırma; her örnek için GVCF (Genomic VCF) adı verilen özel bir varyant çıktısı oluşturmayı, ardından tüm örneklerin GVCF verilerini birleştirmeyi ve bir 'ortak genotipleme' istatistiksel analizi çalıştırmayı kapsar.

Bir örneğin GVCF'sini özel kılan şey, yalnızca programın varyasyon kanıtı bulduğu konumları değil, genomun hedeflenen alanındaki tüm konumlara ilişkin dizi verisi istatistiklerini özetleyen kayıtlar içermesidir. Bu, ortak genotipleme hesaplaması için kritik öneme sahiptir (daha fazla okuma).
GVCF, az önce test ettiğimiz araç olan GATK HaplotypeCaller tarafından ek bir parametre (-ERC GVCF) ile üretilir.
GVCF'lerin birleştirilmesi, örnek başına çağrıları bir veri deposunda (veritabanına benzer şekilde) birleştiren GATK GenomicsDBImport ile yapılır.
Asıl 'ortak genotipleme' analizi ise GATK GenotypeGVCFs ile gerçekleştirilir.
Burada GVCF'leri oluşturmak ve ortak genotipleme çalıştırmak için gereken komutları test ediyoruz. Bu kursta Bölüm 3'te bir Nextflow iş akışına sarmalayacağımız komutlar bunlardır.
- Samtools kullanarak her BAM girdi dosyası için indeks dosyası oluşturma
- Örnek başına genomik varyant çağrılarının GVCF'sini oluşturmak için her BAM girdi dosyası üzerinde GATK HaplotypeCaller'ı çalıştırma
- Tüm GVCF'leri toplayıp bir GenomicsDB veri deposunda birleştirme
- Kohort düzeyinde bir VCF üretmek için birleştirilmiş GVCF veri deposu üzerinde ortak genotipleme çalıştırma
Şimdi üç BAM dosyasının tamamını indekslemekten başlayarak bu komutların hepsini test etmemiz gerekiyor.
2.1. Üç Örneğin Tamamı İçin BAM Dosyalarını İndeksleme¶
Yukarıdaki ilk bölümde yalnızca bir BAM dosyasını indeksledik. Şimdi GATK HaplotypeCaller'ın bunları işleyebilmesi için üç örneğin tamamını indeklememiz gerekiyor.
2.1.1. Samtools Konteynerini Etkileşimli Olarak Başlatma¶
Samtools konteyner imajını zaten indirdik; bu nedenle doğrudan başlatabiliriz:
İstemciniz, daha önce olduğu gibi veri dizini bağlı şekilde konteynerin içinde olduğunuzu gösterecek şekilde değişir.
2.1.2. İndeksleme Komutunu Üç Örneğin Tamamında Çalıştırma¶
İndeksleme komutunu üç BAM dosyasının her birinde çalıştırın:
samtools index /data/bam/reads_mother.bam
samtools index /data/bam/reads_father.bam
samtools index /data/bam/reads_son.bam
Dizin içeriği
Bu işlem, indeks dosyalarını ilgili BAM dosyalarıyla aynı dizinde oluşturmalıdır.
2.1.3. Samtools Konteynerinden Çıkma¶
Konteynerden çıkmak için exit yazın.
İstemciniz normale dönmüş olmalıdır.
2.2. Üç Örneğin Tamamı İçin GVCF Oluşturma¶
Ortak genotipleme adımını çalıştırmak için üç örneğin tamamına ait GVCF'lere ihtiyacımız var.
2.2.1. GATK Konteynerini Etkileşimli Olarak Başlatma¶
GATK konteyner imajını daha önce indirdik; bu nedenle doğrudan başlatabiliriz:
İstemciniz, GATK konteynerinin içinde olduğunuzu gösterecek şekilde değişir.
2.2.2. GVCF Seçeneğiyle Varyant Çağırma Komutunu Çalıştırma¶
Genomik bir VCF (GVCF) üretmek için temel komuta -ERC GVCF seçeneğini ekliyoruz; bu seçenek HaplotypeCaller'ın GVCF modunu etkinleştirir.
Ayrıca çıktı dosyasının uzantısını .vcf'den .g.vcf'ye değiştiriyoruz.
Bu teknik olarak zorunlu değildir; ancak güçlü biçimde önerilen bir kuraldır.
gatk HaplotypeCaller \
-R /data/ref/ref.fasta \
-I /data/bam/reads_mother.bam \
-O reads_mother.g.vcf \
-L /data/ref/intervals.bed \
-ERC GVCF
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar HaplotypeCaller -R /data/ref/ref.fasta -I /data/bam/reads_mother.bam -O reads_mother.g.vcf -L /data/ref/intervals.bed -ERC GVCF
16:51:00.620 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
16:51:00.749 INFO HaplotypeCaller - ------------------------------------------------------------
16:51:00.751 INFO HaplotypeCaller - The Genome Analysis Toolkit (GATK) v4.5.0.0
16:51:00.751 INFO HaplotypeCaller - For support and documentation go to https://software.broadinstitute.org/gatk/
16:51:00.751 INFO HaplotypeCaller - Executing as root@be1a0302f6c7 on Linux v6.8.0-1030-azure amd64
16:51:00.751 INFO HaplotypeCaller - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
16:51:00.752 INFO HaplotypeCaller - Start Date/Time: February 11, 2026 at 4:51:00 PM GMT
16:51:00.752 INFO HaplotypeCaller - ------------------------------------------------------------
16:51:00.752 INFO HaplotypeCaller - ------------------------------------------------------------
16:51:00.752 INFO HaplotypeCaller - HTSJDK Version: 4.1.0
16:51:00.753 INFO HaplotypeCaller - Picard Version: 3.1.1
16:51:00.753 INFO HaplotypeCaller - Built for Spark Version: 3.5.0
16:51:00.753 INFO HaplotypeCaller - HTSJDK Defaults.COMPRESSION_LEVEL : 2
16:51:00.753 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
16:51:00.753 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
16:51:00.754 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
16:51:00.754 INFO HaplotypeCaller - Deflater: IntelDeflater
16:51:00.754 INFO HaplotypeCaller - Inflater: IntelInflater
16:51:00.754 INFO HaplotypeCaller - GCS max retries/reopens: 20
16:51:00.754 INFO HaplotypeCaller - Requester pays: disabled
16:51:00.755 INFO HaplotypeCaller - Initializing engine
16:51:00.893 INFO FeatureManager - Using codec BEDCodec to read file file:///data/ref/intervals.bed
16:51:00.905 INFO IntervalArgumentCollection - Processing 6369 bp from intervals
16:51:00.910 INFO HaplotypeCaller - Done initializing engine
16:51:00.912 INFO HaplotypeCallerEngine - Tool is in reference confidence mode and the annotation, the following changes will be made to any specified annotations: 'StrandBiasBySample' will be enabled. 'ChromosomeCounts', 'FisherStrand', 'StrandOddsRatio' and 'QualByDepth' annotations have been disabled
16:51:00.917 INFO NativeLibraryLoader - Loading libgkl_utils.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_utils.so
16:51:00.919 INFO NativeLibraryLoader - Loading libgkl_smithwaterman.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_smithwaterman.so
16:51:00.919 INFO SmithWatermanAligner - Using AVX accelerated SmithWaterman implementation
16:51:00.923 INFO HaplotypeCallerEngine - Standard Emitting and Calling confidence set to -0.0 for reference-model confidence output
16:51:00.923 INFO HaplotypeCallerEngine - All sites annotated with PLs forced to true for reference-model confidence output
16:51:00.933 INFO NativeLibraryLoader - Loading libgkl_pairhmm_omp.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_pairhmm_omp.so
16:51:00.945 INFO IntelPairHmm - Flush-to-zero (FTZ) is enabled when running PairHMM
16:51:00.945 INFO IntelPairHmm - Available threads: 4
16:51:00.945 INFO IntelPairHmm - Requested threads: 4
16:51:00.945 INFO PairHMM - Using the OpenMP multi-threaded AVX-accelerated native PairHMM implementation
16:51:00.984 INFO ProgressMeter - Starting traversal
16:51:00.985 INFO ProgressMeter - Current Locus Elapsed Minutes Regions Processed Regions/Minute
16:51:01.452 WARN InbreedingCoeff - InbreedingCoeff will not be calculated at position 20_10037292_10066351:3480 and possibly subsequent; at least 10 samples must have called genotypes
16:51:02.358 INFO HaplotypeCaller - 7 read(s) filtered by: MappingQualityReadFilter
0 read(s) filtered by: MappingQualityAvailableReadFilter
0 read(s) filtered by: MappedReadFilter
0 read(s) filtered by: NotSecondaryAlignmentReadFilter
0 read(s) filtered by: NotDuplicateReadFilter
0 read(s) filtered by: PassesVendorQualityCheckReadFilter
0 read(s) filtered by: NonZeroReferenceLengthAlignmentReadFilter
0 read(s) filtered by: GoodCigarReadFilter
0 read(s) filtered by: WellformedReadFilter
7 total reads filtered out of 1867 reads processed
16:51:02.359 INFO ProgressMeter - 20_10037292_10066351:13499 0.0 35 1529.5
16:51:02.359 INFO ProgressMeter - Traversal complete. Processed 35 total regions in 0.0 minutes.
16:51:02.361 INFO VectorLoglessPairHMM - Time spent in setup for JNI call : 0.0022800000000000003
16:51:02.361 INFO PairHMM - Total compute time in PairHMM computeLogLikelihoods() : 0.061637120000000004
16:51:02.361 INFO SmithWatermanAligner - Total compute time in native Smith-Waterman : 0.02 sec
16:51:02.362 INFO HaplotypeCaller - Shutting down engine
[February 11, 2026 at 4:51:02 PM GMT] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.03 minutes.
Runtime.totalMemory()=257949696
Bu işlem, konteynerdeki geçerli çalışma dizininde reads_mother.g.vcf GVCF çıktı dosyasını ve indeks dosyası olan reads_mother.g.vcf.idx'i oluşturur.
Dosya içeriğinin ilk 200 satırını görüntülemek için head -200 reads_mother.g.vcf komutunu çalıştırırsanız, dosyanın ilk bölümde oluşturduğumuz eşdeğer VCF'den çok daha uzun olduğunu ve satırların büyük bölümünün VCF'de gördüklerimizden oldukça farklı göründüğünü fark edeceksiniz.
Dosya içeriği (kısaltılmış)
| reads_mother.g.vcf | |
|---|---|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 | |
Son başlık satırını ve dosyadaki ilk üç 'gerçek' varyant çağrısını bir kez daha vurguladık.
Varyant çağrı satırlarının, varyant çağırıcının varyasyon kanıtı bulmadığı varyant olmayan bölgeleri temsil eden pek çok varyant dışı satırla iç içe geçtiğini fark edeceksiniz. Yukarıda kısaca belirtildiği gibi, GVCF modunun varyant çağırmayı özel kılan özelliği budur: varyant çağırıcı, varyasyon yokluğuna ilişkin güven düzeyini tanımlayan bazı istatistikleri kaydeder. Bu sayede iki çok farklı durum arasında ayrım yapmak mümkün olur: (1) örneğin homozigot-referans olduğunu gösteren kaliteli veriler mevcuttur ve (2) herhangi bir yönde karar vermek için yeterli kaliteli veri bulunmamaktadır.
Bunun gibi bir GVCF'de genellikle çok sayıda varyant dışı satır bulunur; bunların arasına daha az sayıda varyant kaydı serpiştirilmiştir.
2.2.3. İşlemi Diğer İki Örnek İçin Tekrarlama¶
Şimdi aşağıdaki komutları birbiri ardına çalıştırarak kalan iki örnek için GVCF'leri oluşturalım.
gatk HaplotypeCaller \
-R /data/ref/ref.fasta \
-I /data/bam/reads_father.bam \
-O reads_father.g.vcf \
-L /data/ref/intervals.bed \
-ERC GVCF
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar HaplotypeCaller -R /data/ref/ref.fasta -I /data/bam/reads_father.bam -O reads_father.g.vcf -L /data/ref/intervals.bed -ERC GVCF
17:28:30.677 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
17:28:30.801 INFO HaplotypeCaller - ------------------------------------------------------------
17:28:30.803 INFO HaplotypeCaller - The Genome Analysis Toolkit (GATK) v4.5.0.0
17:28:30.804 INFO HaplotypeCaller - For support and documentation go to https://software.broadinstitute.org/gatk/
17:28:30.804 INFO HaplotypeCaller - Executing as root@be1a0302f6c7 on Linux v6.8.0-1030-azure amd64
17:28:30.804 INFO HaplotypeCaller - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
17:28:30.804 INFO HaplotypeCaller - Start Date/Time: February 11, 2026 at 5:28:30 PM GMT
17:28:30.804 INFO HaplotypeCaller - ------------------------------------------------------------
17:28:30.804 INFO HaplotypeCaller - ------------------------------------------------------------
17:28:30.805 INFO HaplotypeCaller - HTSJDK Version: 4.1.0
17:28:30.805 INFO HaplotypeCaller - Picard Version: 3.1.1
17:28:30.805 INFO HaplotypeCaller - Built for Spark Version: 3.5.0
17:28:30.806 INFO HaplotypeCaller - HTSJDK Defaults.COMPRESSION_LEVEL : 2
17:28:30.806 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
17:28:30.806 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
17:28:30.806 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
17:28:30.806 INFO HaplotypeCaller - Deflater: IntelDeflater
17:28:30.807 INFO HaplotypeCaller - Inflater: IntelInflater
17:28:30.807 INFO HaplotypeCaller - GCS max retries/reopens: 20
17:28:30.807 INFO HaplotypeCaller - Requester pays: disabled
17:28:30.807 INFO HaplotypeCaller - Initializing engine
17:28:30.933 INFO FeatureManager - Using codec BEDCodec to read file file:///data/ref/intervals.bed
17:28:30.946 INFO IntervalArgumentCollection - Processing 6369 bp from intervals
17:28:30.951 INFO HaplotypeCaller - Done initializing engine
17:28:30.953 INFO HaplotypeCallerEngine - Tool is in reference confidence mode and the annotation, the following changes will be made to any specified annotations: 'StrandBiasBySample' will be enabled. 'ChromosomeCounts', 'FisherStrand', 'StrandOddsRatio' and 'QualByDepth' annotations have been disabled
17:28:30.957 INFO NativeLibraryLoader - Loading libgkl_utils.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_utils.so
17:28:30.959 INFO NativeLibraryLoader - Loading libgkl_smithwaterman.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_smithwaterman.so
17:28:30.960 INFO SmithWatermanAligner - Using AVX accelerated SmithWaterman implementation
17:28:30.963 INFO HaplotypeCallerEngine - Standard Emitting and Calling confidence set to -0.0 for reference-model confidence output
17:28:30.963 INFO HaplotypeCallerEngine - All sites annotated with PLs forced to true for reference-model confidence output
17:28:30.972 INFO NativeLibraryLoader - Loading libgkl_pairhmm_omp.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_pairhmm_omp.so
17:28:30.987 INFO IntelPairHmm - Flush-to-zero (FTZ) is enabled when running PairHMM
17:28:30.987 INFO IntelPairHmm - Available threads: 4
17:28:30.987 INFO IntelPairHmm - Requested threads: 4
17:28:30.987 INFO PairHMM - Using the OpenMP multi-threaded AVX-accelerated native PairHMM implementation
17:28:31.034 INFO ProgressMeter - Starting traversal
17:28:31.034 INFO ProgressMeter - Current Locus Elapsed Minutes Regions Processed Regions/Minute
17:28:31.570 WARN InbreedingCoeff - InbreedingCoeff will not be calculated at position 20_10037292_10066351:3480 and possibly subsequent; at least 10 samples must have called genotypes
17:28:32.865 INFO HaplotypeCaller - 9 read(s) filtered by: MappingQualityReadFilter
0 read(s) filtered by: MappingQualityAvailableReadFilter
0 read(s) filtered by: MappedReadFilter
0 read(s) filtered by: NotSecondaryAlignmentReadFilter
0 read(s) filtered by: NotDuplicateReadFilter
0 read(s) filtered by: PassesVendorQualityCheckReadFilter
0 read(s) filtered by: NonZeroReferenceLengthAlignmentReadFilter
0 read(s) filtered by: GoodCigarReadFilter
0 read(s) filtered by: WellformedReadFilter
9 total reads filtered out of 2064 reads processed
17:28:32.866 INFO ProgressMeter - 20_10037292_10066351:13338 0.0 38 1245.2
17:28:32.866 INFO ProgressMeter - Traversal complete. Processed 38 total regions in 0.0 minutes.
17:28:32.868 INFO VectorLoglessPairHMM - Time spent in setup for JNI call : 0.0035923200000000004
17:28:32.868 INFO PairHMM - Total compute time in PairHMM computeLogLikelihoods() : 0.10765202500000001
17:28:32.868 INFO SmithWatermanAligner - Total compute time in native Smith-Waterman : 0.03 sec
17:28:32.869 INFO HaplotypeCaller - Shutting down engine
[February 11, 2026 at 5:28:32 PM GMT] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.04 minutes.
Runtime.totalMemory()=299892736
gatk HaplotypeCaller \
-R /data/ref/ref.fasta \
-I /data/bam/reads_son.bam \
-O reads_son.g.vcf \
-L /data/ref/intervals.bed \
-ERC GVCF
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar HaplotypeCaller -R /data/ref/ref.fasta -I /data/bam/reads_son.bam -O reads_son.g.vcf -L /data/ref/intervals.bed -ERC GVCF
17:30:10.017 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
17:30:10.156 INFO HaplotypeCaller - ------------------------------------------------------------
17:30:10.159 INFO HaplotypeCaller - The Genome Analysis Toolkit (GATK) v4.5.0.0
17:30:10.159 INFO HaplotypeCaller - For support and documentation go to https://software.broadinstitute.org/gatk/
17:30:10.159 INFO HaplotypeCaller - Executing as root@be1a0302f6c7 on Linux v6.8.0-1030-azure amd64
17:30:10.159 INFO HaplotypeCaller - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
17:30:10.159 INFO HaplotypeCaller - Start Date/Time: February 11, 2026 at 5:30:09 PM GMT
17:30:10.159 INFO HaplotypeCaller - ------------------------------------------------------------
17:30:10.160 INFO HaplotypeCaller - ------------------------------------------------------------
17:30:10.160 INFO HaplotypeCaller - HTSJDK Version: 4.1.0
17:30:10.160 INFO HaplotypeCaller - Picard Version: 3.1.1
17:30:10.161 INFO HaplotypeCaller - Built for Spark Version: 3.5.0
17:30:10.161 INFO HaplotypeCaller - HTSJDK Defaults.COMPRESSION_LEVEL : 2
17:30:10.161 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
17:30:10.161 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
17:30:10.161 INFO HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
17:30:10.161 INFO HaplotypeCaller - Deflater: IntelDeflater
17:30:10.162 INFO HaplotypeCaller - Inflater: IntelInflater
17:30:10.162 INFO HaplotypeCaller - GCS max retries/reopens: 20
17:30:10.162 INFO HaplotypeCaller - Requester pays: disabled
17:30:10.162 INFO HaplotypeCaller - Initializing engine
17:30:10.277 INFO FeatureManager - Using codec BEDCodec to read file file:///data/ref/intervals.bed
17:30:10.290 INFO IntervalArgumentCollection - Processing 6369 bp from intervals
17:30:10.296 INFO HaplotypeCaller - Done initializing engine
17:30:10.298 INFO HaplotypeCallerEngine - Tool is in reference confidence mode and the annotation, the following changes will be made to any specified annotations: 'StrandBiasBySample' will be enabled. 'ChromosomeCounts', 'FisherStrand', 'StrandOddsRatio' and 'QualByDepth' annotations have been disabled
17:30:10.302 INFO NativeLibraryLoader - Loading libgkl_utils.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_utils.so
17:30:10.303 INFO NativeLibraryLoader - Loading libgkl_smithwaterman.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_smithwaterman.so
17:30:10.304 INFO SmithWatermanAligner - Using AVX accelerated SmithWaterman implementation
17:30:10.307 INFO HaplotypeCallerEngine - Standard Emitting and Calling confidence set to -0.0 for reference-model confidence output
17:30:10.307 INFO HaplotypeCallerEngine - All sites annotated with PLs forced to true for reference-model confidence output
17:30:10.315 INFO NativeLibraryLoader - Loading libgkl_pairhmm_omp.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_pairhmm_omp.so
17:30:10.328 INFO IntelPairHmm - Flush-to-zero (FTZ) is enabled when running PairHMM
17:30:10.329 INFO IntelPairHmm - Available threads: 4
17:30:10.329 INFO IntelPairHmm - Requested threads: 4
17:30:10.329 INFO PairHMM - Using the OpenMP multi-threaded AVX-accelerated native PairHMM implementation
17:30:10.368 INFO ProgressMeter - Starting traversal
17:30:10.369 INFO ProgressMeter - Current Locus Elapsed Minutes Regions Processed Regions/Minute
17:30:10.875 WARN InbreedingCoeff - InbreedingCoeff will not be calculated at position 20_10037292_10066351:3480 and possibly subsequent; at least 10 samples must have called genotypes
17:30:11.980 INFO HaplotypeCaller - 14 read(s) filtered by: MappingQualityReadFilter
0 read(s) filtered by: MappingQualityAvailableReadFilter
0 read(s) filtered by: MappedReadFilter
0 read(s) filtered by: NotSecondaryAlignmentReadFilter
0 read(s) filtered by: NotDuplicateReadFilter
0 read(s) filtered by: PassesVendorQualityCheckReadFilter
0 read(s) filtered by: NonZeroReferenceLengthAlignmentReadFilter
0 read(s) filtered by: GoodCigarReadFilter
0 read(s) filtered by: WellformedReadFilter
14 total reads filtered out of 1981 reads processed
17:30:11.981 INFO ProgressMeter - 20_10037292_10066351:13223 0.0 35 1302.7
17:30:11.981 INFO ProgressMeter - Traversal complete. Processed 35 total regions in 0.0 minutes.
17:30:11.983 INFO VectorLoglessPairHMM - Time spent in setup for JNI call : 0.0034843710000000004
17:30:11.983 INFO PairHMM - Total compute time in PairHMM computeLogLikelihoods() : 0.048108363
17:30:11.983 INFO SmithWatermanAligner - Total compute time in native Smith-Waterman : 0.02 sec
17:30:11.984 INFO HaplotypeCaller - Shutting down engine
[February 11, 2026 at 5:30:11 PM GMT] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.03 minutes.
Runtime.totalMemory()=226492416
Bu işlem tamamlandığında, geçerli dizininizde .g.vcf ile biten üç dosya (her örnek için bir tane) ve .g.vcf.idx ile biten ilgili indeks dosyaları bulunmalıdır.
Dizin içeriği
Bu noktada, girdi örneklerimizin her biri için GVCF modunda varyant çağrısı yaptık. Artık ortak çağırmaya geçme zamanı.
Ancak konteynerden çıkmayın! Bir sonraki adımda aynı konteyneri kullanacağız.
2.3. Ortak Genotiplemeyi Çalıştırma¶
Artık tüm GVCF'lere sahip olduğumuza göre, bir örnek kohortu için varyant çağrıları üretmek amacıyla ortak genotipleme yaklaşımını deneyebiliriz. Bu iki adımlı bir yöntemdir: tüm GVCF'lerin verilerini bir veri deposunda birleştirmek, ardından ortak çağrılan varyantların nihai VCF'ini oluşturmak için ortak genotipleme analizini çalıştırmak.
2.3.1. Örnek Başına Tüm GVCF'leri Birleştirme¶
Bu ilk adım, tüm GVCF'lerin verilerini bir GenomicsDB veri deposunda birleştirmek için GenomicsDBImport adlı başka bir GATK aracını kullanır. GenomicsDB veri deposu, varyant bilgileri için ara depolama görevi gören bir tür veritabanı formatıdır.
gatk GenomicsDBImport \
-V reads_mother.g.vcf \
-V reads_father.g.vcf \
-V reads_son.g.vcf \
-L /data/ref/intervals.bed \
--genomicsdb-workspace-path family_trio_gdb
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar GenomicsDBImport -V reads_mother.g.vcf -V reads_father.g.vcf -V reads_son.g.vcf -L /data/ref/intervals.bed --genomicsdb-workspace-path family_trio_gdb
17:37:07.569 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
17:37:07.699 INFO GenomicsDBImport - ------------------------------------------------------------
17:37:07.702 INFO GenomicsDBImport - The Genome Analysis Toolkit (GATK) v4.5.0.0
17:37:07.702 INFO GenomicsDBImport - For support and documentation go to https://software.broadinstitute.org/gatk/
17:37:07.703 INFO GenomicsDBImport - Executing as root@be1a0302f6c7 on Linux v6.8.0-1030-azure amd64
17:37:07.703 INFO GenomicsDBImport - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
17:37:07.704 INFO GenomicsDBImport - Start Date/Time: February 11, 2026 at 5:37:07 PM GMT
17:37:07.704 INFO GenomicsDBImport - ------------------------------------------------------------
17:37:07.704 INFO GenomicsDBImport - ------------------------------------------------------------
17:37:07.706 INFO GenomicsDBImport - HTSJDK Version: 4.1.0
17:37:07.706 INFO GenomicsDBImport - Picard Version: 3.1.1
17:37:07.707 INFO GenomicsDBImport - Built for Spark Version: 3.5.0
17:37:07.709 INFO GenomicsDBImport - HTSJDK Defaults.COMPRESSION_LEVEL : 2
17:37:07.709 INFO GenomicsDBImport - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
17:37:07.709 INFO GenomicsDBImport - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
17:37:07.710 INFO GenomicsDBImport - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
17:37:07.710 INFO GenomicsDBImport - Deflater: IntelDeflater
17:37:07.711 INFO GenomicsDBImport - Inflater: IntelInflater
17:37:07.711 INFO GenomicsDBImport - GCS max retries/reopens: 20
17:37:07.711 INFO GenomicsDBImport - Requester pays: disabled
17:37:07.712 INFO GenomicsDBImport - Initializing engine
17:37:07.883 INFO FeatureManager - Using codec BEDCodec to read file file:///data/ref/intervals.bed
17:37:07.886 INFO IntervalArgumentCollection - Processing 6369 bp from intervals
17:37:07.889 INFO GenomicsDBImport - Done initializing engine
17:37:08.560 INFO GenomicsDBLibLoader - GenomicsDB native library version : 1.5.1-84e800e
17:37:08.561 INFO GenomicsDBImport - Vid Map JSON file will be written to /tmp/family_trio_gdb/vidmap.json
17:37:08.561 INFO GenomicsDBImport - Callset Map JSON file will be written to /tmp/family_trio_gdb/callset.json
17:37:08.561 INFO GenomicsDBImport - Complete VCF Header will be written to /tmp/family_trio_gdb/vcfheader.vcf
17:37:08.561 INFO GenomicsDBImport - Importing to workspace - /tmp/family_trio_gdb
17:37:08.878 INFO GenomicsDBImport - Importing batch 1 with 3 samples
17:37:09.359 INFO GenomicsDBImport - Importing batch 1 with 3 samples
17:37:09.487 INFO GenomicsDBImport - Importing batch 1 with 3 samples
17:37:09.591 INFO GenomicsDBImport - Done importing batch 1/1
17:37:09.592 INFO GenomicsDBImport - Import completed!
17:37:09.592 INFO GenomicsDBImport - Shutting down engine
[February 11, 2026 at 5:37:09 PM GMT] org.broadinstitute.hellbender.tools.genomicsdb.GenomicsDBImport done. Elapsed time: 0.03 minutes.
Runtime.totalMemory()=113246208
Tool returned:
true
Bu adımın çıktısı, birleştirilmiş varyant verilerini çeşitli farklı dosyalar biçiminde tutan iç içe geçmiş dizinler içeren bir dizinden oluşur. İçeriği inceleyebilirsiniz; ancak bu veri deposu formatının doğrudan insanlar tarafından okunması için tasarlanmadığını hemen göreceksiniz.
İpucu
GATK, gerektiğinde veri deposundan varyant çağrısı verilerini incelemeyi ve çıkarmayı mümkün kılan araçlar içerir.
2.3.2. Ortak Genotipleme Analizini Çalıştırma¶
Bu ikinci adım, kohortdaki tüm örneklerde mevcut verileri göz önünde bulundurarak varyant istatistiklerini ve bireysel genotipleri yeniden hesaplamak için GenotypeGVCFs adlı başka bir GATK aracını kullanır.
Komut çıktısı
Using GATK jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar
Running:
java -Dsamjdk.use_async_io_read_samtools=false -Dsamjdk.use_async_io_write_samtools=true -Dsamjdk.use_async_io_write_tribble=false -Dsamjdk.compression_level=2 -jar /opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar GenotypeGVCFs -R /data/ref/ref.fasta -V gendb://family_trio_gdb -O family_trio.vcf
17:38:45.084 INFO NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/opt/conda/share/gatk4-4.5.0.0-0/gatk-package-4.5.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
17:38:45.217 INFO GenotypeGVCFs - ------------------------------------------------------------
17:38:45.220 INFO GenotypeGVCFs - The Genome Analysis Toolkit (GATK) v4.5.0.0
17:38:45.220 INFO GenotypeGVCFs - For support and documentation go to https://software.broadinstitute.org/gatk/
17:38:45.220 INFO GenotypeGVCFs - Executing as root@be1a0302f6c7 on Linux v6.8.0-1030-azure amd64
17:38:45.220 INFO GenotypeGVCFs - Java runtime: OpenJDK 64-Bit Server VM v17.0.11-internal+0-adhoc..src
17:38:45.221 INFO GenotypeGVCFs - Start Date/Time: February 11, 2026 at 5:38:45 PM GMT
17:38:45.221 INFO GenotypeGVCFs - ------------------------------------------------------------
17:38:45.221 INFO GenotypeGVCFs - ------------------------------------------------------------
17:38:45.221 INFO GenotypeGVCFs - HTSJDK Version: 4.1.0
17:38:45.222 INFO GenotypeGVCFs - Picard Version: 3.1.1
17:38:45.222 INFO GenotypeGVCFs - Built for Spark Version: 3.5.0
17:38:45.222 INFO GenotypeGVCFs - HTSJDK Defaults.COMPRESSION_LEVEL : 2
17:38:45.222 INFO GenotypeGVCFs - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
17:38:45.222 INFO GenotypeGVCFs - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
17:38:45.222 INFO GenotypeGVCFs - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
17:38:45.223 INFO GenotypeGVCFs - Deflater: IntelDeflater
17:38:45.223 INFO GenotypeGVCFs - Inflater: IntelInflater
17:38:45.223 INFO GenotypeGVCFs - GCS max retries/reopens: 20
17:38:45.223 INFO GenotypeGVCFs - Requester pays: disabled
17:38:45.223 INFO GenotypeGVCFs - Initializing engine
17:38:45.544 INFO GenomicsDBLibLoader - GenomicsDB native library version : 1.5.1-84e800e
17:38:45.561 INFO NativeGenomicsDB - pid=221 tid=222 No valid combination operation found for INFO field InbreedingCoeff - the field will NOT be part of INFO fields in the generated VCF records
17:38:45.561 INFO NativeGenomicsDB - pid=221 tid=222 No valid combination operation found for INFO field MLEAC - the field will NOT be part of INFO fields in the generated VCF records
17:38:45.561 INFO NativeGenomicsDB - pid=221 tid=222 No valid combination operation found for INFO field MLEAF - the field will NOT be part of INFO fields in the generated VCF records
17:38:45.577 INFO GenotypeGVCFs - Done initializing engine
17:38:45.615 INFO ProgressMeter - Starting traversal
17:38:45.615 INFO ProgressMeter - Current Locus Elapsed Minutes Variants Processed Variants/Minute
17:38:45.903 WARN InbreedingCoeff - InbreedingCoeff will not be calculated at position 20_10037292_10066351:3480 and possibly subsequent; at least 10 samples must have called genotypes
GENOMICSDB_TIMER,GenomicsDB iterator next() timer,Wall-clock time(s),0.07757032800000006,Cpu time(s),0.07253379200000037
17:38:46.421 INFO ProgressMeter - 20_10037292_10066351:13953 0.0 3390 252357.3
17:38:46.422 INFO ProgressMeter - Traversal complete. Processed 3390 total variants in 0.0 minutes.
17:38:46.423 INFO GenotypeGVCFs - Shutting down engine
[February 11, 2026 at 5:38:46 PM GMT] org.broadinstitute.hellbender.tools.walkers.GenotypeGVCFs done. Elapsed time: 0.02 minutes.
Runtime.totalMemory()=203423744
Bu işlem, konteynerdeki geçerli çalışma dizininde family_trio.vcf VCF çıktı dosyasını ve indeksi olan family_trio.vcf.idx'i oluşturur.
Bu da makul ölçüde küçük bir dosyadır; dolayısıyla dosya içeriğini görüntülemek için cat family_trio.vcf komutunu çalıştırabilir ve ilk birkaç varyant satırını bulmak için aşağı kaydırabilirsiniz.
Dosya içeriği (kısaltılmış)
Varyant çağrısı verilerinin başlangıcını işaretleyen son başlık satırını bir kez daha vurguladık.
Bu, daha önce oluşturduğumuz VCF'e benzer görünür; ancak bu sefer üç örneğin tamamı için genotip düzeyinde bilgiye sahibiz. Dosyadaki son üç sütun, vurgulanan başlık satırında gösterildiği gibi ID alanlarının alfabetik sırasına göre listelenen örneklerin genotip bloklarıdır.
Test aile üçlümüz için çağrılan genotiplere ilk varyant açısından bakarsak, babanın heterozigot-varyant (0/1), anne ve oğulun ise her ikisinin de homozigot-varyant (1/1) olduğunu görürüz.
Veri setinden çıkarmak istediğimiz bilgi türü nihayetinde budur!
2.3.3. Çıktı Dosyalarını Taşıma¶
Daha önce belirtildiği gibi, konteynerin içinde kalan her şey gelecekteki çalışmalar için erişilemez olacaktır. Konteynerden çıkmadan önce, GVCF dosyalarını, nihai çok örnekli VCF'i ve tüm indeks dosyalarını konteynerin dışındaki dosya sistemine manuel olarak taşıyacağız. Bu sayede tüm bu çalışmayı otomatikleştirmek için iş akışımızı oluştururken karşılaştıracak bir şeyimiz olacak.
Dizin içeriği" hl_lines="14-19 22-23
data
├── bam
│ ├── reads_father.bam
│ ├── reads_father.bam.bai
│ ├── reads_mother.bam
│ ├── reads_mother.bam.bai
│ ├── reads_son.bam
│ └── reads_son.bam.bai
├── ref
│ ├── intervals.bed
│ ├── ref.dict
│ ├── ref.fasta
│ └── ref.fasta.fai
├── samplesheet.csv
└── vcf
├── family_trio.vcf
├── family_trio.vcf.idx
├── reads_father.g.vcf
├── reads_father.g.vcf.idx
├── reads_mother.g.vcf
├── reads_mother.g.vcf.idx
├── reads_mother.vcf
├── reads_mother.vcf.idx
├── reads_son.g.vcf
└── reads_son.g.vcf.idx
Bu işlem tamamlandığında tüm dosyalar normal dosya sisteminizde erişilebilir hale gelir.
2.3.4. GATK Konteynerinden Çıkma¶
Konteynerden çıkmak için exit yazın.
İstemciniz normale dönmüş olmalıdır. Bu, ortak varyant çağırma komutlarının manuel testini tamamlar.
Özetle¶
Samtools indeksleme ve GATK varyant çağırma komutlarını ilgili konteynerlerinde nasıl test edeceğinizi öğrendiniz; GVCF'lerin nasıl oluşturulacağını ve birden fazla örnek üzerinde ortak genotiplemenin nasıl çalıştırılacağını da bu kapsamda gördünüz.
Sırada ne var?¶
Kısa bir mola verin, ardından aynı komutları konteynerleri kullanarak çalışmayı yürüten iş akışlarına nasıl saracağınızı öğrenmek için Bölüm 2'ye geçin.