Reference gene count variation

[dhooper1]

Hi TOGA2 team,

Thanks first for a great resource! I was wondering if you could help explain why the number of reference genes might differ between query assemblies run against the same reference. In the summary output files generated by TOG2, two examples attached here, the total number of reference genes considered present plus those considered lost/missing is different across queries.

summary.apus_apus.txt

summary.hemiprocne_comata.txt

I am trying to compare results between queries but am having trouble proceeding when the denominator (i.e., the total number of reference genes) differs between queries. In the examples shared above there were 16,071 reference genes when Apus_apus was used as query and 16,139 when Hemiprocne_comata was used as query. Both TOGA2 runs used the same reference genome, annotation, and isoform files and both query assemblies were aligned to the same reference assembly.

Many thanks for your time and guidance!

[ReverendCasy]

Hello @dhooper1 ,
Thank you! Could you please check how many reference genes the loss_summary.tsv files contain in both cases? Here’s a fast way to do it:

grep -c GENE ${your_output_directory}/loss_summary.tsv

If the number is the same, I suspect a minor hiccup occurred when generating the summary file; otherwise it’s a runtime bug. Could you please specify the TOGA2 version you used?
Also, I do not see the reference gene numbers you specified (16,071 and 16,139) in the respective summary files. Could you please specify where you got them from?

Best,
Yury

[dhooper1]

Hi @ReverendCasy,

Thanks so much for the speedy response! I am using TOGA2 v2.0.6.

Both example cases have the same number of GENE's in the loss_summary.tsv file: 16,799.

The numbers I referred to came from summing the counts for "#reference genes considered present" and "#reference genes considered lost/missing" in the 'Reference genes level' in each summary.txt file as their proportions for reference genes sum to 100%.

[dhooper1]

Hi again @ReverendCasy,

I think I understand better now. The GENE counts for each loss status in the loss_summary.tsv file match what is observed in the summary.txt file but the summary isn't counting the number of GENEs with loss status PG or N in the loss_summary.tsv file. So when these counts differ between queries, the denominator (i.e., the number of reference GENEs) appears to differ.

In practice, what would be the best way to incorporate these PG and N calls for comparative purposes between queries (i.e., for generating the classic TOGA summary plot comparing different assemblies' counts for genes with intact reading frames, missing genes, and genes with inactivating mutations)? I would assume best to count PGs as genes with inactivating mutations but am not clear on how 'N's are defined - as missing? - can you point me to the working definition?

Many thanks!

[ReverendCasy]

Hi @dhooper1 ,

I think I understand better now. The GENE counts for each loss status in the loss_summary.tsv file match what is observed in the summary.txt file but the summary isn't counting the number of GENEs with loss status PG or N in the loss_summary.tsv file. So when these counts differ between queries, the denominator (i.e., the number of reference GENEs) appears to differ.

Good catch, thank you. I am marking the issue as bug, the fix likely going to the v2.0.7b release.

In practice, what would be the best way to incorporate these PG and N calls for comparative purposes between queries (i.e., for generating the classic TOGA summary plot comparing different assemblies' counts for genes with intact reading frames, missing genes, and genes with inactivating mutations)?

For high-level analysis, you can treat both categories as Missing; PG do not have any confirmed orthologs, and N do not have any annotated items at all but in both cases there is not enough evidence for TOGA2 to state that the genes were lost.

but am not clear on how 'N's are defined - as missing? - can you point me to the working definition?

The documentation on loss status classification is on its way; in the meanwhile, please see the rejection file breakdown page for details on which items are classified as N. In brief, N (Discarded) stands for items that were discarded for technical reasons, such as frameshifts in the reference or illegal sequence names.

Best,
Yury

[dhooper1]

Thank you so much for your advice and help with this, Yury. I will keep an eye out for the bug fix in the next release and for more details on loss status classification in the weeks ahead.

Thanks again!
Daniel