Seqanswers Leaderboard Ad

Collapse

Announcement

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

  • Adding features to a multi entry genbank

    I'm attempting to add features using BioPython to a multi-entry Genbank. Each entry in the genbank is a contig, and I'd like to annotate some BLAST hits.

    I've performed a BLAST of my query against a multi-entry fasta, so I know which hits belong to which contig and what the contig is called. I'm able to parse my Genbank and find which contig I'd like to add the feature to. I then append the feature, but I'm problems at the final stage where I attempt to then write out the entire Genbank with the new feature added to one of the entries.

    This is my code so far:


    Code:
    >>> handle = open("D13_multi.gbk", "rU")
    >>> genbank = SeqIO.parse(handle, "genbank")
    >>> for record in genbank:
    ...     if record.id == node:
                    record.features.append(SeqFeature.SeqFeature(SeqFeature.FeatureLocation(start, end), type = "misc_feature"))
    ...             print record.features
    ...
    [SeqFeature(FeatureLocation(ExactPosition(0), ExactPosition(35), strand=1), type='misc_feature'), SeqFeature(FeatureLocation(ExactPosition(652), ExactPosition(145)), type='misc_feature')]
    After this for loop, SeqIO.write() writes 0 features to the new file I specify.
    If I do this - SeqIO.write(genbank, "D13_multi_edited.gbk", "genbank"), it writes out all the entries after the the contig that I've just append a feature to.

    I've also tried this:


    Code:
    >>> for record in genbank:
    ...     if record.id != node:
    ...             SeqIO.write(genbank, "D13_test.gbk", "genbank")
    ...     if record.id == node:
    ...             record.features.append(SeqFeature.SeqFeature(SeqFeature.FeatureLocation(start, end), type = "misc_feature"))
    ...             SeqIO.write(record, "D13_test.gbk", "genbank")
    ...
    109
    In this case it misses the first entry (there should be 110 entries) and writes out the entry that I've added a feature to but without the feature being added.

    What am I doing wrong?
    Last edited by yekwah; 11-11-2013, 01:37 PM. Reason: Attempted to fix indenting for code

  • #2
    When editing a post choose "Edit" --> "Go Advanced" button at the bottom of the edit window. In the subsequent window that opens up use the mose to highlight the text you to designate as "code" then use the "#" button in the edit menu to add
    Code:
    $ code_example
    to properly format your code.

    For example:
    Code:
    >>> handle = open("D13_multi.gbk", "rU")
    >>> genbank = SeqIO.parse(handle, "genbank")
    >>> for record in genbank:
    ... if record.id == node:
    record.features.append(SeqFeature.SeqFeature(SeqFeature.FeatureLocation(start, end), type = "misc_feature"))
    ... print record.features
    ...
    [SeqFeature(FeatureLocation(ExactPosition(0), ExactPosition(35), strand=1), type='misc_feature'), SeqFeature(FeatureLocation(ExactPosition(652), ExactPosition(145)), type='misc_feature')]
    Try it on the post above.
    Last edited by GenoMax; 11-11-2013, 04:12 AM.

    Comment


    • #3
      IF you can add the [ code ] ... [ /code ] tags it would help, one possibility is you keep calling SeqIO.write and each time over-write it with a singe record?
      Last edited by maubp; 11-11-2013, 10:31 AM.

      Comment


      • #4
        Thanks. I've fixed the indentation for the code now, so hopefully that will help.

        I'm pretty certain that if I keep calling SeqIO.write that I overwrite whatever I had already. It seems like when I open the multi entry genbank with SeqIO.parse, it gets to the end of the file and then doesn't go back to the start. So in the for loop, I loop over all the records, get to the record I want and then add the feature, then it loops over the rest of the records. Once the for loop is done, I call SeqIO.write, but there isn't anything to write because its gotten to the end of the file.

        So I'm just not sure how to write out my edited genbank with the features I want added in.

        Comment


        • #5
          OK, keeping this simple (but not memory efficient), load all the records into memory as a list, modify the list, save list of records:

          Code:
          from Bio import SeqIO
          from Bio.SeqFeature import SeqFeature, FeatureLocation
          node = ...
          records = list(SeqIO.parse("D13_multi.gbk", "gb"))
          for record in records:
              if record.id == node:
                  record.features.append(SeqFeature(FeatureLocation(start, end), type = "misc_feature"))
          SeqIO.write(records, "D13_edited.gbk", "gb")
          Getting a bit more clever, avoid loading all the records into memory, loop over the file once - and writing out the modified file one record at a time in the main loop:

          Code:
          from Bio import SeqIO
          from Bio.SeqFeature import SeqFeature, FeatureLocation
          node = ...
          handle = open("D13_edited.gbk", "w")
          for record in SeqIO.parse("D13_multi.gbk", "gb"):
              if record.id == node:
                  record.features.append(SeqFeature(FeatureLocation(start, end), type = "misc_feature"))
              SeqIO.write(record, handle, "gb")
          handle.close()
          If you want to do this in a memory efficient way (like the second version) but also making one call to SeqIO.write, you could define a generator function, or use a generator expression with a custom function to annotate an individual SeqRecord. But I'm guessing that might be a little too complicated for now?

          Comment


          • #6
            Thank you so much for this!

            The second option seems to be the best, though I would like to be as memory efficient as possible. I'm not entirely sure what a generator function is, so that might be a little complicated for the moment.

            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, 02-28-2024, 06:12 AM
            0 responses
            27 views
            0 likes
            Last Post seqadmin  
            Started by seqadmin, 02-23-2024, 04:11 PM
            0 responses
            74 views
            0 likes
            Last Post seqadmin  
            Started by seqadmin, 02-21-2024, 08:52 AM
            0 responses
            81 views
            0 likes
            Last Post seqadmin  
            Started by seqadmin, 02-20-2024, 08:57 AM
            0 responses
            69 views
            0 likes
            Last Post seqadmin  
            Working...
            X