I modified the GFF3 file with the below code:
#!/usr/bin/python3
import click
import gffutils
import gffutils.gffwriter as gffwriter
@click.command()
@click.option('--gff3', help="Provide GFF3 file", required=True)
@click.option('--keep', help="Keep GFF3 file", required=True)
@click.option('--reject', help="Reject GFF3 file", required=True)
def run(gff3, keep, reject):
#db = gffutils.create_db(gff3, dbfn='data/test.db', force=True, keep_order=True, merge_strategy='merge',
# sort_attribute_values=True)
out_keep = gffwriter.GFFWriter(keep)
db = gffutils.FeatureDB('data/test.db', keep_order=True)
for gene in db.features_of_type('gene', order_by='start'):
print(gene)
for i in db.children(gene, featuretype='mRNA', order_by='start'):
print(i)
i.attributes['Note'] = ["Gene description"]
print(i)
out_keep.write_rec(gene)
Unfortunately, the above code saved only this:
#GFF3 file (created by gffutils (v0.9) on 2019-08-28 22:14:13)
NbV1Ch08 AUGUSTUS gene 7015 29794 0.01 - . ID=g1;
NbV1Ch08 AUGUSTUS gene 60876 63944 0.03 + . ID=g2;
NbV1Ch08 AUGUSTUS gene 64722 65524 0.32 - . ID=g3;
NbV1Ch08 AUGUSTUS gene 66706 67484 0.12 - . ID=g4;
NbV1Ch08 AUGUSTUS gene 115451 122754 0.16 - . ID=g5;
NbV1Ch08 AUGUSTUS gene 147786 155025 0.24 + . ID=g6;
NbV1Ch08 AUGUSTUS gene 155152 161104 0.28 - . ID=g7;
NbV1Ch08 AUGUSTUS gene 237126 238098 0.05 + . ID=g8;
NbV1Ch08 AUGUSTUS gene 265796 297830 0.82 + . ID=g9;
NbV1Ch08 AUGUSTUS gene 265796 297830 0.82 + . ID=g9;
Rather than:
NbV1Ch08 AUGUSTUS gene 7015 29794 0.01 - . ID=g1;
NbV1Ch08 AUGUSTUS mRNA 7015 29794 0.01 - . ID=g1.t1;Parent=g1;Note=Gene description;
NbV1Ch08 AUGUSTUS transcription_end_site 7015 7015 . - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS three_prime_utr 7015 8531 0.2 - . ID=g1.t1.3UTR1;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 7015 8747 . - . ID=g1.t1.exon1;Parent=g1.t1;
NbV1Ch08 AUGUSTUS stop_codon 8532 8534 . - 0 Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 8532 8747 0.31 - 0 ID=g1.t1.CDS1;Parent=g1.t1
NbV1Ch08 AUGUSTUS intron 8748 9191 0.49 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 9192 9342 0.66 - 1 ID=g1.t1.CDS2;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 9192 9342 . - . ID=g1.t1.exon2;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 9343 9915 0.58 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 9916 10006 0.71 - 2 ID=g1.t1.CDS3;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 9916 10006 . - . ID=g1.t1.exon3;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 10007 10101 0.74 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 10102 10201 0.78 - 0 ID=g1.t1.CDS4;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 10102 10201 . - . ID=g1.t1.exon4;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 10202 10712 0.8 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 10713 11107 0.11 - 2 ID=g1.t1.CDS5;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 10713 11107 . - . ID=g1.t1.exon5;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 11108 11569 0.07 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 11570 12151 0.09 - 2 ID=g1.t1.CDS6;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 11570 12151 . - . ID=g1.t1.exon6;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 12152 12588 0.34 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 12589 12717 0.39 - 2 ID=g1.t1.CDS7;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 12589 12717 . - . ID=g1.t1.exon7;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 12718 12789 0.42 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 12790 13075 0.39 - 0 ID=g1.t1.CDS8;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 12790 13075 . - . ID=g1.t1.exon8;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 13076 14832 0.51 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 14833 15009 0.39 - 0 ID=g1.t1.CDS9;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 14833 15009 . - . ID=g1.t1.exon9;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 15010 15278 0.59 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 15279 15415 0.56 - 2 ID=g1.t1.CDS10;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 15279 15415 . - . ID=g1.t1.exon10;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 15416 15487 0.58 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 15488 15612 0.96 - 1 ID=g1.t1.CDS11;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 15488 15612 . - . ID=g1.t1.exon11;Parent=g1.t1;
NbV1Ch08 AUGUSTUS intron 15613 15706 0.96 - . Parent=g1.t1;
NbV1Ch08 AUGUSTUS CDS 15707 15957 0.98 - 0 ID=g1.t1.CDS12;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 15707 15958 . - . ID=g1.t1.exon12;Parent=g1.t1;
NbV1Ch08 AUGUSTUS start_codon 15955 15957 . - 0 Parent=g1.t1;
NbV1Ch08 AUGUSTUS five_prime_utr 15958 15958 0.99 - . ID=g1.t1.5UTR1;Parent=g1.t1
NbV1Ch08 AUGUSTUS five_prime_utr 27458 28250 0.37 - . ID=g1.t1.5UTR2;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 27458 28250 . - . ID=g1.t1.exon13;Parent=g1.t1;
NbV1Ch08 AUGUSTUS five_prime_utr 29272 29794 0.08 - . ID=g1.t1.5UTR3;Parent=g1.t1
NbV1Ch08 AUGUSTUS exon 29272 29794 . - . ID=g1.t1.exon14;Parent=g1.t1;
NbV1Ch08 AUGUSTUS transcription_start_site 29794 29794 . - . Parent=g1.t1;
What did I miss?
Thank you in advance,
Update This is base on @daniel-standage code.
import click
import tag
def add_description(instream):
for feature in instream:
for mrna in tag.select.features(feature, type='mRNA', traverse=True):
mrna.add_attribute('Note', 'Gene description')
yield feature
@click.command()
@click.option('--gff3', help="Provide GFF3 file", required=True)
@click.option('--keep', help="Keep GFF3 file", required=True)
@click.option('--reject', help="Reject GFF3 file", required=True)
def run(gff3, keep, reject):
reader = tag.GFF3Reader(infilename=gff3)
desc_stream = add_description(reader)
writer = tag.GFF3Writer(desc_stream)
writer.write()
if __name__ == '__main__':
run()
which caused this error:
mod_braker_gff3.py --gff3 data/NbV1Ch08-augustus.hints_utr.gff3 --keep data/keep.gff3 --reject data/reject.gff3
Traceback (most recent call last):
File "/Users/lorencm/projects/bioinf-scripts/mod_braker_gff3.py", line 26, in <module>
run()
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/click/core.py", line 764, in __call__
return self.main(*args, **kwargs)
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/click/core.py", line 717, in main
rv = self.invoke(ctx)
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/click/core.py", line 956, in invoke
return ctx.invoke(self.callback, **ctx.params)
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/click/core.py", line 555, in invoke
return callback(*args, **kwargs)
File "/Users/lorencm/projects/bioinf-scripts/mod_braker_gff3.py", line 22, in run
writer.write()
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/tag/writer.py", line 57, in write
for entry in self._instream:
File "/Users/lorencm/projects/bioinf-scripts/mod_braker_gff3.py", line 8, in add_description
for mrna in tag.select.features(feature, type='mRNA', traverse=True):
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/tag/select.py", line 24, in features
for feature in entry_type_filter(entrystream, tag.Feature):
File "/anaconda/envs/bioinf-scripts/lib/python3.7/site-packages/tag/select.py", line 99, in entry_type_filter
for entry in entrystream:
TypeError: 'Directive' object is not iterable
Process finished with exit code 1
print db.children(gene, featuretype='mRNA', order_by='start')
, you'd presumably see that it is empty. $\endgroup$