About Genecode version of o3d

[HuZhuocheng]

Dear authors,

Thank you very much for the software developed by you and your team. During use, I encountered an issue: many of the results I expected did not appear. The main reason is the mismatch between my Transcript_ID and O3D_transcript_ID. For example, regarding the TP53 gene, since I used version v47 for annotation, the ID of TP53 in my data is ENST00000504290, while in O3D it is ENST00000269305.

I would like to know which version of Genecode O3D uses, so that I can re-perform the annotation. Alternatively, could you please tell me if a parameter can be added to oncodrive3d run to ignore the mismatch and use hugo_symbol, such that it is considered acceptable as long as the gene names match?

[St3451]

Hi @HuZhuocheng , thanks again for reaching out!

You already nailed the root cause: your annotations tag TP53 as ENST00000504290, while the O3D build expects ENST00000269305. Our default human datasets come from the AlphaFold Protein Database human proteome (v4) mapped to the Ensembl release 111 / GENCODE v45, where TP53 is defined on ENST00000269305. This transcript is still present in GENCODE v47 and it is still flagged as the canonical one, so the mismatch is most likely driven by how your pipeline selects transcripts (which isoform it keeps) rather than by the GENCODE version itself. That said, using the same Ensembl/GENCODE version as O3D is still ideal, because other genes may not align as cleanly and can introduce additional mismatches. With the standard run settings the gene should still be processed, but the output flags it as Transcript_status = Mismatch, which is what you’re seeing.

To resolve this, you can try the following:

1. Re-annotate against the same Ensembl/GENCODE release used by O3D (Ensembl 111 / GENCODE v45), ensuring that transcript IDs in your VEP output match the O3D datasets. This is not strictly required for TP53 alone (since its canonical transcript exists in v47), but it reduces the risk of transcript mismatches for other genes.
2. Keep your current annotations and rerun O3D on the unfiltered VEP output (including all transcripts) while adding
   --o3d_transcripts and --use_input_symbols. These flags tell the pipeline to use its known transcript set while refreshing gene symbols from your input.
3. Check (and, if possible, adjust) your transcript filtering step. After annotation, verify whether your pipeline is filtering transcripts (for example, keeping only the canonical transcript or the first listed one). For TP53, the canonical transcript should be ENST00000269305, so your filtering step should ideally retain this one. If you can modify the transcript-filtering rules, we recommend keeping the canonical transcript for each gene so that the IDs are consistent with the O3D datasets.

Any of these approaches should help clear the mismatch warning and recover the TP53 results you expect. In practice, whichever option you choose, I still recommend double-checking how transcripts are filtered in your pipeline and, if possible, running O3D on unfiltered VEP output with --o3d_transcripts (and --use_input_symbols) so that the expected transcript IDs are not dropped by accident.

Show more

Just to provide more context: that transcript mismatch flag doesn’t filter TP53 out by itself. However, if the transcript change alters the protein sequence (so the residues no longer align), per-position mapping may start dropping rows. O3D cross-checks the wild-type amino acid in your input against the sequence used in the structural annotations and removes mismatching entries by default. Only when the mismatch ratio exceeds `--thr_mapping_issue` (10 % by default) is the entire gene discarded. Setting this threshold to 1 allows keeping mismatched entries as tagged rows instead of removing them, but it is generally safer to keep the default behaviour so discrepancies are flagged and filtered.

By the way, that flag isn’t documented in the README yet, I’m updating it now so this behaviour is easier to trace.

Details

[HuZhuocheng]

Dear Author, thank you very much for your detailed reply. First, I would like to express my gratitude for your thorough response to my questions—such work is highly commendable. With your help, I have successfully resolved this Mismatch issue. Since I used Annovar instead of VEP for annotation, I modified the column names to meet the software's input requirements. However, I have two further questions and would like to seek your advice:

(1) First, I noticed you suggested using unfiltered VEP results. Does this mean I should also input the full, complete results after variants calling here? As my research focuses on somatic mutations, I filtered out variants with excessively high population frequencies, those with insufficient depth and variant allele frequency, and those not marked as "PASS" in the GATK pipeline. After filtering each sample individually, I merged them into a single large cohort MAF file to run o3d. Do you have any valuable suggestions regarding this?

