-
Notifications
You must be signed in to change notification settings - Fork 17
Description
Hi,
Thanks for developing this nice tool.
python leviosam2.py -t 16 -i hifi_normal_minimap2.bam -o hifi_normal_minimap2_lifted_to_hg38 -C chm13v2-grch38.clft -f hg38.fa --sequence_type pb_hifi -a minimap2 --use_preset --lift_realign_config pacbio_all.yaml --keep_tmp_files
When running the full workflow, I encountered the following errors:
terminate called after throwing an instance of 'std::out_of_range'
[W::bam_hdr_read] EOF marker is absent. The input is probably truncated
After the workflow finished, the *-final.bam file contained only 38,311 reads, while the original BAM file had more than 10 million reads.
After running the leviosam2 lift command, I obtained unexpectedly low mapping results: total reads: 10,257,872 mapped: 638,305 (6.22%) unmapped: 9,619,567 (93.78%).
How should I address this issue, or is the result simply expected? I aligned the COLO829 HiFi reads to T2T-CHM13 and then used leviosam2 to convert the BAM coordinates to hg38.
log:
[2025-10-19 19:51:35 /home/ljz/haplotype/data/benchmark/script/leviosam2/workflow/leviosam2.py:940 INFO] `--target_aligner_index` is empty. Set to hg38.fa
[2025-10-19 19:51:35 /home/ljz/haplotype/data/benchmark/script/leviosam2/workflow/leviosam2.py:385 INFO] leviosam2 lift -C chm13v2-grch38.clft -O bam -a hifi_normal_minimap2.bam -p hifi_normal_minimap2_lifted_to_hg38 -t 16 -m -f hg38.fa -S hdist:100 -G 1000 -x /home/ljz/haplotype/data/benchmark/script/leviosam2/configs/pacbio_all.yaml
[I::add_split_rule] Adding rule `hdist:100`
[I::aln::AlnOpts::print_parameters] engine = ksw_extd2_sse
[I::aln::AlnOpts::print_parameters] flag = 8
[I::aln::AlnOpts::print_parameters] nm_threshold = 0
[I::aln::AlnOpts::print_parameters] a = 0
[I::aln::AlnOpts::print_parameters] b = 4
[I::aln::AlnOpts::print_parameters] q = 6
[I::aln::AlnOpts::print_parameters] e = 2
[I::aln::AlnOpts::print_parameters] q2 = 26
[I::aln::AlnOpts::print_parameters] e2 = 1
[I::aln::AlnOpts::print_parameters] w = 601
[I::aln::AlnOpts::print_parameters] zdrop = 400
[I::aln::AlnOpts::print_parameters] end_bonus = -1
[I::update_thread_allocation] --theads=16, --lift_threads=12, --hts_threads=4
[I::lift_run] Loading levioSAM index...done
[I::load_fasta] Loading FASTA...done
[I::WriteDeferred::print_info] Alignments with any of the below features are deferred:
- unlifted
- NM:i (post-liftover) > 100
terminate called after throwing an instance of 'std::out_of_range'
what(): map::at
[2025-10-19 19:52:28 /home/ljz/haplotype/data/benchmark/script/leviosam2/workflow/leviosam2.py:412 INFO] samtools sort -@ 16 -o hifi_normal_minimap2_lifted_to_hg38-committed-sorted.bam hifi_normal_minimap2_lifted_to_hg38-committed.bam
[W::bam_hdr_read] EOF marker is absent. The input is probably truncated
[bam_sort_core] merging from 0 files and 16 in-memory blocks...
[2025-10-19 19:52:31 /home/ljz/haplotype/data/benchmark/script/leviosam2/workflow/leviosam2.py:768 INFO] samtools fastq hifi_normal_minimap2_lifted_to_hg38-deferred.bam | bgzip > hifi_normal_minimap2_lifted_to_hg38-deferred.fq.gz
[W::bam_hdr_read] EOF marker is absent. The input is probably truncated
[M::bam2fq_mainloop] discarded 0 singletons
[M::bam2fq_mainloop] processed 14092 reads
[2025-10-19 19:53:06 /home/ljz/haplotype/data/benchmark/script/leviosam2/workflow/leviosam2.py:662 INFO] minimap2 -a --MD -t 13 hg38.fa hifi_normal_minimap2_lifted_to_hg38-deferred.fq.gz | samtools sort -@ 3 -o hifi_normal_minimap2_lifted_to_hg38-realigned.bam
[M::mm_idx_gen::37.705*1.57] collected minimizers
[M::mm_idx_gen::41.238*2.42] sorted minimizers
[M::main::41.238*2.42] loaded/built the index for 194 target sequence(s)
[M::mm_mapopt_update::42.528*2.38] mid_occ = 706
[M::mm_idx_stat] kmer size: 15; skip: 10; is_hpc: 0; #seq: 194
[M::mm_idx_stat::43.479*2.35] distinct minimizers: 100159079 (38.75% are singletons); average occurrences: 5.545; average spacing: 5.581; total length: 3099750718
[M::worker_pipeline::120.891*8.56] mapped 14092 sequences
[M::main] Version: 2.28-r1209
[M::main] CMD: minimap2 -a --MD -t 13 hg38.fa hifi_normal_minimap2_lifted_to_hg38-deferred.fq.gz
[M::main] Real time: 121.046 sec; CPU: 1035.020 sec; Peak RSS: 13.485 GB
[bam_sort_core] merging from 0 files and 3 in-memory blocks...
Best regards.