New assemblies and annotations for multiple glyma collections

[jd-campbell]

---

name: New genomes+annotations checklist
about: Checklist of tasks associated with new genome and annotation collections

Main steps for adding new genome and annotation collections

Genus/species/collection names:

What are the collection types and names? Example:

• Glycine/max/genomes/Enrei.gnm3.MF2C

• Glycine/max/annotations/Enrei.gnm3.ann1.YRG3

• Glycine/max/genomes/Fukuyutaka.gnm3.424H

• Glycine/max/annotations/Fukuyutaka.gnm3.ann1.RWJM

• Glycine/max/genomes/Harosoy.gnm3.TKJB

• Glycine/max/annotations/Harosoy.gnm3.ann1.R35Y

• Glycine/max/genomes/Himeshirazu.gnm3.61YZ

• Glycine/max/annotations/Himeshirazu.gnm3.ann1.0KJD

• Glycine/max/genomes/Jack.gnm3.05MK

• Glycine/max/annotations/Jack.gnm3.ann1.S50K

• Glycine/max/genomes/Kosuzu.gnm3.9GPW

• Glycine/max/annotations/Kosuzu.gnm3.ann1.4MVB

• Glycine/max/genomes/Miyagishirome.gnm3.5Z0W

• Glycine/max/annotations/Miyagishirome.gnm3.ann1.P8QZ

• Glycine/max/genomes/Ooyachi12.gnm3.3469

• Glycine/max/annotations/xOoyachi12.gnm3.ann1.6929

• Glycine/max/genomes/Peking.gnm3.KKCD

• Glycine/max/annotations/Peking.gnm3.ann1.D1WG

• Glycine/max/genomes/Tanba.gnm3.4HTT

• Glycine/max/annotations/Tanba.gnm3.ann1.17JT

• Glycine/max/genomes/UA4805.gnm3.THQ0

• Glycine/max/annotations/UA4805.gnm3.ann1.DTQV

• Add collection(s) to the Data Store, including commits to datastore-metadata

• Validate the README(s)

• Update about_this_collection.yml

• Calculate AHRD functional annotations

• Calculate gene family assignments (.gfa)

• Add to pan-gene set

• Load relevant mine

• Add BLAST targets

• Incorporate into GCV

• Update the jekyll collections listing

• Update browser configs

• run BUSCO

• Update DSCensor

• Add LINKOUTS to datastore, refresh linkout service

[adf-ncgr]

Just a heads-up that the AHRD (old-school version) for all of these is nearing completion, at which point I'll update the versions in the annex on ceres. I think I saw a message suggesting that someone on the ISU side handle the gene family assignments; let me know if that's not the case, happy to do it if needed.

[StevenCannon-USDA]

Just a heads-up that the AHRD (old-school version) for all of these is nearing completion, at which point I'll update the versions in the annex on ceres. I think I saw a message suggesting that someone on the ISU side handle the gene family assignments; let me know if that's not the case, happy to do it if needed.

The gene family assignments -- I think @jd-campbell will do this -- and also about_this_collection.yml

[jd-campbell]

@StevenCannon-USDA @adf-ncgr -- I noticed I have made a mistake with the cultivar Ooyachi12 -- the name of the cultivar is supposed to be Ooyachi**-2**

I will be correcting my mistake.

DONE.

[adf-ncgr]

I've updated the gff files in the annex with AHRD descriptors and added the BUSCO results for genomes and annotations. I think this is traditionally the point where we'd move things from the annex to the datastore proper and add the metadata files to the github (I think checksums still need to be generated). Let me know if you'd like me to do any of that. Can proceed with downstream steps (as time permits) when they're in place.

[jd-campbell]

@adf-ncgr
You mentioned during our last LIS/PB meeting that there was an issue with one of the files for these glyma genomes. I vaguely remember you mention that Sam validator script did not like something. Can you please refresh my memory?

[adf-ncgr]

Hi @jd-campbell; the immediate problem that the validator code is having is reported as:

