diekhans
commented
6 years ago From: Mark Diekhans markd@soe.ucsc.edu To: NCBI Genomes genomes@ncbi.nlm.nih.gov Cc: Ian Fiddes ian.t.fiddes@gmail.com Subject: RE: GFF input for ape genome (NBAG00000000) Date: Tue, 3 Oct 2017 17:15:57 -0700
Hi Karen,
The file you looked at did have a bug we have just fixed were our code to insert gaps to avoid frameshifts end up creating downstream frameshifts Apologizes.
We shifted from GFF3 to tbl format because table2asn didn't honor the phase and didn't output GBFF files. It just easier to control things with the tbl format. GBFF is much easier for understand the results than the ASN.1
Answers (and questions) to your questions below.
NCBI Genomes genomes@ncbi.nlm.nih.gov writes:
[1] What's the ncRNA gnl|CK280|T0170363 doing there?
GENCODE annotates non-coding transcripts in coding genes. These are consider part of the same gene, even though they don't code for proteins.
[2] You need to include the protein_id on the corresponding mRNA in the .tbl file OR use GFF and rely on parent:child for linkages (will use the locus_tag as the base). If this came from the .sqn that came from the .gff, then apologies. We are having problems with that GFF converter & hoped to have it finally working next week.
We will add this.
[3] don't include ncRNA processed transcripts
While I have mixed feelings about the `processed transcripts' category, they are part of the gene set. One successful mapped to another species has value.
[4] remove the GENCODE_biotype|xxx notes, eg GENCODE_biotype| processed_transcript. Why are they included?
Why do they need to be removed? Aren't they comments? The biotype of the source transcript is very useful in understanding the annotation, especially since the mapping from GENCODE biotype to feature is not one-to-one. I will change the note to make it a bit clearer: CAT source GENCODE transcript biotype: XXX
[5] don't include ncRNA_class qualifiers on mRNAs. They belong just on ncRNAs. ANd they cause WrongQualOnFeature errors.
That is a mistake, we will fix it.
[6] why does one mRNA have "GENCODE_biotype|nonsense_mediated_decay"? Is CDS gnl|CK280|T0170365_prot not actually made? If not, then don’t include that CDS/ mRNA.
There is good evidence that that NMD is an imperfect mechanism (e.g. doi:10.1038/sj.ejhg.5201649) and thus GENCODE annotates the CDS. It's also really useful for understanding why the transcript is considered NMD if you know where the CDS and not have recompute it.
[7] why do they all have the same product name but have different translations? Include 'isoform' or something on the product name, yes?
We just used the gene name due to the lack of a standardized way to name protein isoforms (stable isofrom names are something GENCODE want's to implement).
We can tack on an `_isoformN' suffix, but this will only be for the NCBI submission, not a stable identifier.
[8] in the .tbl format, you have to provide the codon_start, meaning the first base of the first full codon of the CDS. Eg this doesn't have any internal stop codons if it begins with codon_start=2 (so bp519295 is the first base of the first full codon, since this is on the minus strand):
We will add this.
[8] some product names are actually the _ids:
This is because these genes have no HUGO symbol. I will generate better names that don't include the `gnl|'.
[9] in the .tbl format, we need the partial symbols added. (The conversion program will add them when GFF is the input)
We can get this in, at least for the CDS.
[10] internal stop codons are illegal, so you need to do one of these things for those CDS (which is in my earlier email below your reply, below):
our bug fix greatly reduces the number of in-frame stops.
- if you think that the gene really does encode a protein but the nucleotide sequence is poor, then change the CDS location so that it produces the desired translation AND include the "low-quality sequence region" exception on the CDS. That will allow the CDS to have a translation & will quiet the validation errors. The protein definition line will be prepended with "LOW QUALITY:" as a warning to database users.
This is the approach we are taking.
Thank you!
Mark
ifiddes
Here is what NCBI complained about:
I looked over your .tbl file and used it to annotate 200 of the sequences in your submission (the first 100, and then 301-400th). I had to split the fasta file like that because tbl2asn was choking on the fasta file that was >2Gb. Here’s the command line for the first 100 (I didn’t include the discrepancy report because the errors were enough of an issue to deal with):
tbl2asn -i clint_ctg_200bp_exclude0TooBig00000001.fsa -f Clint_Chimp.tbl -M n -j "[organism=Pan troglodytes]" -t SubmissionTemplate.sbt -o clint_ctg_200bp_exclude0TooBig00000001.annot.1003.sqn &
fyi, it spews a warning for each identifier in the .tbl file that isn’t in the .fsa file. But that’s expected with this set-up since there’s only 100 sequences in the fasta file.
Please let us know of any questions.
Thanks, Karen
Here are the errors from those first 100 sequences:
================================================================= 28 REJECT-level messages exist
SEQ_FEAT.Range 28
================================================================= 61568 ERROR-level messages exist
SEQ_FEAT.MissingGeneXref 15 SEQ_FEAT.InternalStop 15041 SEQ_FEAT.MissingCDSproduct 64 SEQ_FEAT.NoStop 3386 SEQ_INST.StopInProtein 15086 SEQ_FEAT.UnindexedFeature 16 SEQ_FEAT.NoProtein 19 SEQ_FEAT.WrongQualOnFeature 2240 = ncRNA_class on mRNAs SEQ_FEAT.CDSwithMultipleMRNAs 231 SEQ_FEAT.AbuttingIntervals 114 SEQ_FEAT.StartCodon 12299 SEQ_FEAT.MissingQualOnFeature 391 = missing ncRNA_class on ncRNA SEQ_FEAT.ShortIntron 77 SEQ_FEAT.BadInternalCharacter 76 = gnl|CK280|T..... as the product name SEQ_FEAT.PartialProblem 20 SEQ_FEAT.FeatureProductInconsistency 30 SEQ_FEAT.TransLen 25 = extends past Stop codon (maybe an internalstop codon?) SEQ_FEAT.PseudoCdsViaGeneHasProduct 101 SEQ_INST.BadProteinStart 12336 SEQ_INST.ShortSeq 1
================================================================= 485337 WARNING-level messages exist
SEQ_FEAT.FeatContentDup 771 SEQ_FEAT.ShortExon 213 SEQ_FEAT.NotSpliceConsensusDonor 222540 SEQ_FEAT.InconsistentPseudogeneValue 177 SEQ_FEAT.MultipleGeneOverlap 141 SEQ_FEAT.PartialProblem 8842 SEQ_FEAT.NotSpliceConsensusAcceptor 231978 SEQ_FEAT.ProductLength 3 SEQ_FEAT.CDSmRNAmismatch 17201 SEQ_FEAT.CDSmRNArange 1286 SEQ_FEAT.CDSwithNoMRNA 271 SEQ_FEAT.DuplicateFeat 1914
================================================================= 14702 NOTE-level messages exist
SEQ_FEAT.RareSpliceConsensusDonor 14682 SEQ_FEAT.PseudoCDSmRNArange 20
=================================================================
These are the kinds of errors that are definitely problems:
SEQ_FEAT.Range = the location of the feature extends beyond the end of the sequence SEQ_FEAT.InternalStop/SEQ_INST.StopInProtein SEQ_FEAT.TransLen = extends past Stop codon (maybe an internalstop codon?) SEQ_FEAT.NoStop = needs a stop codon OR to be partial SEQ_FEAT.StartCodon/SEQ_INST.BadProteinStart = needs a start codon OR to be partial SEQ_FEAT.AbuttingIntervals = needs the “low-quality sequence region” exception SEQ_FEAT.ShortIntron = needs the “low-quality sequence region” exception SEQ_INST.ShortSeq = protein <10 amino acids. The preference is for proteins to be at least 30aa.
SEQ_FEAT.BadInternalCharacter = gnl|CK280|T..... as the product name SEQ_FEAT.WrongQualOnFeature = ncRNA_class on mRNAs SEQ_FEAT.MissingQualOnFeature = missing ncRNA_class on ncRNA SEQ_FEAT.PseudoCdsViaGeneHasProduct = /pseudo not set up correctly
These may get sorted out: SEQ_FEAT.MissingCDSproduct SEQ_FEAT.NoProtein SEQ_FEAT.CDSwithMultipleMRNAs SEQ_FEAT.FeatureProductInconsistency SEQ_FEAT.CDSmRNArange
And here is a closer look at the actual annotation:
In 000400F_1_903666_quiver_pilon, Focus on CK280_G0054588:
.tbl file:
519401 454972 gene locus_tag CK280_G0054588 518622 518370 ncRNA 518271 518095 457283 454972 transcript_id gnl|CK280|T0170363 ncRNA_class other product MKRN3 note GENCODE_biotype|processed_transcript 519362 518992 mRNA 474410 473483 transcript_id gnl|CK280|T0170364 product MKRN3 note GENCODE_biotype|protein_coding 519296 518992 CDS 474410 474398 product MKRN3 protein_id gnl|CK280|T0170364_prot transcript_id gnl|CK280|T0170364 519401 519246 mRNA 518271 518095 509719 509583 transcript_id gnl|CK280|T0170365 ncRNA_class other product MKRN3 note GENCODE_biotype|nonsense_mediated_decay 519296 519246 CDS 518271 518236 product MKRN3 protein_id gnl|CK280|T0170365_prot transcript_id gnl|CK280|T0170365 519306 518992 mRNA 518271 518095 509719 509583 transcript_id gnl|CK280|T0170366 product MKRN3 note GENCODE_biotype|protein_coding 519296 518992 CDS 518271 518095 509719 509653 product MKRN3 protein_id gnl|CK280|T0170366_prot transcript_id gnl|CK280|T0170366 519395 517058 mRNA transcript_id gnl|CK280|T0170368 product MKRN3 note GENCODE_biotype|protein_coding 519296 517773 CDS product MKRN3 protein_id gnl|CK280|T0170368_prot transcript_id gnl|CK280|T0170368
flatfile view from the .sqn file:
Questions & comments:
[1] What's the ncRNA gnl|CK280|T0170363 doing there? [2] You need to include the protein_id on the corresponding mRNA in the .tbl file OR use GFF and rely on parent:child for linkages (will use the locus_tag as the base). If this came from the .sqn that came from the .gff, then apologies. We are having problems with that GFF converter & hoped to have it finally working next week. [3] don't include ncRNA processed transcripts [4] remove the GENCODE_biotype|xxx notes, eg GENCODE_biotype|processed_transcript. Why are they included? [5] don't include ncRNA_class qualifiers on mRNAs. They belong just on ncRNAs. ANd they cause WrongQualOnFeature errors. [6] why does one mRNA have "GENCODE_biotype|nonsense_mediated_decay"? Is CDS gnl|CK280|T0170365_prot not actually made? If not, then don’t include that CDS/mRNA. [7] why do they all have the same product name but have different translations? Include 'isoform' or something on the product name, yes? [8] in the .tbl format, you have to provide the codon_start, meaning the first base of the first full codon of the CDS. Eg this doesn't have any internal stop codons if it begins with codon_start=2 (so bp519295 is the first base of the first full codon, since this is on the minus strand):
->
Although that might not be the problem; it could be indels in the sequence, which would have a different solution (below).
[8] some product names are actually the _ids:
518713 519222 gene locus_tag CK280_G0054590 518713 518741 mRNA 519156 519222 transcript_id gnl|CK280|T0170369 product gnl|CK280|T0170369 note GENCODE_biotype|protein_coding <518713 518741 CDS 519156 >519222 product gnl|CK280|T0170369 protein_id gnl|CK280|T0170369_prot transcript_id gnl|CK280|T0170369
[9] in the .tbl format, we need the partial symbols added. (The conversion program will add them when GFF is the input)
So the .tbl file should be like this, with the desired product name, 5' and 3' partial symbols on the features since the CDS doesn't have a start or stop codon, protein_id & transcript_id on both the mRNA & CDS:
<518713 >519222 gene locus_tag CK280_G0054590 <518713 518741 mRNA 519156 >519222 transcript_id gnl|CK280|T0170369 protein_id gnl|CK280|T0170369_prot product hypothetical protein <518713 518741 CDS 519156 >519222 product hypothetical protein protein_id gnl|CK280|T0170369_prot transcript_id gnl|CK280|T0170369 [10] internal stop codons are illegal, so you need to do one of these things for those CDS (which is in my earlier email below your reply, below): * add /pseudo to the gene if you think that the CDS is 'broken', meaning that it doesn't translate as you expect. (in a .tbl file, you'll change 'product' to 'note'; in a .gff file you don't need to make that change) * add /pseudogene (and the right qualifier) if you think the gene is actually a pseudogene (in a .tbl file, you'll change 'product' to 'note'; in a .gff file you don't need to make that change) * replace the CDS and gene with a misc_feature if you think that this is just a remnant of similarity that was found, but isn't a real gene * if you think that the gene really does encode a protein but the nucleotide sequence is poor, then change the CDS location so that it produces the desired translation AND include the "low-quality sequence region" exception on the CDS. That will allow the CDS to have a translation & will quiet the validation errors. The protein definition line will be prepended with "LOW QUALITY:" as a warning to database users. B. Sequence 000100F_1_9097766_quiver_pilon comments [11] Not sure what this is supposed to do. However, if bp31124 is an insert relative to the human genome, then you'd adjust the CDS location to jump over it & add the "low-quality sequence region" exception. CDS join(30043..30064,31121..31348,32921..32954,36740..36817, 40104..40284,40924..41078,41244..41384,42592..42685, 44775..44855,50299..50472,51704..51793,51950..52118, 54688..54740) /locus_tag="CK280_G0032352" /codon_start=1 /product="GSPT1" /protein_id="CK280:T0101284_prot" /translation="- TFRTYW*KWRDRNVSRRIMGAQRRNK*SRARGWFLGRWKAARGKCP*NDGGGRGNPKT *VCGCTARCS*ERACKCSIHWARRCWQVNHWRTNNVYLTGMVDKRTLEKYEREAKEKN RETWYLSWALDTNQEERDKGKTVEVGRAYFETEKKHFTILDAPGHKSFVPNMIGGASQ ADLAVLVISARKGEFETGFEKGGQTREHAMLAKTAGVKHLIVLINKMDDPTVNWSNER YEECKEKLVPFLKKVGFNPKKDIHFMPCSGLTGANLKEQSDFCPWYIGLPFIPYLDNL PNFNRSVDGPIRLPIVDKYKIWALWSWESWNQDLFVKASSL**CQTRNVEVLGILSDD VETDTVAPGENLKIRLKGIEEEEILPGFILCDPNNLCHSGRTFDAQVVIIEHKSIICP GYNAVLHIHTCIEEVEITVLICLVDKKSGEKSKTRPRFVKQDQVCIARLRTAGTICLE TFKDFPQMGRFTLRDEGKTIAIGKVLKLVPEKD*A" misc_feature 31124 /locus_tag="CK280_G0032352" /note="gap added in CDS to maintain frame, possibly due to error in genome" [12] Here's another one... it has abutting intervals to adjust for an indel BUT it still has internal stop codons: CDS complement(join(7338291..7338508,7339244..7339324, 7343790..7343889,7344083..7344242,7345277..7345369, 7346224..7346286,7346807..7346925,7350256..7350420, 7351494..7351568,7358473..7358583,7359561..7359646, 7361952..7362511,7362524..7362973,7362974..7363346, 7363462..7363548,7363731..7363811,7366863..7366967, 7381325..7381744)) /locus_tag="CK280_G0032515" /note="gaps were added to CDS to maintain frame" /codon_start=1 /product="C16orf96" /protein_id="CK280:T0101780_prot" /translation="- PLHSCLGNRERLCHTHYTHTHTHTHTHTRKWNTEGGLWMRPEMNVVLQILFFFLRRSL ALSPRLECSGTILAH*NLHLPDLSNSLPQPPE*LGLQARATTPG*FLYF**RRGFTVL ARLVLNS*PRGPPAFASQSAGITQPTCFSVTLSVLSVFPSGIL*RALLSACMFIIIIS KKKVGSEGWRGGSRL*SQHFGRLRWADCTPAWATE*DSVSKKKKKKKEKKNNALLASA SLVQVILLPQPPE*LGLQVQATTPS*FLCF*WRWGLHHIGQAGPELLNS*STRLSIPK CWDYRREPLCPALLEILNVLLLIQGSGKSCSTETC*NFIYLFIYLLIFETKSHSVAQA GVQWRHLSSLQPLPPGFK*FSCLSFLIAGITGACHQAQLIFVFLVGTVFHHVGQAAVE LLTSSYPTTLTSQSAGVTGMRHCARPDMLSSLAFPKLVCPWNPFPLLGKHWIWNPFFL LKLLVT*SHRKRQGWDRTWTPSPQTQFL*QVHPERMDIFAAQNWKMVALQREVVRATI PVFPHSPWDPPTSRVLGAEGREVQGGEHSLERFLCRCGDTAPSRLTALIFYPHLAGFS PE*V*NHPQNRGHGALEWPP*CHVHLSE*GRKGRLAGSWGAGGRLTESGRRGWKPAEP LCPLFCRKLVDHRWTCGSL*SSSQRQPWPRPPSTLKLLVPSRSPSPSKTPSYCCSGSR