(2) Second, in the results, I found that some genes (e.g., B2M, CIC, CIITA) and most genes in the IGH family (e.g., IGHE, IGHG1, IGHG2) are marked as O3D_missing. Could you please advise on a good solution to this problem? Since I used panel sequencing, each gene is very valuable. Among the 3172 rows in the results, 560 entries are marked as O3D_missing. Is this caused by the use of v47 for my gene_symbol annotation, or are these genes inherently unavailable in O3D?

I am extremely grateful for your assistance and would deeply appreciate it if you could provide a further reply!

[St3451]

Hi @HuZhuocheng, thanks so much for the kind words! I’m very happy to hear you managed to resolve the mismatch issue.

Regarding your two additional questions:

1. On “unfiltered” VEP results and variant filtering:
   I should have been more specific when I suggested using “unfiltered” VEP results. What I meant there was not to filter at the transcript level, i.e. to keep the mapping of each variant to all possible transcripts rather than restricting to just one per gene. Since your mismatch issue is already solved, you do not need to change what you are currently doing unless you run into a new transcript mismatch issue.

   Also, it’s important to note that a transcript mismatch by itself is not necessarily a problem. It only becomes an issue if it leads to mismatches between the wild-type amino acids in your input mutation data and the residues in the structures used by O3D. When that happens, it will be recorded in the Status column of the O3D gene-level output (e.g. as one of the mismatch-related statuses), so that column is the first place to look when checking if a gene has been processed correctly.

   For the variant-level filters you describe (population frequency, depth, VAF, GATK PASS flag), I strongly recommend keeping them as you have done. These filters are important to clean the somatic mutation data (e.g., to remove likely germline variants and technical artefacts) before running Oncodrive3D. In other words, “unfiltered VEP” in my previous message referred only to transcript selection, not to skipping quality/germline filters on the variants themselves.

2. On genes marked as O3D_missing:
   I assume that the O3D_missing label you see for some genes appears in the Transcript_status column of the O3D gene-level output. This label indicates that the Oncodrive3D dataset is missing transcript information for that gene, not that the gene itself is unavailable or that your annotation with v47 is incorrect.

   Briefly, O3D maps AlphaFold-predicted structures to transcript IDs using the EMBL-EBI Protein API and the Ensembl REST API. For some structures, however, it is not possible to map them unambiguously to a specific transcript, and in those cases the transcript information is missing; this is reported as O3D_missing in Transcript_status.

   As mentioned above, the key column to check is the Status column in the O3D gene-level output. This tells you whether the gene was processed correctly or whether O3D encountered any issues when running the 3D analysis. If a gene appears in the output at all, it means that it is present in the O3D dataset; O3D_missing only refers to the missing transcript mapping, not to the gene being excluded from O3D. For more details, please also see the Output section of the “Input and Output” documentation.

I hope this clarifies both points, and thank you again for the very detailed questions and for your patience while using O3D.

[HuZhuocheng]

Thank you for your reply, which has made my understanding of o3d much clearer. Please allow me to ask one more question. Following your suggestion, I checked the results of the clustering_genes.csv file and found that some genes I am interested in, such as those from the KMT2 family, are all labeled as WT_mismatch, which seems incorrect. Examples are as follows:

Gene,Uniprot_ID,pval,qval,C_gene,C_pos,C_label,Pos_top_vol,Score_obs_sim_top_vol,Mut_in_gene,Clust_mut,Clust_res,Mut_in_top_vol,Mut_in_top_cl_vol,PAE_top_vol,pLDDT_top_vol,pLDDT_top_cl_vol,Ratio_not_in_structure,Ratio_WT_mismatch,Mut_zero_mut_prob,Pos_zero_mut_prob,Cancer,Cohort,F,Transcript_ID,O3D_transcript_ID,Transcript_status,Status
KMT2A,Q03164,,,0,,,,,12,,,,,,,,0.0,0.333,0.0,,,b1,14M,ENST00000389506,ENST00000389506,Match,WT_mismatch
KMT2C,Q8NEZ4,,,0,,,,,25,,,,,,,,0.0,0.96,0.0,,,b1,19M,ENST00000262189,ENST00000262189,Match,WT_mismatch
KMT2D,O14686,,,0,,,,,30,,,,,,,,0.0,0.833,0.0,,,b1,22M,ENST00000301067,ENST00000301067,Match,WT_mismatch

I also saw the following warnings in the output during execution:

