Is there any way of using biopython to write Swissprot files?

Question

I have the following mwe for filtering a Swissprot file based on a certain feature, in this case, transmembrane regions.

from Bio import SeqIO
records=[]
for record in SeqIO.parse("Input.txt", "swiss"):
    transmembrane_protein=False
    print record.id
    for i, feature in enumerate(record.features):
        if feature.type == "TRANSMEM":
            transmembrane_protein=True
    if transmembrane_protein==True:
        records.append(record)
SeqIO.write(records, "Output.txt", "swiss")

The script works when SeqIO.write(records, "Output.txt", "swiss") becomes SeqIO.write(records, "Output.txt", "fasta")

However, such a method is not yet supported.

ValueError: Reading format 'swiss' is supported, but not writing

From the docs, I see that writing is not supported for swiss:

Note that while Bio.SeqIO can read all the above file formats, it cannot write to all of them.

Is there any unofficial way of using biopython/python to write swissprot files from parsed swissprot files?

Answer 1

Using SeqIO.index rather than SeqIO.parse lets you read all the records into a dict, from which you can then extract the raw entry:

from Bio import SeqIO

record_dict = SeqIO.index('Input.txt', 'swiss')
for key in record_dict:
    print(record_dict.get_raw(key).decode())

Now you should be able to apply your test for a transmembrane protein to each entry, and write out only the ones you want to keep.

from Bio import SeqIO

record_dict = SeqIO.index('Input.txt', 'swiss')
out = open('Output.txt', 'wb')
for key in record_dict:
    record = record_dict[key]
    for i, feature in enumerate(record.features):
    if feature.type == "TRANSMEM":
        out.write(record_dict.get_raw(key).decode())

Answer 2

I don't know how to do this with biopython, but it's simple enough in Perl using Swisskinfe:

#!/usr/bin/env perl
use strict;
use warnings;
use SWISS::Entry;

my $flat=$ARGV[0]||die("Need a flat file\n");

## Change the line termination string so we read an 
## entire entry at a time
local $/ = "\n//\n";

my $fh;
## Deal with compressed or uncompressed files
if($flat=~/\.gz$/){
    open($fh,"zcat $flat |")|| die("cannot open flat file $flat : $! \n");
}
else{
    open($fh,"$flat")|| die("cannot open flat file $flat : $!: $@\n");
}
## Parse the file. 
entry: while (my $current_entry = <$fh>) {
    ## Read the current entry as an "entry" object
    my $entry = SWISS::Entry->fromText($current_entry);
    ## Iterate over the entry's features
    foreach my $feature ($entry->FTs->elements()) {
        ## If any of the features are 'TRANSMEM', print the entry
        if ($feature->[0] eq 'TRANSMEM') {
            print "$current_entry\n";
            ## Skip to the next entry
            next entry;
        }
    }
}

You need to install the Swissknife module first:

1. Download Swissknife_1.73.tar.gz from here (SourceForge link).

2. Install

   tar xvzf Swissknife.tar.gz
   cd Swissknife_1.73/
   perl Makefile.PL
   sudo make install
   

Now, make the script above executable and run with:

foo.pl file.flat