J'ai un souci assez agaçant avec haplotypecaller. Je le lance par interval (avec une boucle de srun) et je vais un sbatch array pour faire plusieurs bam en même temps.
Ceci étant dit, il arrive parfois pour une raison qui m'échappe que le job échoue sur un des srun, d'apres le log (voir ci après) il semble qu'il n'arrive pas à trouver le vcf qu'il est censé lui même créer... D'après gatk ça viendrait plutôt du cluster ou du système de fichiers. Il suffit ensuite de relancer le script avec l'array qui a échoué, sans rien changer, et ça marche. Mais du coup je perd du temps (sachant que les autres srun fonctionnent, le job ne vas échouer que lorsque tous les intervals seront finis)
htsjdk.samtools.util.RuntimeIOException: File not found: /shared/projects/gentaumix/processed_data_set5/05_VCF/A_gVCF/G632_DA_C002A5Q_1_H5N3JDSX2.DUAL2.trimmed.align.filtered.dedup.fixed.recal.bam/temp_gVCF/0018.G632_DA_C002A5Q_1_H5N3JDSX2.DUAL2.trimmed.align.filtered.dedup.fixed.recal.g.vcf.gz
....
Caused by: java.nio.file.FileSystemException: /shared/projects/gentaumix/processed_data_set5/05_VCF/A_gVCF/G632_DA_C002A5Q_1_H5N3JDSX2.DUAL2.trimmed.align.filtered.dedup.fixed.recal.bam/temp_gVCF/0018.G632_DA_C002A5Q_1_H5N3JDSX2.DUAL2.trimmed.align.filtered.dedup.fixed.recal.g.vcf.gz: Erreur de protocole
at sun.nio.fs.UnixException.translateToIOException(UnixException.java:91)
at sun.nio.fs.UnixException.rethrowAsIOException(UnixException.java:102)
at sun.nio.fs.UnixException.rethrowAsIOException(UnixException.java:107)
at sun.nio.fs.UnixFileSystemProvider.newByteChannel(UnixFileSystemProvider.java:214)
at java.nio.file.spi.FileSystemProvider.newOutputStream(FileSystemProvider.java:434)
at java.nio.file.Files.newOutputStream(Files.java:216)
at htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder.build(VariantContextWriterBuilder.java:447)
... 11 more
srun: error: cpu-node-61: task 0: Exited with exit code 3
Bien sur. Petite précision sur le script : j'y ai mis des points de contrôle justement parceque parfois un interval echoue, et qu'à la fin ca merge le tout donc si 49 intervals sur 50 fonctionnent j'obtiendrai bien un gvcf à la fin si j'avais pas de controle sans savoir qu'il y a eu un soucis.
#!/usr/bin/env bash
#SBATCH --mail-user='chartreq@igbmc.fr'
#SBATCH --mail-type=END,ARRAY_TASKS
#SBATCH --job-name='Gentaumix_haplotypecaller'
#SBATCH --error=/shared/projects/gentaumix/logs/haplotypecaller/haplotypecaller_%A_%a.err
#SBATCH --output=/shared/projects/gentaumix/logs/haplotypecaller/haplotypecaller_%A_%a.log
#SBATCH --mem-per-cpu=5000
#SBATCH --cpus-per-task=1
#SBATCH --array=1-74
#SBATCH --ntasks=50
#SBATCH -A gentaumix
#Path to analysis ready bam
BAM_INPUT_DIR=/shared/projects/gentaumix/processed_data/03_Alignment/D_recal
#Path to the output gVCF directory:
gVCF_OUTPUT_DIR=/shared/projects/gentaumix/processed_data/05_VCF/A_gVCF
#Path of the reference genome:
REF_Genome=/shared/projects/gentaumix/Ressources/grch38_BWA_2/GCA_000001405.15_GRCh38_no_alt_plus_hs38d1_analysis_set.fa
#Path to file containing the list of intervals:
List_of_interval=/shared/projects/gentaumix/Ressources/wgs_calling_regions.hg38.interval_list
#Number of intervals to be processed in parallel (must be equal to sbatch --ntasks)
SCATTER_COUNT=50
##################### NO CHANGE NEED AFTER #########################################################################
#get the basename of the BAM file corresponding to the number of the job array
BAM_INPUT=$(basename $(ls ${BAM_INPUT_DIR}/*.bam | awk "NR==${SLURM_ARRAY_TASK_ID}"))
#Creation of a temporary directory to store all temporary files of the jobs
mkdir -p ${gVCF_OUTPUT_DIR}/${BAM_INPUT}
#################### List and creation of logs dir for all the command of the jobs #################################
#For haplotypecaller
logs_HC=${gVCF_OUTPUT_DIR}/logs_HC/${SLURM_ARRAY_TASK_ID}
mkdir -p ${logs_HC}
#For IntervalListTools
logs_IntervalListTools=${gVCF_OUTPUT_DIR}/logs_IntervalListTools/${SLURM_ARRAY_TASK_ID}
mkdir -p ${logs_IntervalListTools}
#For MergeVcfs
logs_MergeVcfs=${gVCF_OUTPUT_DIR}/logs_MergeVcfs/${SLURM_ARRAY_TASK_ID}
mkdir -p ${logs_MergeVcfs}
#For check if the number of temporary gVCF is ok:
failed_log_DIR=${gVCF_OUTPUT_DIR}/failled_logs
mkdir -p ${failed_log_DIR}
#####################################################################################################################
module load conda
source activate /shared/projects/gentaumix/conda/envs/gatk_4.2.2.0
######################################################## Split intervals #############################################
#Definition and creation of directory to store intervals:
Scattered_DIR=${gVCF_OUTPUT_DIR}/${BAM_INPUT}/interval
mkdir -p ${Scattered_DIR}
#Main command to split interval:
gatk \
IntervalListTools \
-INPUT ${List_of_interval} \
-OUTPUT ${Scattered_DIR} \
-UNIQUE true \
-SORT true \
-SCATTER_COUNT ${SCATTER_COUNT} \
-BREAK_BANDS_AT_MULTIPLES_OF 1000000 \
-SUBDIVISION_MODE BALANCING_WITHOUT_INTERVAL_SUBDIVISION_WITH_OVERFLOW \
-PADDING 0 \
-ACTION CONCAT \
-INCLUDE_FILTERED false \
-INVERT false \
-VERBOSITY INFO \
-QUIET false \
-VALIDATION_STRINGENCY STRICT \
-COMPRESSION_LEVEL 5 \
-MAX_RECORDS_IN_RAM 500000 \
-CREATE_INDEX false \
-CREATE_MD5_FILE false \
2> ${logs_IntervalListTools}/${BAM_INPUT}.log
wait
#######################################################################################################################
#Creation of the basename of gVCF_output:
GVCF_OUTPUT=${BAM_INPUT/.bam/.g.vcf.gz}
#Creation of temporary output folder for gvcf:
temp_gVCF_OUTPUT_DIR=${gVCF_OUTPUT_DIR}/${BAM_INPUT}/temp_gVCF
mkdir -p ${temp_gVCF_OUTPUT_DIR}
###################################################### HaplotypeCaller ################################################
for i in $(seq -w 0001 00${SCATTER_COUNT})
do
srun --ntasks=1 gatk --java-options "-Xmx${SLURM_MEM_PER_CPU}M" HaplotypeCaller \
-R ${REF_Genome} \
-L ${Scattered_DIR}/temp_${i}_of_${SCATTER_COUNT}/scattered.interval_list \
-I ${BAM_INPUT_DIR}/${BAM_INPUT} \
-O ${temp_gVCF_OUTPUT_DIR}/${i}.${GVCF_OUTPUT} \
-G StandardAnnotation -G AS_StandardAnnotation -G StandardHCAnnotation \
-GQB 10 -GQB 20 -GQB 30 -GQB 40 -GQB 50 -GQB 60 -GQB 70 -GQB 80 -GQB 90 \
-ERC GVCF \
--pcr-indel-model NONE \
2> ${logs_HC}/${i}.${GVCF_OUTPUT}.log &
done
wait
###########################################################################################################################
############################## Check if all intermediate gvcf have been created ###########################################
#The folder should have 2 files per interval so:
number_of_expected_file=$(( 2 * ${SCATTER_COUNT} ))
if (( $(ls ${temp_gVCF_OUTPUT_DIR} | wc -l ) != ${number_of_expected_file} ))
then echo "The file ${BAM_INPUT} run by array number ${SLURM_ARRAY_TASK_ID} had failed" > ${failed_log_DIR}/${BAM_INPUT}.failed.logs && exit 1
fi
wait
###########################################################################################################################
########################################## MergeVcfs ######################################################################
# Create a list with all temporary gvcf
for i in $(seq -w 0001 00${SCATTER_COUNT})
do
echo "${temp_gVCF_OUTPUT_DIR}/${i}.${GVCF_OUTPUT}" >>${gVCF_OUTPUT_DIR}/${BAM_INPUT}/${GVCF_OUTPUT}.list
done
wait
srun --ntasks=1 gatk --java-options "-Xmx${SLURM_MEM_PER_CPU}M" MergeVcfs \
-I ${gVCF_OUTPUT_DIR}/${BAM_INPUT}/${GVCF_OUTPUT}.list \
-O ${gVCF_OUTPUT_DIR}/${GVCF_OUTPUT} \
2> ${logs_MergeVcfs}/${GVCF_OUTPUT}.log
if [ $? -ne 0 ]
then echo "The file ${BAM_INPUT} run by array number ${SLURM_ARRAY_TASK_ID} failed" > ${failed_log_DIR}/${BAM_INPUT}.merged.failed.logs && exit 1
fi
wait
#############################################################################################################################
rm -r ${gVCF_OUTPUT_DIR}/${BAM_INPUT}
Celui ci a du échouer (mais ca date du 21, j'en ai pas de plus récent, si besoin la prochaine fois que ca m'arrive je le met tout de suite ici) :
19453518_9
Ce n'est pas évident de diagnostiquer précisément. Les données d'entrées ne sont plus présentes sur ton espace projet (ou pas en l'état).
Je veux bien que tu nous alertes immédiatement si le problème se reproduit sans nettoyer les données d'entrée. Nous pourrons essayer d'étudier de plus près s'il y a un soucis côté système de fichier.
voila ca vient de se produire avec le job 19531903_1.
Le job est encore en cours (pour les autres interval mais je vais le tuer pour eviter de prendre des ressources inutilement).
log
[27 octobre 2021 09:52:33 CEST] org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller done. Elapsed time: 0.23 minutes.
Runtime.totalMemory()=2144337920
htsjdk.samtools.util.RuntimeIOException: File not found: /shared/projects/gentaumix/HG001/dragmap_gatk/05_VCF/B_gVCF/HG001.trimmed.align.filtered.dedup.bam/temp_gVCF/0021.HG001.trimmed.align.filtered.dedup.g.vcf.gz
at htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder.build(VariantContextWriterBuilder.java:451)
at htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder.build(VariantContextWriterBuilder.java:415)
at org.broadinstitute.hellbender.utils.variant.GATKVariantContextUtils.createVCFWriter(GATKVariantContextUtils.java:123)
at org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCallerEngine.makeVCFWriter(HaplotypeCallerEngine.java:374)
at org.broadinstitute.hellbender.tools.walkers.haplotypecaller.HaplotypeCaller.onTraversalStart(HaplotypeCaller.java:263)
at org.broadinstitute.hellbender.engine.GATKTool.doWork(GATKTool.java:1083)
at org.broadinstitute.hellbender.cmdline.CommandLineProgram.runTool(CommandLineProgram.java:140)
at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMainPostParseArgs(CommandLineProgram.java:192)
at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMain(CommandLineProgram.java:211)
at org.broadinstitute.hellbender.Main.runCommandLineProgram(Main.java:160)
at org.broadinstitute.hellbender.Main.mainEntry(Main.java:203)
at org.broadinstitute.hellbender.Main.main(Main.java:289)
Caused by: java.nio.file.FileSystemException: /shared/projects/gentaumix/HG001/dragmap_gatk/05_VCF/B_gVCF/HG001.trimmed.align.filtered.dedup.bam/temp_gVCF/0021.HG001.trimmed.align.filtered.dedup.g.vcf.gz: Erreur de protocole
at sun.nio.fs.UnixException.translateToIOException(UnixException.java:91)
at sun.nio.fs.UnixException.rethrowAsIOException(UnixException.java:102)
at sun.nio.fs.UnixException.rethrowAsIOException(UnixException.java:107)
at sun.nio.fs.UnixFileSystemProvider.newByteChannel(UnixFileSystemProvider.java:214)
at java.nio.file.spi.FileSystemProvider.newOutputStream(FileSystemProvider.java:434)
at java.nio.file.Files.newOutputStream(Files.java:216)
at htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder.build(VariantContextWriterBuilder.java:447)
... 11 more
srun: error: cpu-node-74: task 0: Exited with exit code 3
Est-ce normal que le fichier /shared/projects/gentaumix/HG001/dragmap_gatk/05_VCF/B_gVCF/HG001.trimmed.align.filtered.dedup.bam/temp_gVCF/0021.HG001.trimmed.align.filtered.dedup.g.vcf.gz soit vide ?
Non. Enfin oui et non, disont qu'haplotypecaller est censé le créer et le remplir au fur et à mesure que le job avance, mais la il a échoué à avant même de commencer.
Ce que j'ai du mal à comprendre c'est que le log sous entend qu'il ne trouve pas ce fichier (qu'il a pourtant lui même créér) alors qu'il est bien là...