Unconfigured Ad

Collapse
X
 
  • Time
  • Show
Clear All
new posts
  • IonTom
    Member
    • Apr 2014
    • 32

    #1

    Python samtools BytesIO

    I currently have a a problem with samtools and python.
    Iḿ trying to write a sam or bam file into a file like object
    which works well using the following code

    import sh
    import io

    #define samtools as a python functions
    samtools = sh.Command("samtools")

    #where is the input sam
    sam_input = "test.sam"

    #generate a file like object
    sam_IO2 = io.BytesIO()

    #let samtools view write in the file like object
    samtools.view(sam_output,S=True,h=True,_out = sam_IO2)

    #show the sam file
    sam_IO2.getvalue()

    or

    bam_IO2 = io.BytesIO()
    samtools.view(bam_input,_out = bam_IO2)

    #show the bam file
    bam_IO2.getvalue()


    but now i try to read from the file like object for both cases:

    samtools.view("/dev/stdin",S=True,_in = sam_IO2)

    samtools.view("/dev/stdin",_in = bam_IO2)

    it fails in both tcases.


    Does anyone have an idea how to do this correctly ?
    Would be really helpful.
    Last edited by IonTom; 05-05-2014, 12:28 PM.
  • dpryan
    Devon Ryan
    • Jul 2011
    • 3478

    #2
    Without actually going through your code, is there a reason you're not just using pysam? That would seem simpler.

    Also, it's generally better to give the actual error message produced instead of just writing "it fails".

    Comment

    • IonTom
      Member
      • Apr 2014
      • 32

      #3
      Sorry forgot tom post this.
      Here is the error message:

      samtools.view("/dev/stdin",_in = bam_IO2)

      Traceback (most recent call last):
      File "<stdin>", line 1, in <module>
      File "/usr/local/lib/python2.7/dist-packages/sh.py", line 769, in __call__
      return RunningCommand(cmd, call_args, stdin, stdout, stderr)
      File "/usr/local/lib/python2.7/dist-packages/sh.py", line 330, in __init__
      self.wait()
      File "/usr/local/lib/python2.7/dist-packages/sh.py", line 334, in wait
      self._handle_exit_code(self.process.wait())
      File "/usr/local/lib/python2.7/dist-packages/sh.py", line 348, in _handle_exit_code
      self.process.stderr
      sh.ErrorReturnCode_1:

      RAN: '/usr/bin/samtools view /dev/stdin'

      STDOUT:

      STDERR:
      [bam_header_read] EOF marker is absent. The input is probably truncated.
      [bam_header_read] invalid BAM binary header (this is not a BAM file).
      [main_samview] fail to read the header from "/dev/stdin".

      Isn't pysam more or less doing the same a the sh solution ?
      Basically calling the c implementation ?
      So would it help to interact with the file like objects ?

      Comment

      • dpryan
        Devon Ryan
        • Jul 2011
        • 3478

        #4
        Pysam actually wraps the C API, so it's a bit more useful. Why bother reinventing the wheel.

        Comment

        • maubp
          Peter (Biopython etc)
          • Jul 2009
          • 1544

          #5
          I would also recommend trying pysam instead, but it could be your code is breaking due to a simple user error:

          Where is the BytesIO handle position within the data stream? Does a
          Code:
          .seek(0)
          to ensure it is at the start help? My guess is the handle position is at the end of the data, and thus when samtools tries to read from it there is no (more) data.

          Comment

          • IonTom
            Member
            • Apr 2014
            • 32

            #6
            Here is a solution to the problem above using pipes in python

            bam_input = "input.bam"
            bam_output = "output.bam"

            outputter = open(bam_output,mode = "w")

            p1 = subprocess.Popen(['samtools',"view","-u",bam_output],stdout=subprocess.PIPE,bufsize=1) #Set up the samtools view command to output to pipe
            p2 = subprocess.Popen(["samtools","view","-b","-"], stdin=p1.stdout,stdout=outputter,bufsize = 1) #send p1's output to p2 and p2 to bam file
            p1.stdout.close() #make sure we close the output so p2 doesn't hang waiting for more input
            output = p2.communicate()[0] #run our commands

            outputter.close()

            This example is just a toy example but the principle can also be used for real tasks.
            Last edited by IonTom; 05-13-2014, 03:19 PM.

            Comment

            Latest Articles

            Collapse

            • kinpict.com
              Reply to 16S/ITS Profiling vs. Shotgun Metagenomics: What Changes in Your Workflow?
              by kinpict.com
              The part about this topic feels surprisingly grounded, and the details on the page is the bit that made me slow down. I’d keep <a href="https://kinpict.com/">ai family photo generator</a> bookmarked for later.
              Today, 06:37 PM
            • SEQadmin2
              16S/ITS Profiling vs. Shotgun Metagenomics: What Changes in Your Workflow?
              by SEQadmin2





              NGS has revolutionized clinical microbiology by offering new ways to identify and surveil pathogens without prior knowledge of the target. One popular approach is 16S rRNA and ITS amplicon sequencing, which targets specific microbial markers. But studies are increasingly incorporating shotgun metagenomics, a technique that sequences all DNA from a sample without targeting specific microbial markers. This approach allows more detailed taxonomic classification,...
              10-30-2026, 06:46 AM

            ad_right_rmr

            Collapse

            News

            Collapse

            Topics Statistics Last Post
            Started by SEQadmin2, 11-12-2026, 07:59 AM
            0 responses
            26 views
            0 reactions
            Last Post SEQadmin2  
            Started by SEQadmin2, 11-05-2026, 09:25 AM
            0 responses
            35 views
            0 reactions
            Last Post SEQadmin2  
            Started by SEQadmin2, 11-03-2026, 11:54 AM
            2 responses
            56 views
            0 reactions
            Last Post kinpict.com  
            Started by SEQadmin2, 10-28-2026, 09:57 AM
            0 responses
            175 views
            0 reactions
            Last Post SEQadmin2  
            Working...