GATK生成GVCF文件报错

Question

使用GATK为每个样本生成GVCF文件时，
GATK运行出现报错，未能正确输出结果

gatk HaplotypeCaller -R GWHBKHI00000000.genome.fasta --emit-ref-confidence GVCF -I s1.markdup.bam -O s1.g.vcf.gz

运行结果如下

Using GATK jar /home/tx/miniconda3/envs/gatk4/share/gatk4-4.3.0.0-0/gatk-package-4.3.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 /home/tx/miniconda3/envs/gatk4/share/gatk4-4.3.0.0-0/gatk-package-4.3.0.0-local.jar HaplotypeCaller -R GWHBKHI00000000.genome.fasta --emit-ref-confidence GVCF -I s1.markdup.bam -O s1.g.vcf.gz
10:01:58.146 INFO  NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/home/tx/miniconda3/envs/gatk4/share/gatk4-4.3.0.0-0/gatk-package-4.3.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
10:01:58.308 INFO  HaplotypeCaller - ------------------------------------------------------------
10:01:58.308 INFO  HaplotypeCaller - The Genome Analysis Toolkit (GATK) v4.3.0.0
10:01:58.308 INFO  HaplotypeCaller - For support and documentation go to https://software.broadinstitute.org/gatk/
10:01:58.308 INFO  HaplotypeCaller - Executing as tx@utx on Linux v6.8.0-45-generic amd64
10:01:58.309 INFO  HaplotypeCaller - Java runtime: OpenJDK 64-Bit Server VM v11.0.13+7-b1751.21
10:01:58.309 INFO  HaplotypeCaller - Start Date/Time: 2024年9月22日 CST 上午10:01:58
10:01:58.309 INFO  HaplotypeCaller - ------------------------------------------------------------
10:01:58.309 INFO  HaplotypeCaller - ------------------------------------------------------------
10:01:58.310 INFO  HaplotypeCaller - HTSJDK Version: 3.0.1
10:01:58.310 INFO  HaplotypeCaller - Picard Version: 2.27.5
10:01:58.310 INFO  HaplotypeCaller - Built for Spark Version: 2.4.5
10:01:58.310 INFO  HaplotypeCaller - HTSJDK Defaults.COMPRESSION_LEVEL : 2
10:01:58.310 INFO  HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
10:01:58.310 INFO  HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
10:01:58.310 INFO  HaplotypeCaller - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
10:01:58.310 INFO  HaplotypeCaller - Deflater: IntelDeflater
10:01:58.310 INFO  HaplotypeCaller - Inflater: IntelInflater
10:01:58.310 INFO  HaplotypeCaller - GCS max retries/reopens: 20
10:01:58.310 INFO  HaplotypeCaller - Requester pays: disabled
10:01:58.311 INFO  HaplotypeCaller - Initializing engine
10:01:58.313 INFO  HaplotypeCaller - Shutting down engine
[2024年9月22日 CST 上午10:01:58] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.00 minutes.
Runtime.totalMemory()=620756992
***********************************************************************

A USER ERROR has occurred: Fasta dict file file:///media/tx/TX/S/922/GWHBKHI00000000.genome.dict for reference file:///media/tx/TX/S/922/GWHBKHI00000000.genome.fasta does not exist. Please see http://gatkforums.broadinstitute.org/discussion/1601/how-can-i-prepare-a-fasta-file-to-use-as-reference for help creating it.

***********************************************************************
Set the system property GATK_STACKTRACE_ON_USER_EXCEPTION (--java-options '-DGATK_STACKTRACE_ON_USER_EXCEPTION=true') to print the stack trace.

我文件夹下的文件如下：

-rwxr-xr-x 1 tx tx 216M  9月 22 08:17 gindex.1.ht2
-rwxr-xr-x 1 tx tx 159M  9月 22 08:17 gindex.2.ht2
-rwxr-xr-x 1 tx tx  314  9月 22 08:11 gindex.3.ht2
-rwxr-xr-x 1 tx tx 159M  9月 22 08:11 gindex.4.ht2
-rwxr-xr-x 1 tx tx 278M  9月 22 08:18 gindex.5.ht2
-rwxr-xr-x 1 tx tx 162M  9月 22 08:18 gindex.6.ht2
-rwxr-xr-x 1 tx tx   12  9月 22 08:11 gindex.7.ht2
-rwxr-xr-x 1 tx tx    8  9月 22 08:11 gindex.8.ht2
-rwxr-xr-x 1 tx tx 640M  9月 12 07:49 GWHBKHI00000000.genome.fasta
-rwxr-xr-x 1 tx tx 3.3K  9月 22 08:21 GWHBKHI00000000.genome.fasta.dict
-rwxr-xr-x 1 tx tx 1.1K  9月 22 09:59 GWHBKHI00000000.genome.fasta.fai
-rwxr-xr-x 1 tx tx 2.6G  9月 12 11:49 s1_1.fq.gz
-rwxr-xr-x 1 tx tx 2.5G  9月 12 11:49 s1_2.fq.gz
-rwxr-xr-x 1 tx tx 5.7G  9月 22 09:03 s1.bam
-rwxr-xr-x 1 tx tx 5.8G  9月 22 09:21 s1.markdup.bam
-rwxr-xr-x 1 tx tx  22G  9月 22 08:57 s1.sam
-rwxr-xr-x 1 tx tx 4.2G  9月 22 09:09 s1.sorted.bam
-rwxr-xr-x 1 tx tx 8.6K  9月 22 09:21 s1.sorted.markdup_metrics.txt    

Answer 1

以下回复参考：皆我百晓生、券券喵儿等免费微信小程序作答:

根据你提供的错误信息，问题似乎出在缺少 FASTA 的字典文件（.dict文件）。GATK 的 HaplotypeCaller 工具需要这个字典文件来正确地读取和解析你的参考基因组 FASTA 文件。这个字典文件通常是通过 CreateSequenceDictionary 工具从 FASTA 文件生成的。

在你的情况下，你已经有了一个 .fasta.dict 文件，但是它似乎没有被 GATK 正确识别。这可能是由于文件路径的问题，或者文件权限问题，或者文件本身的问题（例如格式不正确或者损坏）。

你可以尝试以下步骤来解决这个问题：

1. 确保你的 .dict 文件路径是正确的。在你的 GATK 命令中，你需要使用绝对路径来指定 .dict 文件。在你的错误信息中，看起来你正在尝试使用相对路径（file:///media/tx/...），这可能会导致问题。你可以尝试使用绝对路径来替换这个路径。例如，如果你的 .dict 文件实际上在 /home/tx/... 下，你需要修改你的 GATK 命令来反映这一点。
2. 检查文件的权限。确保你的 .dict 文件对运行 GATK 的用户是可读的。你可以使用 ls -l 命令来检查文件的权限。如果文件权限不正确，你可以使用 chmod 命令来更改它。
3. 检查 .dict 文件是否损坏或格式不正确。你可以尝试使用其他工具（如 samtools）来读取这个文件，看看是否有问题。如果文件损坏或格式不正确，你可能需要重新生成它。你可以使用 GATK 的 CreateSequenceDictionary 工具来从你的 .fasta 文件生成一个新的 .dict 文件。

如果以上步骤都不能解决问题，你可能需要更详细地检查你的环境设置和文件结构，或者寻求更专业的帮助。