Hi, this is quite a basic question I feel, but I'm not exactly sure of what the best practice would be. I am currently using IGV for an interactive view, but I would rather have this automated on the command line with any bespoke elements put in place by myself. What tools (command line, not interactive) would I need for the following:
I have indexed BAM data of specific reads along with their corresponding vcf file. I would like to go back all the way to my reference genome, so I can take each read with an identified SNP and pull out the sequence from the reference genome plus an extra N bp e.g., 50 bps, on either side from where the read aligns on to the reference. In short, I want to know what exists to the left and to the right of my aligned read. Could this be done using the #CHROM and POS data in the vcf?
Hi Anthony,
simply do a blast in the UCSC website (http://genome.ucsc.edu), the results will lead you to the sequences and also to a genome browser where you can zoom in or out, what you're searching for.
fred
Once you are acquainted with VCF files columns and headers, you can manipulate similar file in MS Excel as well.
If you are non-programmer then I would suggest to open your VCF files on MS Excel or any other similar programme. Identify SNPs for which you need flanking region, modify Start and End co-ordinate for them by adding and subtracting 50 and save this data as separate VCF file. Feed this VCF file to getfasta module of bedtools (http://bedtools.readthedocs.io/en/latest/content/tools/getfasta.html) and you would flanking nucleotide sequence.
- Badges
- Science method