[convertFile] ## INVALID: glyma.Hwangkeum.gnm1.ann1.1G4F.gene_models_main.gff3.gz 898 gene records are missing the Note attribute.
...
[convertFile] ## INVALID: glyma.Zh13.gnm2.ann1.FJ3G.gene_models_main.gff3.gz 3740 gene records are missing the Note attribute.

The reason these genes are Note-less is that they did not have proteins in the files used by AHRD, which is in most cases because they are actually non-coding genes. Sam had set a hard-coded value of 100 as the maximum number of Note-less genes to tolerate, and I am now remembering that Sven updated the validation code so that we could be selectively more permissive if desired: ncgr/java@064ae85
So, I could in theory use that selective permissiveness as a workaround, but I think we more or less agreed that it would be better for several reasons to separate non-coding genes into their own gff file.

As much for my reference as yours, here's what appears to be the full set of files currently with non-zero numbers of Note-less genes; it's likely that some of them are noteless for other reasons than being non-coding, e.g. the Arachis hypogaea ones are things Sudhansu made via liftoff but we did not run through AHRD (yet).

[adf@dal /falafel/legumeinfo/rsyncd_datastore/v2]$for f in */*/annotations/*/*gene_models_main.gff3.gz; do echo -n "$f "; zgrep '   gene    ' $f | grep -vc 'Note='; done | grep -v ' 0$'
Arachis/duranensis/annotations/V14167.gnm1.ann1.cxSM/aradu.V14167.gnm1.ann1.cxSM.gene_models_main.gff3.gz 1
Arachis/hypogaea/annotations/C76-16.gnm1.ann1.5MP1/arahy.C76-16.gnm1.ann1.5MP1.gene_models_main.gff3.gz 78394
Arachis/hypogaea/annotations/F07.gnm1.ann1.1DJN/arahy.F07.gnm1.ann1.1DJN.gene_models_main.gff3.gz 80079
Arachis/hypogaea/annotations/ICG1471.gnm1.ann1.BK1B/arahy.ICG1471.gnm1.ann1.BK1B.gene_models_main.gff3.gz 75344
Arachis/hypogaea/annotations/Tifrunner.gnm1.ann1.CCJH/arahy.Tifrunner.gnm1.ann1.CCJH.gene_models_main.gff3.gz 4
Cicer/arietinum/annotations/ICC4958.gnm2.ann1.LCVX/cicar.ICC4958.gnm2.ann1.LCVX.gene_models_main.gff3.gz 153
Faidherbia/albida/annotations/WAFC.gnm1.ann1.RTP9/faial.WAFC.gnm1.ann1.RTP9.gene_models_main.gff3.gz 28979
Glycine/max/annotations/Enrei.gnm3.ann1.YRG3/glyma.Enrei.gnm3.ann1.YRG3.gene_models_main.gff3.gz 4735
Glycine/max/annotations/Fukuyutaka.gnm3.ann1.RWJM/glyma.Fukuyutaka.gnm3.ann1.RWJM.gene_models_main.gff3.gz 3030
Glycine/max/annotations/Harosoy.gnm3.ann1.R35Y/glyma.Harosoy.gnm3.ann1.R35Y.gene_models_main.gff3.gz 3000
Glycine/max/annotations/Himeshirazu.gnm3.ann1.0KJD/glyma.Himeshirazu.gnm3.ann1.0KJD.gene_models_main.gff3.gz 2913
Glycine/max/annotations/Hwangkeum.gnm1.ann1.1G4F/glyma.Hwangkeum.gnm1.ann1.1G4F.gene_models_main.gff3.gz 898
Glycine/max/annotations/Jack.gnm3.ann1.S50K/glyma.Jack.gnm3.ann1.S50K.gene_models_main.gff3.gz 2930
Glycine/max/annotations/JD17.gnm1.ann1.CLFP/glyma.JD17.gnm1.ann1.CLFP.gene_models_main.gff3.gz 82
Glycine/max/annotations/Kosuzu.gnm3.ann1.4MVB/glyma.Kosuzu.gnm3.ann1.4MVB.gene_models_main.gff3.gz 2148
Glycine/max/annotations/Miyagishirome.gnm3.ann1.P8QZ/glyma.Miyagishirome.gnm3.ann1.P8QZ.gene_models_main.gff3.gz 2237
Glycine/max/annotations/Ooyachi2.gnm3.ann1.6929/glyma.Ooyachi2.gnm3.ann1.6929.gene_models_main.gff3.gz 2308
Glycine/max/annotations/Peking.gnm3.ann1.D1WG/glyma.Peking.gnm3.ann1.D1WG.gene_models_main.gff3.gz 3130
Glycine/max/annotations/Tanba.gnm3.ann1.17JT/glyma.Tanba.gnm3.ann1.17JT.gene_models_main.gff3.gz 2187
Glycine/max/annotations/UA4805.gnm3.ann1.DTQV/glyma.UA4805.gnm3.ann1.DTQV.gene_models_main.gff3.gz 3219
Glycine/max/annotations/Wm82_NJAU.gnm1.ann1.KM71/glyma.Wm82_NJAU.gnm1.ann1.KM71.gene_models_main.gff3.gz 1
Glycine/max/annotations/Zh13.gnm1.ann1.8VV3/glyma.Zh13.gnm1.ann1.8VV3.gene_models_main.gff3.gz 127
Glycine/max/annotations/Zh13.gnm2.ann1.FJ3G/glyma.Zh13.gnm2.ann1.FJ3G.gene_models_main.gff3.gz 3740
Lens/culinaris/annotations/CDC_Redberry.gnm2.ann1.5FB4/lencu.CDC_Redberry.gnm2.ann1.5FB4.gene_models_main.gff3.gz 637
Lupinus/albus/annotations/Amiga.gnm1.ann1.3GKS/lupal.Amiga.gnm1.ann1.3GKS.gene_models_main.gff3.gz 3213
Medicago/truncatula/annotations/A17.gnm5.ann1_6.L2RX/medtr.A17.gnm5.ann1_6.L2RX.gene_models_main.gff3.gz 6631
Medicago/truncatula/annotations/HM004.gnm1.ann1.2XTB/medtr.HM004.gnm1.ann1.2XTB.gene_models_main.gff3.gz 1274
Medicago/truncatula/annotations/HM010.gnm1.ann1.WV9J/medtr.HM010.gnm1.ann1.WV9J.gene_models_main.gff3.gz 1298
Medicago/truncatula/annotations/HM022.gnm1.ann1.6C8N/medtr.HM022.gnm1.ann1.6C8N.gene_models_main.gff3.gz 1292
Medicago/truncatula/annotations/HM023.gnm1.ann1.WZN8/medtr.HM023.gnm1.ann1.WZN8.gene_models_main.gff3.gz 1275
Medicago/truncatula/annotations/HM034.gnm1.ann1.YR6S/medtr.HM034.gnm1.ann1.YR6S.gene_models_main.gff3.gz 1261
Medicago/truncatula/annotations/HM050.gnm1.ann1.GWRX/medtr.HM050.gnm1.ann1.GWRX.gene_models_main.gff3.gz 1284
Medicago/truncatula/annotations/HM056.gnm1.ann1.CHP6/medtr.HM056.gnm1.ann1.CHP6.gene_models_main.gff3.gz 1298
Medicago/truncatula/annotations/HM058.gnm1.ann1.LXPZ/medtr.HM058.gnm1.ann1.LXPZ.gene_models_main.gff3.gz 1279
Medicago/truncatula/annotations/HM060.gnm1.ann1.H41P/medtr.HM060.gnm1.ann1.H41P.gene_models_main.gff3.gz 1270
Medicago/truncatula/annotations/HM095.gnm1.ann1.55W4/medtr.HM095.gnm1.ann1.55W4.gene_models_main.gff3.gz 1281
Medicago/truncatula/annotations/HM125.gnm1.ann1.KY5W/medtr.HM125.gnm1.ann1.KY5W.gene_models_main.gff3.gz 1278
Medicago/truncatula/annotations/HM129.gnm1.ann1.7FTD/medtr.HM129.gnm1.ann1.7FTD.gene_models_main.gff3.gz 1277
Medicago/truncatula/annotations/HM185.gnm1.ann1.GB3D/medtr.HM185.gnm1.ann1.GB3D.gene_models_main.gff3.gz 1271
Medicago/truncatula/annotations/HM324.gnm1.ann1.SQH2/medtr.HM324.gnm1.ann1.SQH2.gene_models_main.gff3.gz 1212
Phaseolus/acutifolius/annotations/W6_15578.gnm2.ann1.LVZ1/phaac.W6_15578.gnm2.ann1.LVZ1.gene_models_main.gff3.gz 1
Trifolium/pratense/annotations/HEN17-A07.gnm1.ann1.LHG2/tripr.HEN17-A07.gnm1.ann1.LHG2.gene_models_main.gff3.gz 2
Vigna/angularis/annotations/Gyeongwon.gnm3.ann1.3Nz5/vigan.Gyeongwon.gnm3.ann1.3Nz5.gene_models_main.gff3.gz 3401

