Hello everyone,
I'm trying to figure out a way to determine whether, in cases where a single read (or read pair) maps across two or more SNP positions in a reference genome, whether or not the SNPs are in agreement within that read.
My data are mouse paired end RNA-seq reads from an F1 cross between the reference strain and an alternate strain. I've mapped the reads to the mm9 reference genome.
I also have a vcf file documenting the positions of all the SNPs between mm9 and the alternate genome.
So,for example, if a single read maps across two SNP positions, I want to figure out if the read matches the reference at both positions, the alternate genome at both positions, or the reference at one position, and the alternate genome at the other.
I started off using bedtools to get the intersection of the mapped reads and the SNPs, but I'm not sure if/how I can get from there to the information I'm after. I think I might be able to write a python script to figure it out from the files I have, but before I do that I want to make sure there isn't an existing tool to do this sort of thing. Any help or advice would be appreciated.
Thanks
I'm trying to figure out a way to determine whether, in cases where a single read (or read pair) maps across two or more SNP positions in a reference genome, whether or not the SNPs are in agreement within that read.
My data are mouse paired end RNA-seq reads from an F1 cross between the reference strain and an alternate strain. I've mapped the reads to the mm9 reference genome.
I also have a vcf file documenting the positions of all the SNPs between mm9 and the alternate genome.
So,for example, if a single read maps across two SNP positions, I want to figure out if the read matches the reference at both positions, the alternate genome at both positions, or the reference at one position, and the alternate genome at the other.
I started off using bedtools to get the intersection of the mapped reads and the SNPs, but I'm not sure if/how I can get from there to the information I'm after. I think I might be able to write a python script to figure it out from the files I have, but before I do that I want to make sure there isn't an existing tool to do this sort of thing. Any help or advice would be appreciated.
Thanks