Hello,
I have several BAC sequences each of which contains 2-7 contigs. I need to align these BACs to a pseudomolecule and scaffolds. Basically, I need to identify a region were BAC maps on the pseudomolecule. At the end I need to get statistics such as alignment identity, aln lenth, coverage, start end. I can get all of these information from BLAST megablast. However, since megablast is local, it results in many HSP's obviously. Do you guys have an idea how can I calculate this and find the region on a pseudomolecule where each of my BAC aligns? if yes, then how? I can use bwa as an alternative, however, I dont know how to extract information from .bam file..
1 answer
You may try to use LAST: http://last.cbrc.jp/
First create database from your assembly, then query it with your BAC fasta. You will get MAF output you can parse.
Caveat: repetitive contigs may be dropped from the output. But then it may be the correct behavior: Imagine that your BAC insert in reality is 100kb long and maps to scaffold/chromosome from assembly starting from position 1, with the last 10k (90-100kb BC & genomic) is full of repeats. On the top of it the next 50k of your assembly scaffold is also repetitive and full of stretches of Ns denoting gaps. You do not want to come to conclusion that your BAC is 150kb because something from BAC maps to these parts of scaffold.
Here we go (Shameless self-plug mode on):
http://openwetware.org/wiki/Wikiomics:WinterSchool_day2#last
One can do this MAF to BAM coversion not just with mapped short reads but also with contigs/BAC sequences or other genomes. And verify some tricky spots by eye if necessary.
Log in to answer this question.
-
Understanding Best-Hits filtering algorithm in BLAST+ applications
written by Denis
-
Converting BACs to genomic locations(coordinates) and vice versa
written by kmark.gr •
-
Pairwise Alignment; sequence identity and query coverage
written by fareehakanwal90 •
-
Blastn/megablast, obtain alignments of more than 100 bp and 99% identity
written by jlerga •
-
Better statergy for Gap Closing
written by second_exon
-
BAC sequencing issues
written by int11ap1
-
Finding repetitive elements from BAC sequences and Illumina sequences
written by int11ap1
-
What Would Be The Next Analysis Steps After Producing Alignments With Blastx
written by User000
-
Extracting Fasta Alignments From Parsed Blastxml File
written by Eric •
-
Genome Assembly With Physical Map
written by Plantae
did you try BLAT instead of BLAST ?
I don't know if it solves my problem..the output is more or less similar to what I get with megablast...Do you think taking the best hit solves the problem?
in the ucsc, the BAC end sequences are "placed on the assembled sequence using BLAT".
Yes but they use only the ends not the full sequence and they have a scoring system to filter hits.