We probably need a better way of extracting non-coding genes than merely looking at Note-less-ness, but wanted to give you a response about this sooner rather than later (and refreshed my memory about a few things in the process). Let me know if you have any immediate questions, but I'll try to propose something more concrete as far as the separation method is concerned....

[adf-ncgr]

@jd-campbell just wanted to let you know I started working on a script to split the problematic gffs into coding and non-coding genes. The strategy I used initially was to look for genes with mRNA children and put those in one pile, while those without mRNA children would go into the other.

Although that seems to do the needful in some cases, it looks like the situation is a little more complex than I'd imagined. For example, applying that strategy to the glyma.Enrei.gnm3.ann1.YRG3.gene_models_main.gff3.gz file which has 4735 Note-less genes, we only get 2518 of them which don't have mRNA subfeatures (and which look like they are probably noncoding, given that they have an attribute Method=cmscan, although they don't actually have any subfeatures like tRNA or ncRNA, which is a little odd).

The remainder of the Note-less genes in this case appear to be genes with mRNA subfeatures but without CDS like:

glyma.Enrei.gnm3.chr01  NARO    gene    796822  796996  .       -       .       ID=glyma.Enrei.gnm3.ann1.GlymaJER.01G300001.1
glyma.Enrei.gnm3.chr01  NARO    mRNA    796822  796996  .       -       .       ID=glyma.Enrei.gnm3.ann1.GlymaJER.01G300001.1.t1;Parent=glyma.Enrei.gnm3.ann1.GlymaJER.01G300001.1
glyma.Enrei.gnm3.chr01  NARO    exon    796822  796996  .       -       .       ID=glyma.Enrei.gnm3.ann1.GlymaJER.01G300001.1.t1.exon;Parent=glyma.Enrei.gnm3.ann1.GlymaJER.01G300001.1.t1

In fact this particular mRNA doesn't appear in the fasta for mRNA, though I've seen other cases (e,g. lentil annotations) where such things do appear in the mRNA fasta but not the cds/protein fastas. Any thoughts from you or @StevenCannon-USDA about what we'd want to do with such non-canonical genes? I'm starting to wonder if separating things into separate gff files is better than just adding attribute tags at the gene level to indicate to other applications whether or not they should include them (for example GCV could easily only import genes that were tagged as gene_biotype=protein_coding, using the Ensembl/NCBI convention, and intermine could only require Notes on such genes). An oddball gene like the one above could either not get a gene_biotype tag at all or be given something indicating its questionable status.

The remainder of this is mostly a note to self:
In other cases from the list of annotations with Note-less genes, the lack of a Note feature seems to have resulted from an issue with converting protein ids back to gene ids in cases where it looked like a standard suffixing convention was being followed, but there were some non-compliant genes with odd transcript/protein ids associated, e.g. one of the medicago sets has most transcripts with ids like medtr.HM129.gnm1.ann1.g299.t1 whereas some look like medtr.HM129.gnm1.ann1.crp4820_scf0000_0M_1. I guess many of these types of issues got swept under the rug before Sam was being quite so rigorous about this specific validation, and I'm likely to just wait for the newer AHRD to get run to take care of those ones.