Seqanswers Leaderboard Ad

Collapse

Announcement

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

  • BWA - number of sequences aligned

    Hello everyone,

    I am quite a newbie in bioinformatics, and I have some difficulties using bwa software.

    I successfully aligned my dataset to my reference genome. Now, I want to see the impact of the settings on my results, but I don't know how to obtain the number of reads aligned in the .sam file.

    Is there a useful trick to obtain the number of reads aligned ?

    I have found no clue on the bwa manual, neither in the samtool package and I begin to feel quite helpless ...

    Thank you for your answers !

  • #2
    Originally posted by Cbon View Post
    Hello everyone,

    Is there a useful trick to obtain the number of reads aligned ?
    Try

    Code:
    samtools flagstat [I]file.bam[/I]
    Although I don't know if it works on sam files (in case you should convert in BAM, which I believe should be the default...).

    d

    Comment


    • #3
      If you are using UNIX and GNU awk,
      Code:
      gawk '!/^[@#]/ && (!and($2,4){print}' bwa.sam| wc -l
      should give you the number of reads aligned - providing BWA is setting the flag for unaligned reads. If you have many reads and BWA writes out only the aligned reads,
      Code:
      wc -l bwa.sam
      is a very good simple approximate of aligned reads (it counts the header of course, but the ~5 lines of the header compared to the millions of reads are negligible).

      Comment


      • #4
        Thank you both for your answers !

        My files are already converted in .bam, so I hope I can use the flagstat command. If it doesn't work, I will try the code.

        thank you again !

        Comment


        • #5
          Originally posted by racetrout View Post
          If you are using UNIX and GNU awk,
          Code:
          gawk '!/^[@#]/ && (!and($2,4){print}' bwa.sam| wc -l
          should give you the number of reads aligned - providing BWA is setting the flag for unaligned reads. If you have many reads and BWA writes out only the aligned reads,
          Code:
          wc -l bwa.sam
          is a very good simple approximate of aligned reads (it counts the header of course, but the ~5 lines of the header compared to the millions of reads are negligible).
          Yes, but removing the lines is easy too. Skipping 2 header lines (think I had that), and starting on the 3rd line, can be done with the tail command and adding a '+' before the number of lines, to denote which line to start from:

          Code:
          tail -n +3 bwa.sam | gawk '!/^[@#]/ && (!and($2,4){print}' | wc -l
          ... or for the case with only aligned reads coming from bwa:
          Code:
          tail -n +3 bwa.sam | wc -l

          Comment


          • #6
            bwa by default includes all reads, even unaligned ones. You can tell that they are unmapped, from the binary tag.

            So the flagstat command will distinguish between the mapped an unmapped reads.

            You can also use samtools view -bF0x4 all.bam > mapped .bam to get a bam which has all the unmapped reads filtered out.

            Comment


            • #7
              Use "samtools view -c -S <in.sam>" and add the "-f/-F" options as desired.

              Comment


              • #8
                I personally find it---at the best case---unacceptable the fact that a short-read aligner does not have an extra line of code to report the total number of aligned reads... Especially when it outputs (on the standard error) things like:

                `[infer_isize] (25, 50, 75) percentile: (137, 179, 227)'

                Anyhow, isn't there a closing parenthesis missing from this one?

                Originally posted by saml View Post
                Code:
                tail -n +3 bwa.sam | gawk '!/^[@#]/ && (!and($2,4){print}' | wc -l
                Thanks
                Last edited by solonas13; 05-13-2012, 03:47 PM.

                Comment

                Latest Articles

                Collapse

                • seqadmin
                  Strategies for Sequencing Challenging Samples
                  by seqadmin


                  Despite advancements in sequencing platforms and related sample preparation technologies, certain sample types continue to present significant challenges that can compromise sequencing results. Pedro Echave, Senior Manager of the Global Business Segment at Revvity, explained that the success of a sequencing experiment ultimately depends on the amount and integrity of the nucleic acid template (RNA or DNA) obtained from a sample. “The better the quality of the nucleic acid isolated...
                  03-22-2024, 06:39 AM
                • seqadmin
                  Techniques and Challenges in Conservation Genomics
                  by seqadmin



                  The field of conservation genomics centers on applying genomics technologies in support of conservation efforts and the preservation of biodiversity. This article features interviews with two researchers who showcase their innovative work and highlight the current state and future of conservation genomics.

                  Avian Conservation
                  Matthew DeSaix, a recent doctoral graduate from Kristen Ruegg’s lab at The University of Colorado, shared that most of his research...
                  03-08-2024, 10:41 AM

                ad_right_rmr

                Collapse

                News

                Collapse

                Topics Statistics Last Post
                Started by seqadmin, Yesterday, 06:37 PM
                0 responses
                11 views
                0 likes
                Last Post seqadmin  
                Started by seqadmin, Yesterday, 06:07 PM
                0 responses
                10 views
                0 likes
                Last Post seqadmin  
                Started by seqadmin, 03-22-2024, 10:03 AM
                0 responses
                51 views
                0 likes
                Last Post seqadmin  
                Started by seqadmin, 03-21-2024, 07:32 AM
                0 responses
                68 views
                0 likes
                Last Post seqadmin  
                Working...
                X