2025-11-06 13:04:32,967 - WARNING | oncodrive3d.run.clustering - Detected 4 (33.3%) mut with ref-structure WT AA mismatch of KMT2A (Q03164-F14M, transcript status = Match): Filtering the gene
2025-11-06 13:04:32,973 - WARNING | oncodrive3d.run.clustering - Detected 24 (96.0%) mut with ref-structure WT AA mismatch of KMT2C (Q8NEZ4-F19M, transcript status = Match): Filtering the gene
2025-11-06 13:04:32,979 - WARNING | oncodrive3d.run.clustering - Detected 25 (83.3%) mut with ref-structure WT AA mismatch of KMT2D (O14686-F22M, transcript status = Match): Filtering the gene

I have 350 genes here, among which 97 are labeled with WT_mismatch. Could you please give some suggestions on how to check and resolve this issue?
Thank you again!

[St3451]

Hi @HuZhuocheng!, thanks again for reaching out.

Let me unpack what’s going on with those WT_mismatch genes and how you can check / fix it.

1. What WT_mismatch actually means

For each missense mutation, Oncodrive3D compares:

• the wild-type amino acid and position coming from your input (MAF / annotated table; using the HGVSp_Short field)
• with the wild-type amino acid at the same position in the protein sequence used internally (from the AlphaFold structures; by default settings: AlphaFold Protein Database, human proteome download v4)

If too many mutations in a gene fail this check, the gene is filtered and marked in the gene-level file as:

• Status = WT_mismatch
• Ratio_WT_mismatch = fraction of mutations in the gene whose WT AA didn’t match the reference sequence.

That’s exactly what your log messages are reporting, e.g.:

• KMT2A – 4 / 12 mutations (33.3%) WT mismatch
• KMT2C – 24 / 25 mutations (96.0%) WT mismatch
• KMT2D – 25 / 30 mutations (83.3%) WT mismatch

Because these fractions are > the default --thr_mapping_issue 0.1, the genes are filtered from the clustering step and end up with Status = WT_mismatch in clustering_genes.csv.

So:

• Transcript_status = Match: the transcript ID is in the O3D datasets and aligns cleanly.
• Status = WT_mismatch: the residues encoded in your input don’t match the residues in the used structures, for a large fraction of mutations.

2. Why this can happen (even when transcripts “Match”)

Given your description (ANNOVAR, GENCODE v47), the high number of WT mismatches in 97/350 genes strongly suggests a systematic mismatch between:

• the protein sequence that your annotation tool used to produce the protein changes, and
• the protein sequence used by Oncodrive3D (AlphaFold human proteome v4, i.e., Ensembl 111 / GENCODE v45).

The most likely cause is that you annotated using a different Ensembl/GENCODE release where the same Ensembl transcript ID exists but its protein sequence has changed (e.g. updated CDS, different start codon, small indels, etc.). Ensembl tracks these updates via the version suffix (e.g. ENST00000262189.8 → ENST00000262189.9), but many annotation tools — and Oncodrive3D itself — do not keep track of the .x suffix. As a result, different transcript versions from different GENCODE releases are collapsed to the same ID (ENST00000262189), even though their underlying protein sequences differ. Since O3D uses the sequence from Ensembl 111 / GENCODE v45, annotations generated from a newer release (e.g. v47) can produce wild-type residues that no longer match the O3D sequence, and those genes are then flagged as WT_mismatch.

3. How to check what’s going wrong (KMT2 as a test case)

Take one KMT2 gene, e.g. KMT2C:

1. Pick one mutation
   From your original input (before O3D), choose a KMT2C variant and note:

   • protein position
   • WT amino acid (from the HGVSp_Short field)

2. Look up the sequence O3D used
   In your O3D build datasets folder, look for seq_for_mut_prob.tsv and find the row with:

   • Gene = KMT2C
   • The Seq column is the exact amino acid sequence O3D used.

3. Compare WT residues
   Compare the WT AA in your annotation vs the residue at that position in Seq
   If they differ for many mutations, that explains the high Ratio_WT_mismatch and Status = WT_mismatch.

4. Check your annotation source
   In that case, your AA changes are likely coming from a different transcript/protein definition (e.g. refGene vs Ensembl, or mixed fields in ANNOVAR). Aligning the annotation source to the Ensembl/Ensembl-like definitions used in O3D should reduce these mismatches.

By the way, since you mentioned that the affected genes include several long ones (like the KMT2 family), could you check whether the WT_mismatch flag also appears for shorter genes in your results? Long genes are inherently trickier, because AlphaFold does not provide a single full-length structure for proteins with more than 2700 residues; instead, it releases them as fragments that O3D reassembles into a single protein. It would therefore be very helpful to know whether the mismatches are confined mainly to these long, fragment-based proteins or whether they also affect shorter genes.

4. How to reduce / fix the WT mismatches

Ideally, I would suggest to make the upstream annotation and O3D use the same protein definitions. The most robust way is to ensure that the protein consequences are computed from the same Ensembl/GENCODE release that Oncodrive3D uses:

• Annotate with VEP using Ensembl release 111 / GENCODE v45, and feed O3D either:

  • the resulting MAF, or
  • the VEP output, including all mapped transcripts, together with --o3d_transcripts --use_input_symbols.

• If you prefer to stay with ANNOVAR:

  • make sure you’re using Ensembl-based annotations for the protein consequences (HGVSp_Short)
  • and that the Ensembl sequence version in ANNOVAR matches (as closely as possible) the one used by O3D. Even small differences in CDS can inflate WT_mismatch for long genes like KMT2s.

After re-annotating, I would re-run O3D and check if the number of genes labeled as WT_mismatch Status reduced.

[HuZhuocheng]

Thank you very much for your reply. Thanks to your detailed explanation, I successfully resolved this issue over the weekend! Now, most of the remaining genes labeled as WT_mismatch are due to O3D_missing（from96/300->18/300）, and I believe there is no need for me to address this. I would like to share the attempts I made here to help others who may encounter the same problem when using ANNOVAR for annotation in the future:

Main Solution Adopted
I replaced the original RefGene annotation with ANNOVAR's ENS annotation.
First, I extracted the TranscriptIDs from the output list, then searched for matching ones in my annotations. For example, in ABCB1:ENST00000543898.5:exon12:c.C1135T:p.L379F, ABCB1:ENST00000622132.5:exon12:c.C1327T:p.L443F, ABCB1:ENST00000265724.8:exon13:c.C1327T:p.L443F, the second ID matched, so I used it for replacement.
If this method failed, I retained the amino acid mutation label in aachange. Special attention should be paid to the accuracy of replacements for Missense_Mutation sites.

Additionally, I would like to ask a few questions, and I hope you don't mind. Thank you very much for your work:
(1) My data is panel sequencing data. Can I consider that the level of Score_obs_sim_top_vol indicates whether a gene is under positive selection?
(2) For my panel containing only 300 genes, are there any additional running parameters that need to be adjusted?

In short, thank you very much for your help this week, which has allowed me to resolve one issue after another!

[St3451]

Hi @HuZhuocheng! It is great that you managed to solve the mapping issue. And thanks a lot for sharing your solution. We can keep this issue open so that it can be useful for other users using ANNOVAR as annotation tool.
Finally, please feel free to ask any additional questions, they are always welcome and I am very happy to help!

Question 1

First let me provide you a bit of context regarding the issue of using panel sequencing data in the detection of signals of positive selection. One of the main problem of panel sequencing data is that (depending on the size of the panel) usually there isn't enough information to accurate model neutral mutagenesis (e.g., the mutational profile of the cohort, for example the frequencies of mutations (SBS) in threenucleotide context). However, while O3D can take as input the mutation profile of the cohort, it is not essential and it works quite well even without it (i.e., assuming an uniform distribution of the mutations on the DNA sequence). Then, of course, if you sequence strictly only the exact regions where you expect to find clusters, you don't have much resolution.

Regarding your specific question, the Score_obs_sim_top_vol is the score that we assign to the gene. It represents how anomalous is to observed a given number of mutations in the volume of the residue with the highest clustering signal of the gene. It also have a direct interpretation, for example, a score of 2 literally means that the clustering score in the target volume is 2 times higher than what we would expect by chance (under neutral mutagenesis). It could be used to compare signal across residues and genes or other downstream analysis. If you are interested to know more about the clustering score, please check out the Supplementary Notes of the O3D paper published on NAR.

However, if you want to know which genes are under positive selection I suggest you to look at the qval column (p-value corrected for FDR) of the gene (which is the FDR corrected p-value of the most significant residue). For example, we used a selection threshold of 0.01. Please check out the Output section of the “Input and Output” documentation.

Question 2

Not really, if you’re not encountering specific issues and don’t have particular requirements (for example, data with highly heterogeneous sequencing depth or variants that must be mapped to predefined transcripts), you can safely run Oncodrive3D using the default configuration.

Thanks again for using Oncodrive3D and for asking interesting questions :)

[HuZhuocheng]

Dear Authors,

Thank you for your response to my questions. I would like to follow up with a few additional inquiries regarding the results. First, I greatly appreciate your impactful work.

Question 1:
In the generated clustering_genes file, I noticed a column named Mut_in_gene, which I understand records the number of mutations per gene. However, there is an inconsistency between these values and the mutation counts from my input MAF file. For example, one gene showed 17 mutations across all samples in my MAF file, but Mut_in_gene indicated 30. Conversely, another gene had over 60 mutations in my analysis, yet Mut_in_gene reported only 30. Could you please explain the reason for this discrepancy?

Question 2:
Among the genes identified as under positive selection in the results, I observed several IG (Immunoglobulin) family genes such as IGHG1 and IGLV3-19. Given that IG genes are known for high population diversity, should this factor be considered in the analysis? Would it be appropriate to exclude these genes from the final results?

Question 3:
I am unclear about the exact meaning of C_pos in the output. Additionally, I would like to obtain detailed information on the protein impact associated with this parameter. For instance, the Pos_top_vol value for the STAT6 gene is 419—what resources or tools would you recommend for further analyzing this position?

Once again, thank you for your valuable work. I would be extremely grateful for your insights into these questions.

[St3451]

Hi @HuZhuocheng, thanks for your questions and sorry for the late reply.

Question 1 – Mut_in_gene vs. your original MAF

It is surprising that the Mut_in_gene shows more mutations than the original mutations across all samples in the MAF file.
Could you please share the oncodrive3d run command you used?

Mut_in_gene is simply the number of rows that survive all internal preprocessing for that gene. During O3D run we:

1. filter to missense mutations,
2. re-map each variant to the corresponding structure (a HUGO symbol always map to a unique structure, while a structure can map to multiple HUGO symbols), and
3. drop rows that cannot be placed on the structure (outside chain, WT mismatch, etc.).

The count is therefore taken after those steps. If a variant appears twice in your input (e.g. because the MAF lists it once per transcript or per sample) it will be counted twice; if a variant is discarded because it cannot be mapped, the count goes down. You can inspect exactly what O3D used by looking at the <cohort>.mutations.processed.tsv (written in the output folder). The sums per gene in that file match the Mut_in_gene column.

Question 2 – IG genes

As you mentioned, immunoglobulin loci are notoriously hyper-variable and this could result in false positive hits just because of their underlying biology. We intentionally leave them in: some users study immune repertoires and want those signals. If your application isn’t focused on IG genes, it’s perfectly fine to filter them after the run. This doesn’t affect the rest of the output.

Question 3 – C_pos, Pos_top_vol, and interpretation

In the gene-level CSV, C_pos lists the residues whose mutation density within the residue and its neighboring residues (“volume”) was called significant (C = 1). Pos_top_vol points to the single residue whose cluster shows the strongest observed vs. expected enrichment, e.g., Pos_top_vol = 419 in your STAT6 example.

To understand what’s happening at that position you could:

• Inspect the residue-level file <cohort>.3d_clustering_pos.csv and look up STAT6, Pos = 419. You’ll see how many mutations land there (Mut_in_res, Mut_in_vol) compared to the total mutations in the gene (Mut_in_gene), the clustering score, and the p-value.
• Also, you could run oncodrive3d plot (after you have run oncodrive3d build-annotations). The per-gene plots overlay mutation counts, missense probabilities, clustering scores, domains, PTMs, etc., so you can visualize where residue 419 sits on the protein, visually see how many mutations contributed to that cluster, and overlap the clustering signal with annotations (e.g., domains). The annotated CSV that oncodrive3d plot generates contains the same information in tabular form. Those visuals make it easier to judge whether the signal is a plausible structural cluster or more likely an artifact. Please check out the documentation of the annotation and plotting features of Oncodrive3D and let us know if you have any problem to set it up.
• It’s also good practice to do a literature search around that residue to explore whether the residue overlaps known functional elements.

Thanks again, and please let us know if you have any updates on the issues, further questions or run into any other issues.