Seqanswers Leaderboard Ad

Collapse

Announcement

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

  • gsnap produces more reads than fastq

    hello, i am running gsnap with rna seq data against a reference genome:

    gsnap -m 10 -B 5 -t 8 -A sam -d GENOME seqs.fastq > res.sam

    but strangely, the output sam file contains more reads than the original fastq files:

    $ cat seqs.fastq | echo $((`wc -l`/4))
    3776979

    $ cat res.sam | grep -v '^ *@' | wc -l
    6009141



    and I don't understand what that means....

    thank you

  • #2
    I'm not familiar with gsnap but it seems like there could be a few factors, the most likely one being your grep command. When I try to find out how many alignments have been reported I use the command:
    Code:
    cut -f1 out.sam | sort | uniq | wc -l
    That way you don't count nonspecific alignments. That being said, your SAM file may contain multiple alignments for many reads. There should be a flag to control how many alignments per read are reported but if not there are ways around it. I hop this helps!

    Comment


    • #3
      Originally posted by twaddlac View Post
      I'm not familiar with gsnap but it seems like there could be a few factors, the most likely one being your grep command. When I try to find out how many alignments have been reported I use the command:
      Code:
      cut -f1 out.sam | sort | uniq | wc -l
      That way you don't count nonspecific alignments. That being said, your SAM file may contain multiple alignments for many reads. There should be a flag to control how many alignments per read are reported but if not there are ways around it. I hop this helps!
      Besides mapping to multiple locations, the sam file also has a header that can contain many lines and throw off your count.

      Comment


      • #4
        Originally posted by severin View Post
        Besides mapping to multiple locations, the sam file also has a header that can contain many lines and throw off your count.
        I forgot to mention, to omit the header lines you should do
        Code:
        grep -v '^@' out.sam | cut -f1 | sort | uniq | wc -l
        If you convert your SAM file to BAM you will:
        A) reduce the size of the file (BAM < SAM)
        B) view the alignment output without the headers
        Code:
        samtools view out.bam
        To convert the sam file to bam, just do
        Code:
        samtools view -bS out.sam > out.bam

        Comment


        • #5
          If a read is aligned to multiple location, it appears multiple times in the SAm file. Sort the file by read name, cut to the first coulmn (read name), pipe it through 'uniq', then count again.

          Comment

          Latest Articles

          Collapse

          • seqadmin
            Current Approaches to Protein Sequencing
            by seqadmin


            Proteins are often described as the workhorses of the cell, and identifying their sequences is key to understanding their role in biological processes and disease. Currently, the most common technique used to determine protein sequences is mass spectrometry. While still a valuable tool, mass spectrometry faces several limitations and requires a highly experienced scientist familiar with the equipment to operate it. Additionally, other proteomic methods, like affinity assays, are constrained...
            04-04-2024, 04:25 PM
          • 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

          ad_right_rmr

          Collapse

          News

          Collapse

          Topics Statistics Last Post
          Started by seqadmin, 04-11-2024, 12:08 PM
          0 responses
          18 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 04-10-2024, 10:19 PM
          0 responses
          22 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 04-10-2024, 09:21 AM
          0 responses
          17 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 04-04-2024, 09:00 AM
          0 responses
          49 views
          0 likes
          Last Post seqadmin  
          Working...
          X