Problem with downloading CDS sequences using Biopython
I'm trying to download CDS sequences for a given genome using Biopython. My script looks like this:
from Bio import Entrez
from Bio import SeqIO
Entrez.email = "c...@gmail.com"
genomeAccessions = ['NC_021353.1', 'NC_020913.1']
handle = Entrez.efetch(db="nucleotide", id=genomeAccessions, rettype="gb")
records = SeqIO.parse(handle, "gb")
for i,record in enumerate(records):
print(len(record.features))
for feature in record.features:
if feature.type == "CDS":
print feature.location
print feature.qualifiers["protein_id"]
print feature.location.extract(record).seq
But using this code I get only one feature (for example for genome NC_021353) even though there are many features http://www.ncbi.nlm.nih.gov/nuccore/NC_021353.
I would be grateful for any suggestion what I'm doing wrong.
1 answer
Your rettype must be "gbwithparts", this should work:
handle = Entrez.efetch(db="nucleotide", id=genomeAccessions, rettype="gbwithparts")
Running your program:
>>>
======================== RESTART: D:/Desktop/test.py ========================
3823
[580:1822](+)
['WP_020447944.1']
ATGGACAACATCAGCGGCAACATCTTTG...
...
[1818:2619](-)
['WP_020447945.1']
GTGAATAAAGAGAACAAATCTACTTTAA...
...
Your code will break later, possibly there are entries without a protein_id:
Traceback (most recent call last):
File "D:/Desktop/test.py", line 20, in <module>
print feature.qualifiers["protein_id"]
KeyError: 'protein_id'
Log in to answer this question.
More posts like this
-
Getting cds files from ncbi using accession numbers
written by lily •
-
Empty in genebank extractor python code
written by MSRS
-
Failure to extract CDS from embl file using BioPython
written by rjqmantaring •
-
parsing gbk file to get fna file
written by ulises.rodriguez •
-
Intron extraction with Biopython.
written by Evergreen •
-
Compare 2 local protein FASTA using NCBI BLAST Biopython
written by ahtmatrix •
-
Downloading Gene Annotation protein.faa files from Genome Accession numbers from NCBI via Entrez
written by ijc2 •
-
Parse the Xml response from Entrez Db=bioproject Using Biopython eFetch
written by Prasad •
-
Problem with writing edited feature qualifier to a new gbk file in Biopython
written by A •
-
How to get .faa files from a lot of genome from ncbi
written by CrLs •
I think something similar has been already discussed here:
But this solution is not in Biopython.
Download cds region coordinates
This one is.
How To Extract Just 'Cds' From Genbank File Into Another Genbank File?
There are some other suggestions here:
BioPython error parsing standard GenBank file