Samtools mpileup for RNA Editing levels generating empty output file

Question

[haydenshinn]

I am working on a bash script that takes in bam files (respective paths listed out as required by mpileup documentation), a reference fasta (that bams were aligned to) and a tab separated "known sites" tsv that contains known editing sites for RNA seq data (chromosome, position 1, position 2).

My script runs with no errors but is giving an output file with size 0. I have looked at other examples, made sure to sort my bams, the positions list, and index the reference file but still have no luck. Has anyone encountered a similar issue or understand what's going on under the hood with mpileup where it would produce a size 0 output? Thanks.

# Samtools mpileup command to get RNA A to I editing levels with reference and known sites, followed by listing the BAM files
samtools mpileup -f ${reference} -l ${known_sites} -b ${bam_list} > ${output_dir}/mpileup_output.txt

[Pierre Lindenbaum]

check known_sites have the same chromosome notation than reference ('chr1' vs '1'). Also the position file MUST be a BED file 0-based, half-open .

[Juke34]

I have no clue about your problem sorry. But I wanted to mention that it is possible to avoid to pass by the extra step of mpileup. As a side project I have written SAPIN a python tool to look at editing sites directly from bam without passing by the mpileup format. It might interest you. I do not filter by specific position but this is something that can be implemented easily. Or this filter can be donne as a downstream step.

[haydenshinn]

Thanks, I'll check it out