Removing the Mitochondrial genome from contigs

[BenAawf]

Hi,

Thanks for sharing this cool work!

I'm having trouble using gemma_los.py to remove mitochondrial contigs from my HiFiasm assembly. I adapted the provided bash script for my animal genome, but the gemma_los.py part doesn't seem to be working.

Could you clarify what gemma_los.py does with the PAF file? My understanding is that mitochondrion.txt should list coordinates (like a BED file) so seqkit grep -v -f can remove the mitochondrial sequences, leaving only nuclear contigs. Is this correct?

Here's the command I'm using:

#!/bin/bash . /home/abdul/mambaforge/etc/profile.d/conda.sh conda activate notebooklab mt=$1 ref=$2 threads=58 minimap2 -t $threads -x asm5 $mt $ref> mitochondrion.paf python3 gemma_los.py mitochondrion.paf > mitochondrion.txt conda deactivate conda activate seqkit_env seqkit grep -v -f mitochondrion.txt $ref > Hifiasm_remove_mt.fasta

[liushoucheng]

try:
awk '{print $1}' mito_vs_assembly.paf | sort | uniq > mito_contigs.txt
seqkit grep -v -f mito_contigs.txt assembly.fasta > nuclear_only.fasta

or use tools like FCS-GX kraken2

[lixp0528]

@BenAawf
Hello, we have uploaded the following files: mitochondrion.paf, chloroplast.paf, mitochondrion.txt, and chloroplast.txt. The gemma_los.py script is a custom tool we developed to filter contigs based on their alignment with mitochondrial and chloroplast genome sequences. Specifically, this script removes contigs where the alignment length exceeds 50% of the total contig length. The absence of generated .txt files may indicate either issues with the input file format or the lack of contigs meeting the specified filtering criteria.

[BenAawf]

The above awk removes the entire contigs where the mito had a blast hit, but I just want to remove the sequence that belongsto the mitochondrial my PAF output looks like this ptg000005l 39299640 12656062 12656122 - ptg000001l_rc_rotated 16694 9701 9761 60 60 27 tp:A:P cm:i:5 s1:i:60 s2:i:0 dv:f:0.0176 rl:i:0 ptg000006l 195265338 45886461 45886760 - ptg000001l_rc_rotated 16694 5069 5373 57 304 13 tp:A:P cm:i:3 s1:i:55 s2:i:0 dv:f:0.0265 rl:i:0 ptg000013l 94777610 15720152 15721044 + ptg000001l_rc_rotated 16694 13395 14298 108 903 60 tp:A:P cm:i:7 s1:i:105 s2:i:0 dv:f:0.0131 rl:i:0

[liushoucheng]

The above awk removes the entire contigs where the mito had a blast hit, but I just want to remove the sequence that belongsto the mitochondrial my PAF output looks like this ptg000005l 39299640 12656062 12656122 - ptg000001l_rc_rotated 16694 9701 9761 60 60 27 tp:A:P cm:i:5 s1:i:60 s2:i:0 dv:f:0.0176 rl:i:0 ptg000006l 195265338 45886461 45886760 - ptg000001l_rc_rotated 16694 5069 5373 57 304 13 tp:A:P cm:i:3 s1:i:55 s2:i:0 dv:f:0.0265 rl:i:0 ptg000013l 94777610 15720152 15721044 + ptg000001l_rc_rotated 16694 13395 14298 108 903 60 tp:A:P cm:i:7 s1:i:105 s2:i:0 dv:f:0.0131 rl:i:0

we updated the script,you can download then
try : python gemma_los.py mitochondrion.paf --length 100 > mitochondrion.txt
--length is length threshold alignment length exceeds 100bp,you can control depend your case.