SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
samtools: "sort by name" and fast random access cycomatto Bioinformatics 7 01-22-2013 03:42 AM
"allele balance ratio" and "quality by depth" in VCF files efoss Bioinformatics 2 10-25-2011 11:13 AM
Relatively large proportion of "LOWDATA", "FAIL" of FPKM_status running cufflink ruben6um Bioinformatics 3 10-12-2011 12:39 AM
The position file formats ".clocs" and "_pos.txt"? Ist there any difference? elgor Illumina/Solexa 0 06-27-2011 07:55 AM
"Systems biology and administration" & "Genome generation: no engineering allowed" seb567 Bioinformatics 0 05-25-2010 12:19 PM

Reply
 
Thread Tools
Old 03-23-2010, 12:59 PM   #1
popto
Junior Member
 
Location: London

Join Date: Mar 2010
Posts: 3
Default Haplotype and "random" chromosomes

Hi all!

I'm trying to make sense of some sequence data (50 base reads, Illumina) and I noticed an area of coverage gap on chromosome 6 - around the region that aligns very well to the haplotype chromosomes. I have some reads that mapped to the haplotype chromosomes (e.g. chr6_cox_hap1), not enough to explain the dip in coverage. I am worried that because of the high homology between the chromosomes, the "missing" reads might be "hiding" as ##:##:## (i.e., as not mappable due to the fact they map equally well to >1 locus).

So basically what I am wondering - and I apologize if this is a very basic question - Do you align your reads to all the available chromosomes or do you omit the "haplotype" and "random" ones from your build? And if you are using all of the chromosomes, do you observe the same dip in coverage?

I would be very grateful for any advice you might have...
Thanks!!
Popto
popto is offline   Reply With Quote
Old 03-23-2010, 01:32 PM   #2
Simon Anders
Senior Member
 
Location: Heidelberg, Germany

Join Date: Feb 2010
Posts: 994
Default

I always omit the haplotype sequences from the reference index, for precisely the reason you mention.

Simon
Simon Anders is offline   Reply With Quote
Old 03-23-2010, 01:47 PM   #3
popto
Junior Member
 
Location: London

Join Date: Mar 2010
Posts: 3
Default

Thank you, Simon, this is very helpful.
popto is offline   Reply With Quote
Old 03-23-2010, 02:23 PM   #4
thinkRNA
Member
 
Location: Carlsbad,CA

Join Date: Jan 2010
Posts: 94
Default

Quote:
Originally Posted by Simon Anders View Post
I always omit the haplotype sequences from the reference index, for precisely the reason you mention.

Simon
How do you determine which region is haplotype sequence?
thinkRNA is offline   Reply With Quote
Old 03-24-2010, 01:58 AM   #5
Simon Anders
Senior Member
 
Location: Heidelberg, Germany

Join Date: Feb 2010
Posts: 994
Default

I took my reference from Ensembl: ftp://ftp.ensembl.org/pub/current_fa...o_sapiens/dna/

All the files with "HSCHR" in the file name are haplotype variants, e.g., the "HSCHR6_MHC" files contain variants to the the MHC region of chromosome 6. I suggest to simply not include these files when building the reference (unless, of course, you are specifically interested in them, but then you need to do some additional tweaking).

The "nonchromosomal" file contains the "random" contigs. I usually include them, but these contigs are so short that it does not really matter.

Do not take, by the way, the repeat masked ("rm" in the filename) sequences. You should leave checking for repeats to the aligner.

Simon
Simon Anders is offline   Reply With Quote
Old 02-23-2011, 07:34 PM   #6
pcg
Junior Member
 
Location: USA

Join Date: Jan 2010
Posts: 8
Default

Simon,

I presume that if you do exclude the haplotypes in the index then you remove those chromosomes from the GTF annotation file aswell? Right?

So basically if I am understanding correctly the reason then Simon, you remove these haplotypes because there is going to be an alignment problem due to the high similarity between the two chromosomes and you may get false mapping to a chromosome?

Thanks,
pcg is offline   Reply With Quote
Old 02-23-2011, 11:58 PM   #7
Simon Anders
Senior Member
 
Location: Heidelberg, Germany

Join Date: Feb 2010
Posts: 994
Default

Quote:
Originally Posted by pcg View Post
I presume that if you do exclude the haplotypes in the index then you remove those chromosomes from the GTF annotation file aswell? Right?
Actually, no. The aligner does not need a GTF file, and when counting later (e.g. with my htseq-count script), a feature in the GTF file with a chromosome name that does not appear in the SAM file will not collect any counts anyway.

Quote:
So basically if I am understanding correctly the reason then Simon, you remove these haplotypes because there is going to be an alignment problem due to the high similarity between the two chromosomes and you may get false mapping to a chromosome?
Especially when looking for differential expression, it is a good idea to discount all non-unique alignments. Now, if the aligner sees several version of, e.g., the MHC, it does not know that these are all variants of the same region but rather treats them as paralogs at different places. So. if a read maps there, the aligner will think that there are multiple mappings, flag the read accordingly, and you will exclude it, ending up with no signal at all at the variant regions, even (or: especially) at the parts of the variant region that are actually conserved and would hence have posed no problem for mapping.

Simon
Simon Anders is offline   Reply With Quote
Old 02-24-2011, 11:55 AM   #8
pcg
Junior Member
 
Location: USA

Join Date: Jan 2010
Posts: 8
Default

Thanks Simon for your reply.

As you rightly point out you do not need a GTF for alignment but if you want to run a cufflinks analysis on the alignment and only want expression for what is currently annotated (in the GTF) then unless you remove those haplotypes from the GTF file you will still see hits to them and expression values?

Thanks in advance,
pcg is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 01:56 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO