Seqanswers Leaderboard Ad

Collapse

Announcement

Collapse
No announcement yet.
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • Mismatch between BAM and GFF file

    I am trying to use htseq-count to count the occurences of features (taken from my GFF file) in my BAM file:

    I use the std htseq-count command viz;
    Code:
    samtools view file.bam | htseq-count [options] - file.gff
    I get a series of error messages telling me that my reads are all skipped as they don't occur in my GFF files:

    Code:
    Warning: Skipping read 'HWI-ST1085:118:C1ALWACXX:1:2213:19460:44494', because chromosome 'chromosome:AGPv2:1:1:301354135:1', to which it has been aligned, did not appear in the GFF file'.
    The following is the first 10 lines of the BAM file:
    Code:
    HWI-ST1085:118:C1ALWACXX:1:1105:4877:50025	0	chromosome:AGPv2:1:1:301354135:1	3412	50	100M	*	0	0	TGAAGTATAAGGCAACCCAAGTCTGCCATCATCTCTTTCTCGTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAG	&&&(((((*****++++++++*++++++++++++++++++++++++++++++++++++++++++++++++++****((((((''&''&&&&&&&&&&$&$	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1104:1183:70980	0	chromosome:AGPv2:1:1:301354135:1	3444	50	100M	*	0	0	CTCTTTCTCGTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAG	&&&((((&*)***)+++++++++++*+++))++&*+)(&$)++)+))++*++&*+)+*+++++++&***)*((((&%%%&&&&&&&&&&$%"$!""&&&&	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1111:2149:77593	0	chromosome:AGPv2:1:1:301354135:1	3444	50	100M	*	0	0	CTCTTTCTCGTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAG	&&&(((((*****++++++++++++++++++++++++++++++++++++++++++++++++++++******((((&&&&&&&&&&&&&&&&&&$%&&&'&	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:2313:11306:75258	0	chromosome:AGPv2:1:1:301354135:1	3446	50	100M	*	0	0	CTTTCTCGTTGAGCTGAAGCGTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTC	$$$&&&&&&&&&!!"''+(+!!$&''++!$'&'++++'+++++'+++'+++++++'&++++&&&&&&&&&&%%$$$$$$$%$$$%$%#$$$###%%&%%#	AS:i:-4	XN:i:0	XM:i:2	XO:i:0	XG:i:0	NM:i:2	MD:Z:12A7C79	YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1214:18295:59994	0	chromosome:AGPv2:1:1:301354135:1	3450	50	100M	*	0	0	CTCGTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTCCGCT	$$$&&&&%(%(((%"&"%"&"('"!"'&)(+!$&%(("&""((&!%))+'(%!$((""#!!#&!$!!"%&$""!!!"!""#!#$""""%"%#%%$!!!!!	AS:i:-2	XN:i:0	XM:i:1	XO:i:0	XG:i:0	NM:i:1	MD:Z:96GYT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1312:13622:34439	0	chromosome:AGPv2:1:1:301354135:1	3453	50	100M	*	0	0	GTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTCGGCTAGT	&&&(((((****)++++++++++++++++++++++++++++++++++++++++++++)++++++++***(((&&&&&&&&&&&&$%&&&'&&&&&&&&&$	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:2112:2965:72104	0	chromosome:AGPv2:1:1:301354135:1	3453	50	100M	*	0	0	GTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTCGGCTAGT	&%&(((((*****+++++++++++++++++++++++*+++++++++++++++++++++++++++++***(((&&&&&&&&&&&&#%%&&'&&&&&&$%&"	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:2314:17254:77725	0	chromosome:AGPv2:1:1:301354135:1	3453	50	100M	*	0	0	GTTGAACTGAAGCCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTCGGCTGGT	&&&(((((*****+++++++++++++++++++++++)+++++++++++++++++++))++++*+++***(((&&&&&&&&&&&&%%&&&&&&&&&&%!!!	AS:i:-2	XN:i:0	XM:i:1	XO:i:0	XG:i:0	NM:i:1	MD:Z:97AYT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1212:1945:32241	0	chromosome:AGPv2:1:1:301354135:1	3465	50	100M	*	0	0	CCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCGGCGGTAGCTAGTCGGCTAGTCCATTGACTGGC	%&&((&&&*(*))++++*+++++++++++)+++'++')++)+)+%&(*'*+++)++))'***(%&''&"$#$!#%&&&&&&$&&&$%#&&"%&&&&#!"%	AS:i:-2	XN:i:0	XM:i:1	XO:i:0	XG:i:0	NM:i:1	MD:Z:67A32	YT:Z:UU	NH:i:1
    HWI-ST1085:118:C1ALWACXX:1:1109:10979:8471	16	chromosome:AGPv2:1:1:301354135:1	3465	50	100M	*	0	0	CCTGTCCATGCCTCCATGGCCCAGTCCAGCATCATCGCCAATCAGAGCTGAGGGCAGCCGCAGAGCCAGCGGTAGCTAGTCGGCTAGTCCATTGACTGGC	&&&&&&&&&&&&&&&&&&%&&&&&&&&&&&(&&&&&&''''''(((((*****++++++++++++++++++++++++++++++++++*****(((((&&&	AS:i:0	XN:i:0	XM:i:0	XO:i:0	XG:i:0	NM:i:0	MD:Z:100YT:Z:UU	NH:i:1
    10 lines of GFF file:

    Code:
    1	ensembl	chromosome	1	301354135	.	.	.	ID=1;Name=chromosome:AGPv2:1:1:301354135:1
    1	ensembl	gene	4854	9652	.	-	.	ID=GRMZM2G059865;Name=GRMZM2G059865;biotype=protein_coding
    1	ensembl	mRNA	4854	9652	.	-	.	ID=GRMZM2G059865_T01;Parent=GRMZM2G059865;Name=GRMZM2G059865_T01;biotype=protein_coding
    1	ensembl	intron	7904	9192	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71462
    1	ensembl	intron	7121	7593	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71463
    1	ensembl	intron	6798	6917	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71464
    1	ensembl	intron	6518	6638	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71465
    1	ensembl	intron	6266	6361	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71466
    1	ensembl	intron	5976	6107	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71467
    1	ensembl	intron	5408	5856	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71468
    1	ensembl	intron	5189	5341	.	-	.	Parent=GRMZM2G059865_T01;Name=intron.71469
    All ID attributes in my GFF file seem to match the BAM headers. Not sure if my untrained eye is missing something. Any suggestions? How should I modify the GFF file?

    Thanks
    Siva

  • #2
    Originally posted by Siva View Post
    All ID attributes in my GFF file seem to match the BAM headers. Not sure if my untrained eye is missing something. Any suggestions? How should I modify the GFF file?
    No, they don't match. The ID for the chromosome in your GFF file is '1'; column 1 of GFF file.

    Comment


    • #3
      Can you also check, if the reference fasta you created to align the bam file has same chromosome names as that of GFF.

      If you see only few reads are ignored the way you showed them here, its possible, your reference fasta has more chromosomes than in gff. For example, 'chromosome:AGPv2:1:1:301354135:1"

      Comment


      • #4
        Originally posted by kmcarr View Post
        No, they don't match. The ID for the chromosome in your GFF file is '1'; column 1 of GFF file.
        Whew! Thanks I changed the ID ("1")in the GFF to the rather unwieldy string ("chromosome:AGPv2:1:1:301354135:1") as given in BAM and it worked!

        Siva

        Comment


        • #5
          Originally posted by ragowthaman View Post
          Can you also check, if the reference fasta you created to align the bam file has same chromosome names as that of GFF.

          If you see only few reads are ignored the way you showed them here, its possible, your reference fasta has more chromosomes than in gff. For example, 'chromosome:AGPv2:1:1:301354135:1"
          Thanks! The reference fasta has the same number of chromosomes as in the GFF. It was the ID issue that I failed to see though it was staring at me.

          Comment

          Latest Articles

          Collapse

          • seqadmin
            The Impact of AI in Genomic Medicine
            by seqadmin



            Artificial intelligence (AI) has evolved from a futuristic vision to a mainstream technology, highlighted by the introduction of tools like OpenAI's ChatGPT and Google's Gemini. In recent years, AI has become increasingly integrated into the field of genomics. This integration has enabled new scientific discoveries while simultaneously raising important ethical questions1. Interviews with two researchers at the center of this intersection provide insightful perspectives into...
            02-26-2024, 02:07 PM
          • seqadmin
            Multiomics Techniques Advancing Disease Research
            by seqadmin


            New and advanced multiomics tools and technologies have opened new avenues of research and markedly enhanced various disciplines such as disease research and precision medicine1. The practice of merging diverse data from various ‘omes increasingly provides a more holistic understanding of biological systems. As Maddison Masaeli, Co-Founder and CEO at Deepcell, aptly noted, “You can't explain biology in its complex form with one modality.”

            A major leap in the field has
            ...
            02-08-2024, 06:33 AM

          ad_right_rmr

          Collapse

          News

          Collapse

          Topics Statistics Last Post
          Started by seqadmin, Today, 06:12 AM
          0 responses
          13 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 02-23-2024, 04:11 PM
          0 responses
          64 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 02-21-2024, 08:52 AM
          0 responses
          70 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 02-20-2024, 08:57 AM
          0 responses
          60 views
          0 likes
          Last Post seqadmin  
          Working...
          X