VCF проблема с использованием Mitty
Я смоделировал файлы VCF, используя Mitty "simulate-варианты", и теперь я пытаюсь смоделировать чтения (60X и средняя длина 50, если это возможно), снова используя Mitty. Вот мой скрипт, который я использую:
#!/usr/bin/env bash
FASTA=/fast/users/me/data/human_g1k_v37.fasta
SAMPLEVCF=/fast/users/me/sim_VCF/VCF1/simtest.vcf.gz
SAMPLENAME=sample1
REGION_BED=/fast/users/me/bedfile/all_chr1.bed
FILTVCF=sample1-filt.vcf.gz
READMODEL=old-Garvan.pkl
COVERAGE=60
READ_GEN_SEED=7
FASTQ_PREFIX=sample1-reads
READ_CORRUPT_SEED=8
mitty -v4 filter-variants \
${SAMPLEVCF} \
${SAMPLENAME} \
${REGION_BED} \
- \
2> vcf-filter.log | bgzip -c > ${FILTVCF}
tabix -p vcf ${FILTVCF}
mitty -v4 generate-reads \
${FASTA} \
${FILTVCF} \
${SAMPLENAME} \
${REGION_BED} \
${READMODEL} \
${COVERAGE} \
${READ_GEN_SEED} \
>(gzip > ${FASTQ_PREFIX}1.fq.gz) \
${FASTQ_PREFIX}_lq.txt \
--fastq2 >(gzip > ${FASTQ_PREFIX}2.fq.gz) \
--threads 2
mitty -v4 corrupt-reads \
${READMODEL} \
${FASTQ_PREFIX}1.fq.gz >(gzip > ${FASTQ_PREFIX}_corrupt1.fq.gz) \
${FASTQ_PREFIX}_lq.txt \
${FASTQ_PREFIX}_corrupt-lq.txt \
${READ_CORRUPT_SEED} \
--fastq2-in ${FASTQ_PREFIX}2.fq.gz \
--fastq2-out >(gzip > ${FASTQ_PREFIX}_corrupt2.fq.gz) \
--threads 2
И вот вывод ошибок множественного числа, которые я получаю, когда я пытаюсь запустить это:
DEBUG:mitty.cli:Mitty version 2.28.3
DEBUG:root:Found model old-Garvan.pkl in builtins
Traceback (most recent call last):
File "/home/me/.local/bin/mitty", line 11, in <module>
load_entry_point('mitty==2.28.3', 'console_scripts', 'mitty')()
File "/home/me/.local/lib/python3.5/site-packages/click/core.py", line 764, in __call__
return self.main(*args, **kwargs)
File "/home/me/.local/lib/python3.5/site-packages/click/core.py", line 717, in main
rv = self.invoke(ctx)
File "/home/me/.local/lib/python3.5/site-packages/click/core.py", line 1137, in invoke
return _process_result(sub_ctx.command.invoke(sub_ctx))
File "/home/me/.local/lib/python3.5/site-packages/click/core.py", line 956, in invoke
return ctx.invoke(self.callback, **ctx.params)
File "/home/me/.local/lib/python3.5/site-packages/click/core.py", line 555, in invoke
return callback(*args, **kwargs)
File "/home/me/.local/lib/python3.5/site-packages/mitty/cli.py", line 162, in generate_reads
threads=threads, seed=seed)
File "/home/me/.local/lib/python3.5/site-packages/mitty/simulation/readgenerate.py", line 74, in process_multi_threaded
vcf_df = vio.load_variant_file(vcf_fname, sample_name, bed_fname)
File "/home/me/.local/lib/python3.5/site-packages/mitty/lib/vcfio.py", line 47, in load_variant_file
vcf_fp = pysam.VariantFile(fname, mode)
File "pysam/libcbcf.pyx", line 3479, in pysam.libcbcf.VariantFile.__init__ (pysam/libcbcf.c:51398)
File "pysam/libcbcf.pyx", line 3669, in pysam.libcbcf.VariantFile.open (pysam/libcbcf.c:54088)
ValueError: invalid file `b'sample1-filt.vcf.gz'` (mode='b'r'') - is it VCF/BCF format?
DEBUG:mitty.cli:Mitty version 2.28.3
Usage: mitty corrupt-reads [OPTIONS] MODELFILE FASTQ1_IN FASTQ1_OUT SIDECAR_IN
SIDECAR_OUT SEED
Try "mitty corrupt-reads --help" for help.
Error: Invalid value for "SIDECAR_IN": Path "sample1-reads_lq.txt" does not exist.
Я думаю, у меня есть, как минимум, две ошибки в моем коде. Когда я избавляюсь от "части с поврежденным чтением", у меня все еще остается первая проблема, и я не знаю, что делать с этим файлом VCF.
Кто-нибудь знает, почему я их получаю?
Большое спасибо, ура