fixing syndna filtering #16
Conversation
lucaspatel
left a comment
lucaspatel
left a comment
There was a problem hiding this comment.
Some room for improvement but fine for now.
| coverm filter --bam-files ${out_folder}/${sample_name}_GCF_000184185_sorted.sam --min-read-percent-identity 99.9 --min-read-aligned-percent 95 --threads {{nprocs}} -o ${out_folder}/${sample_name}_GCF_000184185.bam | ||
| samtools view -O SAM -o ${out_folder}/${sample_name}_no_GCF_000184185_sorted.sam ${out_folder}/${sample_name}_no_inserts.bam | ||
| minimap2 -x map-hifi -t {{nprocs}} -a --MD --eqx -o ${out_folder}/${sample_name}_GCF_000184185.sam ${db_folder}/GCF_000184185.1_ASM18418v1_genomic_chroso.fna ${out_folder}/${sample_name}_no_plasmid.fastq | ||
| samtools view -bS -@ {{ nprocs/2 | int }} ${out_folder}/${sample_name}_no_plasmid.fastq | samtools sort -@ {{ nprocs/2 | int }} -O bam -o ${out_folder}/${sample_name}_GCF_000184185_sorted.bam |
There was a problem hiding this comment.
This seems wrong. You're running samtools view on a FASTQ file. Shouldn't it be samtools view -bS -@ {{ nprocs/2 | int }} ${out_folder}/${sample_name}_GCF_000184185.sam?
There was a problem hiding this comment.
@jianshu93, can you comment? I'm not actually sure.
There was a problem hiding this comment.
should be the SAM file. I think I have the SAM as input in the example commands?
There was a problem hiding this comment.
This is what @jianshu93 sent:
minimap2 -x map-hifi -t 8 -a --MD --eqx -o reads.sam ecoli_genome.fna reads.fastq
samtools view -bS -@ 8 reads.fastq | samtools sort -@ 24 -O bam -o reads.sorted.bam
coverm filter --bam-files reads.sorted.bam --min-read-percent-identity 99.9 --min-read-aligned-percent 95 --threads 8 -o reads_filtered.sorted.bam
samtools view -O SAM -o reads_filtered.sam ./reads_filtered.sorted.bam
awk '{print $1}' reads_filtered.sam > reads_filtered.txt
seqkit grep -v -f reads_filtered.txt reads.fastq > reads_no_ecoli.fastq
@lucaspatel, could you also check?
There was a problem hiding this comment.
I still think is OK based on the commands sent by @jianshu93. Additionally, I have added as comments the original commands to make sure they are the same and remove checking them in the tests to avoid confusion or extra lines.
There was a problem hiding this comment.
Interesting, just tested it on barnacle2 and this really does seem to work:
head -n 100 15814.PSQ.0044_no_plasmid.fastq | samtools view -bS -@ 4 - | samtools sort -@ 4 -O sam
Good with me then
There was a problem hiding this comment.
In general this file is fine, but I think there's considerable room for optimization via piping. There are several SAM and BAM files that appear to be intermediates and do not really need to be written to disk unless they are needed for a downstream process.
There was a problem hiding this comment.
I think you are right; however, trying to simulate the current available commands and processes. Now, if you see some easy, low hanging fruits, please let me know so I can add them.
antgonza
left a comment
antgonza
left a comment
There was a problem hiding this comment.
@lucaspatel, thank you for your review.
There was a problem hiding this comment.
I think you are right; however, trying to simulate the current available commands and processes. Now, if you see some easy, low hanging fruits, please let me know so I can add them.
| coverm filter --bam-files ${out_folder}/${sample_name}_GCF_000184185_sorted.sam --min-read-percent-identity 99.9 --min-read-aligned-percent 95 --threads {{nprocs}} -o ${out_folder}/${sample_name}_GCF_000184185.bam | ||
| samtools view -O SAM -o ${out_folder}/${sample_name}_no_GCF_000184185_sorted.sam ${out_folder}/${sample_name}_no_inserts.bam | ||
| minimap2 -x map-hifi -t {{nprocs}} -a --MD --eqx -o ${out_folder}/${sample_name}_GCF_000184185.sam ${db_folder}/GCF_000184185.1_ASM18418v1_genomic_chroso.fna ${out_folder}/${sample_name}_no_plasmid.fastq | ||
| samtools view -bS -@ {{ nprocs/2 | int }} ${out_folder}/${sample_name}_no_plasmid.fastq | samtools sort -@ {{ nprocs/2 | int }} -O bam -o ${out_folder}/${sample_name}_GCF_000184185_sorted.bam |
There was a problem hiding this comment.
@jianshu93, can you comment? I'm not actually sure.
lucaspatel
left a comment
lucaspatel
left a comment
There was a problem hiding this comment.
Looks good to me, suggest some IO optimizations in the future but important to translate the working code as a start, which you've done very nicely
| coverm filter --bam-files ${out_folder}/${sample_name}_GCF_000184185_sorted.sam --min-read-percent-identity 99.9 --min-read-aligned-percent 95 --threads {{nprocs}} -o ${out_folder}/${sample_name}_GCF_000184185.bam | ||
| samtools view -O SAM -o ${out_folder}/${sample_name}_no_GCF_000184185_sorted.sam ${out_folder}/${sample_name}_no_inserts.bam | ||
| minimap2 -x map-hifi -t {{nprocs}} -a --MD --eqx -o ${out_folder}/${sample_name}_GCF_000184185.sam ${db_folder}/GCF_000184185.1_ASM18418v1_genomic_chroso.fna ${out_folder}/${sample_name}_no_plasmid.fastq | ||
| samtools view -bS -@ {{ nprocs/2 | int }} ${out_folder}/${sample_name}_no_plasmid.fastq | samtools sort -@ {{ nprocs/2 | int }} -O bam -o ${out_folder}/${sample_name}_GCF_000184185_sorted.bam |
There was a problem hiding this comment.
Interesting, just tested it on barnacle2 and this really does seem to work:
head -n 100 15814.PSQ.0044_no_plasmid.fastq | samtools view -bS -@ 4 - | samtools sort -@ 4 -O sam
Good with me then
No description provided.