Oronasal mucosal melanoma is defined by two transcriptional subtypes in humans and dogs with implications for diagnosis and therapy
Mucosal melanoma (MM) is a rare and aggressive form of melanoma arising from melanocytes in sun‐protected mucosal surfaces, including oral, nasal, genital, and anorectal regions . Despite accounting for a minority (1.4%) of melanoma cases, MM carries a disproportionately higher morbidity and mortality risk than cutaneous melanoma (CM) . MM is often diagnosed at an advanced stage, exhibits a higher likelihood of metastasis , and a reduced 5‐year survival . Unlike CM, MM lacks well‐defined precursor lesions and has few oncogenic driver mutations or effective therapeutics . MM from different anatomical sites is mutationally heterogeneous, although interestingly, there is some evidence for ‘upper’ and ‘lower’ body site‐specific mutational profiles . We reasoned that exploring the nuances of anatomical site‐specific MM, with a specific focus on oronasal mucosal melanoma (OMM), had the potential to disentangle the inherent heterogeneity of this rare melanoma subtype. There are few experimental models for OMM; however, OMM is a common and naturally occurring oral tumour in dogs . Previous research identified two molecular subgroups in OMM in dogs based on RNA expression, but the clinical significance of these subgroups for veterinary oncology is unknown . Despite dogs being proposed as a model for human OMM, it is unknown whether these subtypes even exist in human OMM and whether human OMM is analogous to the canine disease. These outstanding questions are critical because the lack of diagnostic subtypes precludes subtype classification that could be pivotal for patient‐centric precision medicine and could inform clinical trial design to account for potential subtype‐specific responses . Here, for the first time, we present transcriptomic and 16S sequencing data side by side for two new cohorts of human and canine OMM. Our analysis reveals OMM tumour stratification for both human and canine patients based on a CTLA4‐high or MET‐high transcriptional signature, defined by 41 genes, shared between species. This patient stratification is facilitated by an immunohistochemical tool, pivotal not only for individualised patient care but also for informing clinical trial design.   Ethics statement and sample collection This study was conducted with the approval of the University of Edinburgh's Veterinary Ethical Review Committee (Reference no. 10.21) and Lothian BioResource Access Committee (REC reference 20/ES/0061 SR1725). Formalin‐fixed, paraffin‐embedded (FFPE) human and canine OMM tissue samples and associated clinical data were collected from patients with naturally occurring tumours who had undergone surgical resection or biopsy of the primary tumour. Tissue was only included in the study with the informed, written consent of the human patient or the caregiver of the dog who bore the tumour. The treatment that a human or canine patient received was unaffected by their inclusion in this study. FFPE tissue samples of treatment‐naïve human OMMs were collected between 2006 and 2021. Clinicopathological data collected included patient gender, ethnicity, age (years), primary tumour location, presence of ulceration, tumour depth (mm), World Health Organization (WHO) stage, metastatic status, provision of adjuvant therapy, and duration of patient survival or lost to follow‐up (in months following initial diagnosis of OMM). FFPE tissue samples of treatment‐naïve canine OMMs were collected between 2011 and 2021 for histopathology from dogs attending the Clinical Veterinary Oncology departments at the Universities of Edinburgh, Bristol, Liverpool, Glasgow and Ohio Small Animal Teaching Hospitals, Vets Now (Glasgow) Referrals, Paragon Referrals, The Ralph Veterinary Referral Centre, Anderson Moores Veterinary Specialists, and North Downs Specialist Referrals. Clinicopathological data collected included sex, breed, age (years), primary tumour anatomical location, melanoma type (melanotic/amelanotic, pleomorphic/composite/nevocytoid/animal, ulcerated, exophytic), depth of invasion (mm), stroma:tumour ratio , mitotic count , metastatic status (confirmed by abdominal ultrasound, thoracic computed tomography or radiography, and/or histological examination of ≥1 regional lymph node), WHO stage, provision of adjuvant therapy, and duration of patient survival or lost to follow‐up (in months following initial diagnosis of OMM). The diagnosis of OMM was confirmed by a consultant histopathologist (MM) or board‐certified veterinary anatomic pathologist (JdP), and the regions of the FFPE blocks that represented OMM tissue tumour identified. RNA and DNA extraction Using a manual tissue arrayer (Beecher MTA‐1, Estigen OÜ, Tartu, Estonia) machine, two cores (1 mm outer diameter, 5 mm length) were harvested from each canine FFPE block: one from an area of OMM and one from histopathologically normal tissue. RNA and DNA were isolated from the tissue cores using a Covaris E220 Evolution Focused Ultrasonicator (Covaris Ltd, Woddington, Brighton, UK, Catalogue No. 500429) and the truXTRAC® FFPE total NA (Nucleic Acid) Kit – Column (Covaris, Catalogue No. 520220). Inhibitors (e.g. melanin) were removed using the OneStep PCR Inhibitor Removal kit (Zymo Research, Irvine, CA, USA, Catalogue No. D6030). Quality control DNA (canine only) and total RNA (human and canine) were assessed for quality and integrity using a Fragment Analyser Automated Capillary Electrophoresis System (Agilent Technologies, Santa Clara, CA, USA, Catalogue No. 5300) with the Genomic DNA 50 kb Kit (Agilent Technologies, Catalogue No. DNF‐467‐0500) and Standard Sensitivity RNA Analysis Kit (Agilent Technologies, Catalogue No. DNF‐471‐0500) respectively, and then assessed using a Qubit 2.0 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA, Catalogue No. Q32866) and the Qubit DNA (Thermo Fisher Scientific, Catalogue No. Q32853) or RNA (Thermo Fisher Scientific, Catalogue No. Q10210) broad range assay kits. In the RNA samples, DNA contamination was quantified using a Qubit dsDNA HS assay kit (Thermo Fisher Scientific, Catalogue No. Q32854). The Fragment Analyser data showed a high degree of RNA degradation, so the fragmentation step was omitted from the library preparation protocol. RNA library preparation First‐strand cDNA was generated from 50 ng of each total RNA sample using a SMARTer® Stranded Total RNA‐Seq Kit version 2 – Pico Input Mammalian kit (Clontech Laboratories, Mountain View, CA, USA, Catalogue No. 634411). Illumina‐compatible adapters and indexes were then added via five cycles of PCR. AMPure XP beads (Beckman Coulter, Brea, CA, USA, Catalogue No. A63881) were then used to purify the cDNA library. Depletion of ribosomal cDNA (cDNA fragments originating from highly abundant rRNA molecules) was performed using ZapR version 2 and R‐probes version 2 specific to mammalian ribosomal RNA and human mitochondrial rRNA. Uncleaved fragments were then enriched by 13 cycles of PCR before a final library purification using AMPure XP beads. Libraries were quantified with the Qubit dsDNA HS assay and assessed for quality and size distribution of library fragments using the Fragment Analyser and a NGS Fragment Kit (Agilent Technologies, Catalogue No. DNF‐473‐0500). DNA library preparation Libraries were prepared from up to 300 ng of each DNA sample using the Sure Select XT Target Enrichment for Illumina System according to the provided protocol for FFPE samples. DNA samples were sheared to an average fragment size of 150–200 bp using a Covaris E220 Evolution and purified using Agencourt AMPure XP beads (Beckman Coulter, Catalogue No. A6881). Purified DNA fragments were end‐repaired to remove 3’‐ and 5’ overhangs before further purification with AMPure XP beads. A single ‘A' nucleotide was added to the 3’ ends of the blunt fragments to prevent them from ligating to another during the subsequent adapter ligation reaction, and a corresponding single ‘T' nucleotide on the 3’ end of the adapter provided a complementary overhang for ligating the adapter to the fragment. dA‐tailed DNA was again purified using AMPure XP beads before paired‐end adapters were ligated to the ends of the dA‐tailed DNA fragments to prepare them for hybridisation in a flow cell. After purification using AMPure XP beads, 13 cycles of PCR were used to selectively enrich those DNA fragments that had adapter molecules on both ends and amplify the amount of DNA in the library suitable for target enrichment. Amplified gDNA libraries were purified using AMPure XP beads. gDNA libraries were quantified by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser using a DNA HS Kit (Agilent Technologies, Catalogue No. 5067‐4626). Seven hundred fifty ng of each gDNA library (where available) was hybridised to target‐specific probes (Canine All‐Exon) before target molecules were captured on streptavidin beads and non‐target molecules removed with a series of washes. Individual hybridised and captured libraries were then amplified using 12 cycles of PCR (with unique indexing primers to allow multiplexed sequencing) and purified with AMPure XP beads. Indexed libraries were measured by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser with a DNA HS Kit. DNA concentration and fragment size data were used to calculate molarity for sequencing. Sequencing Sequencing was performed using the NextSeq 2000 platform (Illumina, San Diego, CA, USA, Catalogue No. 20038897) with NextSeq 2000 P3 Reagents (200 Cycles) (Illumina, Catalogue No. 20040560). Libraries were combined in equimolar pools based on Qubit and bioanalyser assay results, and each pool was sequenced in a P3 flow cell. For RNA sequencing (RNA‐seq), PhiX Control version 3 (Illumina, Catalogue No. FC‐110‐3001) was spiked into each run at a concentration of 1% to allow troubleshooting in the event of any issues. Base call data produced by the NextSeq 1000/2000 Control Software (version 1.4.1.39716) were automatically converted into FASTQ files and uploaded to BaseSpace. Immunohistochemistry Human OMM Immunohistochemistry was conducted on 3‐μm‐thick FFPE human OMM tissue using a Leica Bond immunostainer (Leica Bond III Autostainer, Leica, Sheffield, UK, Catalogue No. 21.2201). A standard protocol was followed in accordance with the manufacturer's instructions. Anti‐CD3 was used to identify T cells and natural killer (NK) cells (1:500 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Sheffield, UK, Catalogue No. NCL‐L‐CD3‐565), Anti‐CD19 was employed to stain B‐lymphocytes (1:50 Novocastra Liquid Mouse Monoclonal Antibody CD19, Leica Biosystems, Catalogue No. NCL‐L‐CD19‐163), Anti‐CD68 was used to label monocytes and macrophages (1:50 Monoclonal Mouse Anti‐Human CD68, Dako, Santa Clara, CA, USA, Catalogue No. M0876). Tissue was also stained for CD163 (1:50, CD163 Monoclonal Antibody, Thermo Fisher Scientific, Catalogue No. MA5‐11458), CTLA4 (1:50 CTLA4 Polyclonal Antibody, Invitrogen, Darmstadt, Germany, Catalogue No. PA5‐115060), and MET (1:400 HGFR/c‐MET Antibody, 1G7NB, Novus Biologicals, Abingdon, Oxford, UK, Catalogue No. 44306SS) using the same protocol. Canine OMM Immunohistochemistry was conducted on 4‐μm‐thick FFPE canine OMM tissue, cut using a microtome (Thermo Fisher Scientific Rotary Microtome Microm HM 340 E). All tissue sections, including relevant controls, were deparaffinised in xylene and rehydrated in graded ethanol through to distilled H 2 O prior to staining. Heat‐induced antigen retrieval was performed using sodium citrate (pH 6.0; 110 °C; 5–12 min overall) in a microwave (HistoS5, Milestone, Sorisole, Italy). Once retrieved, all tissue sections were stained in relevant batches in an autostainer (Autostainer 360, Epredia, Runcorn, UK) following a standard protocol that incorporates a 30‐min incubation with the primary antibody at room temperature. Anti‐CD3 was used to identify T cells and NK cells (1:200 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Catalogue No. NCL‐L‐CD3‐565), anti‐PAX5 was employed to identify B lymphocytes (1:50 Mouse Monoclonal Antibody Pax5, Becton and Dickinson, Workingham, UK, Catalogue No. P67320‐050), and anti‐IBA1 was used to label monocytes and macrophages (1:500 Rabbit polyclonal anti‐Iba1, FUJIFILM Wako, Richmond, VA, USA, Catalogue No. 019‐19741). All sections were treated with Dako Real Peroxidase Blocking solution (Dako S2023) to block endogenous peroxidases. The next step was a 40‐min incubation with either EnVision anti‐Mouse (Dako, EnVision+ System HRP‐labelled polymer; Catalogue No. K400111‐2) or EnVision anti‐Rabbit (Dako, EnVision+ System HRP‐labelled polymer; K400311‐2). This system is based on an HRP‐labelled polymer conjugated with secondary antibodies. The labelled polymer does not contain avidin or biotin. As such, nonspecific staining resulting from endogenous avidin‐biotin activity is eliminated or significantly reduced. Envision mouse was used to stain tissues for CD3 and PAX5, whereas EnVision Anti‐Rabbit with DAB enhancer protocol was used for IBA1 [Dako EnVision® + Dual Link System‐HRP (DAB+), Catalogue No. PD04048_02/K4065]. All sections were then counterstained with Harris' haematoxylin (Varistain Gemini Slide Stainer, Thermo Fisher Scientific), dehydrated with graded ethanol through to xylene, and coverslipped. Canine tissue was stained for CTLA4 and MET using the same immunohistochemical protocol as for human OMM. Image acquisition Slides were digitally scanned (NanoZoomer‐XR, Hamamatsu Photonics K.K, Shizuoka, Japan) and raw image data were saved in ndpi format and handled by the software ndp.view 2 (Hamamatsu Photonics) to save details of the whole image in JPEG format. Images were imported to ImageJ (FiJi Project, 2.11.0, GPLv3+) for labelling and uploaded to QuPath (version 0.5.0) . 16S sequencing DNA was isolated, purified, and quantified from FFPE core tissue samples from human and canine OMM as described earlier. DNA was similarly isolated, purified, and quantified from four FFPE core samples of normal canine oral tissue. Libraries were prepared from 2 μl of each purified DNA sample using a Quick‐16S™ Plus NGS Library Prep Kit (V3‐V4, UDI) (Zymo Research, Catalogue Nos. D6421‐PS1 and D6421‐PS2) following the manufacturer's protocol. Alongside the experimental samples, control libraries were prepared from the ZymoBIOMICS™ Microbial Community DNA Standard (positive), a ‘blank’ sample taken through the DNA extraction and purification processes without tissue (negative) and nuclease‐free water (nontemplate control). The source material for both the positive control and NTC were provided in the Quick‐16S™ Plus kit. The pool was purified using AMPure XP beads (Beckman Coulter, Catalogue No. A63881) to remove primer‐dimer sequences, and DNA concentration was assessed using fluorometry (Qubit dsDNA High Sensitivity Assay, Thermo Fisher Scientific, Catalogue No. Q32851) and assessed for quality and fragment size using an Agilent bioanalyser with the DNA High Sensitivity Kit (Agilent, Catalogue No. 5067‐4626). Fragment size and quantity measurements were used to calculate molarity for sequencing. Sequencing was performed using the NextSeq 2000 platform (Illumina, Catalogue No. SY‐415‐1002) using NextSeq 1000/2000 P1 Reagents (600 cycles) version 3 Kit (Catalogue No. 20075294). Loading concentration was 750 p m as recommended in the Quick‐16S Plus Kit User Guide. PhiX Control version 3 (Catalogue No. FC‐110‐3001) library was spiked in at 40% to help with cluster resolution and facilitate troubleshooting in case of any problems with the run. Expression data analysis pipeline Alignment and gene‐level counts RNA‐seq data were processed using the nf‐core ‘rnaseq’ pipeline version 3.8.1 . In brief, samples were aligned using STAR version 2.6.1d and gene‐based counts produced using Salmon version 1.5.2. Canine data were aligned to the ROS_Cfam_1.0 reference genome, annotated using the corresponding GTF file for build accession GCA_014441545.1. Human sequence data were aligned to hg38 and annotated with GENCODE version 41 . Unsupervised clustering Unsupervised consensus clustering of expression data was performed using the R/Bioconductor package cola . Prior to clustering, gene‐level count data were subjected to a variance‐stabilising transformation using DESeq2 , and the resulting matrix was used as input. Five different clustering algorithms were tested, with an evaluation of two to six clusters in each case. The cola algorithm resamples count data a fixed number of times, repeating the clustering process on each iteration. Upon completion, cola may present multiple stable solutions, with partitions varying in size and number. When determining which of these outputs to investigate, the best solution was selected based on the highest performance measured by the silhouette score, concordance, and proportion of ambiguous clustering (PAC) score. Differential expression analysis All downstream differential expression analysis was performed using DESeq2. Models were fitted treating the clusters identified by cola, sex, and, where appropriate, batch as factors. Log fold‐change (logFC) estimates were produced using the apeglm shrinkage method , which is intended to provide more robust estimates in the event of high within‐group variability. Shrunken logFC estimates were accompanied by s ‐values , an aggregate false sign rate, which are broadly analogous to q values. p values were also computed and adjusted using the independent hypothesis weighting (IHW) method . Results were annotated using biomaRt , and volcano plots were generated using EnhancedVolcano . Over‐representation analysis Over‐representation analysis (ORA) was performed using the R package ClusterProfiler . ORA was performed for up‐ and downregulated genes independently, taking the top differentially expressed genes with an s ‐value <0.01 in each case. These gene lists were then tested for over‐representation against the Biological Process (BP), Cellular Component (CC), and Molecular Function GO ontologies , in addition to known pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG) . Minimum and maximum gene set sizes were set to 10 and 1,000 respectively for all analyses. Immunocyte deconvolution The computational framework CIBERSORTx was used to deconvolute the tumour‐infiltrating immune cells from bulk RNA‐sequencing data . Twenty‐two immune cell subtypes were parsed from the annotated gene signature matrix LM22 and 100 permutations of the CIBERSORTx web portal . After running, only samples with CIBERSORT p values <0.05 were included in subsequent analyses. B‐mode batch correction was applied, and quantile normalisation was disabled to generate absolute scores. CIBERSORT fractions for naïve B cells, memory B cells, and plasma cells were aggregated into B cells; monocytes, M0, M1, M2 macrophages into ‘monocytes and macrophages’; and CD8, CD4, T follicular helper cells, regulatory T cells, and gamma delta T cells into T cells; and resting and activated NK cells into NK cells. Wilcoxon tests were used to compare the immune cell fractions between species, metastatic status, and transcriptomic consensus clustering groups. A Cox regression model was used to describe the time‐to‐event outcomes, with the adjustedCurves package  used to show the effect of a continuous variable adjusted for age and gender/sex. Receiver operating characteristic (ROC) curves were plotted, and optimal cut‐off points were defined using Youden's index. Statistical significance was set at p < 0.05. Image data analysis Using the HE slides to identify the tumour margins, the tumour areas were defined in QuPath . Positive cell detection was performed using the default parameters, except for sigma, which was set at 1.0 μm. Wilcoxon tests were used to compare the ‘number of positive cells detected per mm 2 ’ between transcriptomic consensus clustering groups. Optimal cut‐off points were defined as above. Statistical significance was set at p < 0.05. 16S analysis Cutadapt (version 4.5) was used to remove primers from paired‐end 16S rRNA gene reads . Both paired‐end reads were removed from further analyses if one or both of the pairs did not contain the 16S primers using the –discard‐untrimmed option. Mothur (version 1.48.0) was used for quality control, taxonomic assignment, and operational taxonomic unit (OTU) clustering , following an adjusted version of the mothur Miseq pipeline . Sequences were removed if they were >626 bp in length, <322 bp in length, contained ambiguous bases, had homopolymers of >9 bp in length, did not align to the correct variable region of the 16S rRNA gene, or did not originate from bacteria. Chimeras were identified via the ‘vsearch’ command, then removed. Sequence alignment and taxonomic assignment were conducted using the SILVA database  (version 138.1) formatted for use in mothur. OTUs were clustered using the ‘cluster.split’ command in mothur. The OTU file, phylogenetic tree file, taxonomy file, and sample metadata file were packaged into a phyloseq object . Data were decontaminated using decontam (version 1.22.0) , whereby contaminants were identified by comparing the prevalence (presence/absence across samples) of each sequence feature in the true positive samples with the prevalence in negative controls. The threshold of 0.5 was established to designate as contaminants exhibiting a higher prevalence in negative controls than positive controls. All contaminant sequence features were then removed from the entire dataset. Thirty‐four contaminants were removed from the dataset, resulting in 33,929 OTUs. Samples with <3,000 reads were discarded. Two human and one canine sample did not have >3,000 reads and were also discarded, resulting in 11 human OMM samples, 35 canine OMM samples, and four canine normal oral tissue samples. To adjust for sequencing depth across samples, raw counts for each OTU were converted into relative abundance values by dividing by the total amount of counts for each sample. Goods coverage was calculated using phyloseq‐coverage , and samples with a coverage value below 0.9 were considered to have insufficient sequencing depth and were discarded. To explore the diversity of microbial community composition within samples, we generated alpha‐diversity indices using both Shannon diversity and Chao1 richness. Differences in alpha diversity between cohorts were assessed by applying two‐sided Wilcoxon's rank‐sum tests. p values obtained from pairwise correlations were adjusted using the Benjamini–Hochberg procedure to control the false discovery rate. We presented significant findings as box‐and‐whisker plots. To evaluate whether samples clustered by group according to their overall microbiota compositions, we conducted a beta‐diversity analysis using permutational multivariate analysis of variance (PERMANOVA). This analysis was conducted with the ‘adonis2’ function and Bray–Curtis dissimilarity matrix, implemented in the vegan package in R (version 2.6‐4) . We presented significant findings as principal component analysis (PCA) plots. Statistical significance was defined as p < 0.05. To determine the OTUs most likely to explain differences between groups, we conducted linear discriminant analysis effect size (LefSE) analyses on variables exhibiting statistically significant alpha or beta diversity results. We used the R package microbiomeMarker (1.8.0)  to identify markers generated using all phylogenetic levels and presented results as tables and histographs.   This study was conducted with the approval of the University of Edinburgh's Veterinary Ethical Review Committee (Reference no. 10.21) and Lothian BioResource Access Committee (REC reference 20/ES/0061 SR1725). Formalin‐fixed, paraffin‐embedded (FFPE) human and canine OMM tissue samples and associated clinical data were collected from patients with naturally occurring tumours who had undergone surgical resection or biopsy of the primary tumour. Tissue was only included in the study with the informed, written consent of the human patient or the caregiver of the dog who bore the tumour. The treatment that a human or canine patient received was unaffected by their inclusion in this study. FFPE tissue samples of treatment‐naïve human OMMs were collected between 2006 and 2021. Clinicopathological data collected included patient gender, ethnicity, age (years), primary tumour location, presence of ulceration, tumour depth (mm), World Health Organization (WHO) stage, metastatic status, provision of adjuvant therapy, and duration of patient survival or lost to follow‐up (in months following initial diagnosis of OMM). FFPE tissue samples of treatment‐naïve canine OMMs were collected between 2011 and 2021 for histopathology from dogs attending the Clinical Veterinary Oncology departments at the Universities of Edinburgh, Bristol, Liverpool, Glasgow and Ohio Small Animal Teaching Hospitals, Vets Now (Glasgow) Referrals, Paragon Referrals, The Ralph Veterinary Referral Centre, Anderson Moores Veterinary Specialists, and North Downs Specialist Referrals. Clinicopathological data collected included sex, breed, age (years), primary tumour anatomical location, melanoma type (melanotic/amelanotic, pleomorphic/composite/nevocytoid/animal, ulcerated, exophytic), depth of invasion (mm), stroma:tumour ratio , mitotic count , metastatic status (confirmed by abdominal ultrasound, thoracic computed tomography or radiography, and/or histological examination of ≥1 regional lymph node), WHO stage, provision of adjuvant therapy, and duration of patient survival or lost to follow‐up (in months following initial diagnosis of OMM). The diagnosis of OMM was confirmed by a consultant histopathologist (MM) or board‐certified veterinary anatomic pathologist (JdP), and the regions of the FFPE blocks that represented OMM tissue tumour identified.   and DNA extraction Using a manual tissue arrayer (Beecher MTA‐1, Estigen OÜ, Tartu, Estonia) machine, two cores (1 mm outer diameter, 5 mm length) were harvested from each canine FFPE block: one from an area of OMM and one from histopathologically normal tissue. RNA and DNA were isolated from the tissue cores using a Covaris E220 Evolution Focused Ultrasonicator (Covaris Ltd, Woddington, Brighton, UK, Catalogue No. 500429) and the truXTRAC® FFPE total NA (Nucleic Acid) Kit – Column (Covaris, Catalogue No. 520220). Inhibitors (e.g. melanin) were removed using the OneStep PCR Inhibitor Removal kit (Zymo Research, Irvine, CA, USA, Catalogue No. D6030). Quality control DNA (canine only) and total RNA (human and canine) were assessed for quality and integrity using a Fragment Analyser Automated Capillary Electrophoresis System (Agilent Technologies, Santa Clara, CA, USA, Catalogue No. 5300) with the Genomic DNA 50 kb Kit (Agilent Technologies, Catalogue No. DNF‐467‐0500) and Standard Sensitivity RNA Analysis Kit (Agilent Technologies, Catalogue No. DNF‐471‐0500) respectively, and then assessed using a Qubit 2.0 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA, Catalogue No. Q32866) and the Qubit DNA (Thermo Fisher Scientific, Catalogue No. Q32853) or RNA (Thermo Fisher Scientific, Catalogue No. Q10210) broad range assay kits. In the RNA samples, DNA contamination was quantified using a Qubit dsDNA HS assay kit (Thermo Fisher Scientific, Catalogue No. Q32854). The Fragment Analyser data showed a high degree of RNA degradation, so the fragmentation step was omitted from the library preparation protocol. RNA library preparation First‐strand cDNA was generated from 50 ng of each total RNA sample using a SMARTer® Stranded Total RNA‐Seq Kit version 2 – Pico Input Mammalian kit (Clontech Laboratories, Mountain View, CA, USA, Catalogue No. 634411). Illumina‐compatible adapters and indexes were then added via five cycles of PCR. AMPure XP beads (Beckman Coulter, Brea, CA, USA, Catalogue No. A63881) were then used to purify the cDNA library. Depletion of ribosomal cDNA (cDNA fragments originating from highly abundant rRNA molecules) was performed using ZapR version 2 and R‐probes version 2 specific to mammalian ribosomal RNA and human mitochondrial rRNA. Uncleaved fragments were then enriched by 13 cycles of PCR before a final library purification using AMPure XP beads. Libraries were quantified with the Qubit dsDNA HS assay and assessed for quality and size distribution of library fragments using the Fragment Analyser and a NGS Fragment Kit (Agilent Technologies, Catalogue No. DNF‐473‐0500). DNA library preparation Libraries were prepared from up to 300 ng of each DNA sample using the Sure Select XT Target Enrichment for Illumina System according to the provided protocol for FFPE samples. DNA samples were sheared to an average fragment size of 150–200 bp using a Covaris E220 Evolution and purified using Agencourt AMPure XP beads (Beckman Coulter, Catalogue No. A6881). Purified DNA fragments were end‐repaired to remove 3’‐ and 5’ overhangs before further purification with AMPure XP beads. A single ‘A' nucleotide was added to the 3’ ends of the blunt fragments to prevent them from ligating to another during the subsequent adapter ligation reaction, and a corresponding single ‘T' nucleotide on the 3’ end of the adapter provided a complementary overhang for ligating the adapter to the fragment. dA‐tailed DNA was again purified using AMPure XP beads before paired‐end adapters were ligated to the ends of the dA‐tailed DNA fragments to prepare them for hybridisation in a flow cell. After purification using AMPure XP beads, 13 cycles of PCR were used to selectively enrich those DNA fragments that had adapter molecules on both ends and amplify the amount of DNA in the library suitable for target enrichment. Amplified gDNA libraries were purified using AMPure XP beads. gDNA libraries were quantified by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser using a DNA HS Kit (Agilent Technologies, Catalogue No. 5067‐4626). Seven hundred fifty ng of each gDNA library (where available) was hybridised to target‐specific probes (Canine All‐Exon) before target molecules were captured on streptavidin beads and non‐target molecules removed with a series of washes. Individual hybridised and captured libraries were then amplified using 12 cycles of PCR (with unique indexing primers to allow multiplexed sequencing) and purified with AMPure XP beads. Indexed libraries were measured by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser with a DNA HS Kit. DNA concentration and fragment size data were used to calculate molarity for sequencing. Sequencing Sequencing was performed using the NextSeq 2000 platform (Illumina, San Diego, CA, USA, Catalogue No. 20038897) with NextSeq 2000 P3 Reagents (200 Cycles) (Illumina, Catalogue No. 20040560). Libraries were combined in equimolar pools based on Qubit and bioanalyser assay results, and each pool was sequenced in a P3 flow cell. For RNA sequencing (RNA‐seq), PhiX Control version 3 (Illumina, Catalogue No. FC‐110‐3001) was spiked into each run at a concentration of 1% to allow troubleshooting in the event of any issues. Base call data produced by the NextSeq 1000/2000 Control Software (version 1.4.1.39716) were automatically converted into FASTQ files and uploaded to BaseSpace.   DNA (canine only) and total RNA (human and canine) were assessed for quality and integrity using a Fragment Analyser Automated Capillary Electrophoresis System (Agilent Technologies, Santa Clara, CA, USA, Catalogue No. 5300) with the Genomic DNA 50 kb Kit (Agilent Technologies, Catalogue No. DNF‐467‐0500) and Standard Sensitivity RNA Analysis Kit (Agilent Technologies, Catalogue No. DNF‐471‐0500) respectively, and then assessed using a Qubit 2.0 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA, Catalogue No. Q32866) and the Qubit DNA (Thermo Fisher Scientific, Catalogue No. Q32853) or RNA (Thermo Fisher Scientific, Catalogue No. Q10210) broad range assay kits. In the RNA samples, DNA contamination was quantified using a Qubit dsDNA HS assay kit (Thermo Fisher Scientific, Catalogue No. Q32854). The Fragment Analyser data showed a high degree of RNA degradation, so the fragmentation step was omitted from the library preparation protocol.   First‐strand cDNA was generated from 50 ng of each total RNA sample using a SMARTer® Stranded Total RNA‐Seq Kit version 2 – Pico Input Mammalian kit (Clontech Laboratories, Mountain View, CA, USA, Catalogue No. 634411). Illumina‐compatible adapters and indexes were then added via five cycles of PCR. AMPure XP beads (Beckman Coulter, Brea, CA, USA, Catalogue No. A63881) were then used to purify the cDNA library. Depletion of ribosomal cDNA (cDNA fragments originating from highly abundant rRNA molecules) was performed using ZapR version 2 and R‐probes version 2 specific to mammalian ribosomal RNA and human mitochondrial rRNA. Uncleaved fragments were then enriched by 13 cycles of PCR before a final library purification using AMPure XP beads. Libraries were quantified with the Qubit dsDNA HS assay and assessed for quality and size distribution of library fragments using the Fragment Analyser and a NGS Fragment Kit (Agilent Technologies, Catalogue No. DNF‐473‐0500).   Libraries were prepared from up to 300 ng of each DNA sample using the Sure Select XT Target Enrichment for Illumina System according to the provided protocol for FFPE samples. DNA samples were sheared to an average fragment size of 150–200 bp using a Covaris E220 Evolution and purified using Agencourt AMPure XP beads (Beckman Coulter, Catalogue No. A6881). Purified DNA fragments were end‐repaired to remove 3’‐ and 5’ overhangs before further purification with AMPure XP beads. A single ‘A' nucleotide was added to the 3’ ends of the blunt fragments to prevent them from ligating to another during the subsequent adapter ligation reaction, and a corresponding single ‘T' nucleotide on the 3’ end of the adapter provided a complementary overhang for ligating the adapter to the fragment. dA‐tailed DNA was again purified using AMPure XP beads before paired‐end adapters were ligated to the ends of the dA‐tailed DNA fragments to prepare them for hybridisation in a flow cell. After purification using AMPure XP beads, 13 cycles of PCR were used to selectively enrich those DNA fragments that had adapter molecules on both ends and amplify the amount of DNA in the library suitable for target enrichment. Amplified gDNA libraries were purified using AMPure XP beads. gDNA libraries were quantified by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser using a DNA HS Kit (Agilent Technologies, Catalogue No. 5067‐4626). Seven hundred fifty ng of each gDNA library (where available) was hybridised to target‐specific probes (Canine All‐Exon) before target molecules were captured on streptavidin beads and non‐target molecules removed with a series of washes. Individual hybridised and captured libraries were then amplified using 12 cycles of PCR (with unique indexing primers to allow multiplexed sequencing) and purified with AMPure XP beads. Indexed libraries were measured by fluorometry using a Qubit dsDNA HS assay and assessed for quality and fragment size using the Agilent bioanalyser with a DNA HS Kit. DNA concentration and fragment size data were used to calculate molarity for sequencing.   Sequencing was performed using the NextSeq 2000 platform (Illumina, San Diego, CA, USA, Catalogue No. 20038897) with NextSeq 2000 P3 Reagents (200 Cycles) (Illumina, Catalogue No. 20040560). Libraries were combined in equimolar pools based on Qubit and bioanalyser assay results, and each pool was sequenced in a P3 flow cell. For RNA sequencing (RNA‐seq), PhiX Control version 3 (Illumina, Catalogue No. FC‐110‐3001) was spiked into each run at a concentration of 1% to allow troubleshooting in the event of any issues. Base call data produced by the NextSeq 1000/2000 Control Software (version 1.4.1.39716) were automatically converted into FASTQ files and uploaded to BaseSpace.   Human OMM Immunohistochemistry was conducted on 3‐μm‐thick FFPE human OMM tissue using a Leica Bond immunostainer (Leica Bond III Autostainer, Leica, Sheffield, UK, Catalogue No. 21.2201). A standard protocol was followed in accordance with the manufacturer's instructions. Anti‐CD3 was used to identify T cells and natural killer (NK) cells (1:500 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Sheffield, UK, Catalogue No. NCL‐L‐CD3‐565), Anti‐CD19 was employed to stain B‐lymphocytes (1:50 Novocastra Liquid Mouse Monoclonal Antibody CD19, Leica Biosystems, Catalogue No. NCL‐L‐CD19‐163), Anti‐CD68 was used to label monocytes and macrophages (1:50 Monoclonal Mouse Anti‐Human CD68, Dako, Santa Clara, CA, USA, Catalogue No. M0876). Tissue was also stained for CD163 (1:50, CD163 Monoclonal Antibody, Thermo Fisher Scientific, Catalogue No. MA5‐11458), CTLA4 (1:50 CTLA4 Polyclonal Antibody, Invitrogen, Darmstadt, Germany, Catalogue No. PA5‐115060), and MET (1:400 HGFR/c‐MET Antibody, 1G7NB, Novus Biologicals, Abingdon, Oxford, UK, Catalogue No. 44306SS) using the same protocol. Canine OMM Immunohistochemistry was conducted on 4‐μm‐thick FFPE canine OMM tissue, cut using a microtome (Thermo Fisher Scientific Rotary Microtome Microm HM 340 E). All tissue sections, including relevant controls, were deparaffinised in xylene and rehydrated in graded ethanol through to distilled H 2 O prior to staining. Heat‐induced antigen retrieval was performed using sodium citrate (pH 6.0; 110 °C; 5–12 min overall) in a microwave (HistoS5, Milestone, Sorisole, Italy). Once retrieved, all tissue sections were stained in relevant batches in an autostainer (Autostainer 360, Epredia, Runcorn, UK) following a standard protocol that incorporates a 30‐min incubation with the primary antibody at room temperature. Anti‐CD3 was used to identify T cells and NK cells (1:200 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Catalogue No. NCL‐L‐CD3‐565), anti‐PAX5 was employed to identify B lymphocytes (1:50 Mouse Monoclonal Antibody Pax5, Becton and Dickinson, Workingham, UK, Catalogue No. P67320‐050), and anti‐IBA1 was used to label monocytes and macrophages (1:500 Rabbit polyclonal anti‐Iba1, FUJIFILM Wako, Richmond, VA, USA, Catalogue No. 019‐19741). All sections were treated with Dako Real Peroxidase Blocking solution (Dako S2023) to block endogenous peroxidases. The next step was a 40‐min incubation with either EnVision anti‐Mouse (Dako, EnVision+ System HRP‐labelled polymer; Catalogue No. K400111‐2) or EnVision anti‐Rabbit (Dako, EnVision+ System HRP‐labelled polymer; K400311‐2). This system is based on an HRP‐labelled polymer conjugated with secondary antibodies. The labelled polymer does not contain avidin or biotin. As such, nonspecific staining resulting from endogenous avidin‐biotin activity is eliminated or significantly reduced. Envision mouse was used to stain tissues for CD3 and PAX5, whereas EnVision Anti‐Rabbit with DAB enhancer protocol was used for IBA1 [Dako EnVision® + Dual Link System‐HRP (DAB+), Catalogue No. PD04048_02/K4065]. All sections were then counterstained with Harris' haematoxylin (Varistain Gemini Slide Stainer, Thermo Fisher Scientific), dehydrated with graded ethanol through to xylene, and coverslipped. Canine tissue was stained for CTLA4 and MET using the same immunohistochemical protocol as for human OMM.   Immunohistochemistry was conducted on 3‐μm‐thick FFPE human OMM tissue using a Leica Bond immunostainer (Leica Bond III Autostainer, Leica, Sheffield, UK, Catalogue No. 21.2201). A standard protocol was followed in accordance with the manufacturer's instructions. Anti‐CD3 was used to identify T cells and natural killer (NK) cells (1:500 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Sheffield, UK, Catalogue No. NCL‐L‐CD3‐565), Anti‐CD19 was employed to stain B‐lymphocytes (1:50 Novocastra Liquid Mouse Monoclonal Antibody CD19, Leica Biosystems, Catalogue No. NCL‐L‐CD19‐163), Anti‐CD68 was used to label monocytes and macrophages (1:50 Monoclonal Mouse Anti‐Human CD68, Dako, Santa Clara, CA, USA, Catalogue No. M0876). Tissue was also stained for CD163 (1:50, CD163 Monoclonal Antibody, Thermo Fisher Scientific, Catalogue No. MA5‐11458), CTLA4 (1:50 CTLA4 Polyclonal Antibody, Invitrogen, Darmstadt, Germany, Catalogue No. PA5‐115060), and MET (1:400 HGFR/c‐MET Antibody, 1G7NB, Novus Biologicals, Abingdon, Oxford, UK, Catalogue No. 44306SS) using the same protocol.   Immunohistochemistry was conducted on 4‐μm‐thick FFPE canine OMM tissue, cut using a microtome (Thermo Fisher Scientific Rotary Microtome Microm HM 340 E). All tissue sections, including relevant controls, were deparaffinised in xylene and rehydrated in graded ethanol through to distilled H 2 O prior to staining. Heat‐induced antigen retrieval was performed using sodium citrate (pH 6.0; 110 °C; 5–12 min overall) in a microwave (HistoS5, Milestone, Sorisole, Italy). Once retrieved, all tissue sections were stained in relevant batches in an autostainer (Autostainer 360, Epredia, Runcorn, UK) following a standard protocol that incorporates a 30‐min incubation with the primary antibody at room temperature. Anti‐CD3 was used to identify T cells and NK cells (1:200 Novocastra Liquid Mouse Monoclonal Antibody CD3, Leica Biosystems, Catalogue No. NCL‐L‐CD3‐565), anti‐PAX5 was employed to identify B lymphocytes (1:50 Mouse Monoclonal Antibody Pax5, Becton and Dickinson, Workingham, UK, Catalogue No. P67320‐050), and anti‐IBA1 was used to label monocytes and macrophages (1:500 Rabbit polyclonal anti‐Iba1, FUJIFILM Wako, Richmond, VA, USA, Catalogue No. 019‐19741). All sections were treated with Dako Real Peroxidase Blocking solution (Dako S2023) to block endogenous peroxidases. The next step was a 40‐min incubation with either EnVision anti‐Mouse (Dako, EnVision+ System HRP‐labelled polymer; Catalogue No. K400111‐2) or EnVision anti‐Rabbit (Dako, EnVision+ System HRP‐labelled polymer; K400311‐2). This system is based on an HRP‐labelled polymer conjugated with secondary antibodies. The labelled polymer does not contain avidin or biotin. As such, nonspecific staining resulting from endogenous avidin‐biotin activity is eliminated or significantly reduced. Envision mouse was used to stain tissues for CD3 and PAX5, whereas EnVision Anti‐Rabbit with DAB enhancer protocol was used for IBA1 [Dako EnVision® + Dual Link System‐HRP (DAB+), Catalogue No. PD04048_02/K4065]. All sections were then counterstained with Harris' haematoxylin (Varistain Gemini Slide Stainer, Thermo Fisher Scientific), dehydrated with graded ethanol through to xylene, and coverslipped. Canine tissue was stained for CTLA4 and MET using the same immunohistochemical protocol as for human OMM.   Slides were digitally scanned (NanoZoomer‐XR, Hamamatsu Photonics K.K, Shizuoka, Japan) and raw image data were saved in ndpi format and handled by the software ndp.view 2 (Hamamatsu Photonics) to save details of the whole image in JPEG format. Images were imported to ImageJ (FiJi Project, 2.11.0, GPLv3+) for labelling and uploaded to QuPath (version 0.5.0) .   sequencing DNA was isolated, purified, and quantified from FFPE core tissue samples from human and canine OMM as described earlier. DNA was similarly isolated, purified, and quantified from four FFPE core samples of normal canine oral tissue. Libraries were prepared from 2 μl of each purified DNA sample using a Quick‐16S™ Plus NGS Library Prep Kit (V3‐V4, UDI) (Zymo Research, Catalogue Nos. D6421‐PS1 and D6421‐PS2) following the manufacturer's protocol. Alongside the experimental samples, control libraries were prepared from the ZymoBIOMICS™ Microbial Community DNA Standard (positive), a ‘blank’ sample taken through the DNA extraction and purification processes without tissue (negative) and nuclease‐free water (nontemplate control). The source material for both the positive control and NTC were provided in the Quick‐16S™ Plus kit. The pool was purified using AMPure XP beads (Beckman Coulter, Catalogue No. A63881) to remove primer‐dimer sequences, and DNA concentration was assessed using fluorometry (Qubit dsDNA High Sensitivity Assay, Thermo Fisher Scientific, Catalogue No. Q32851) and assessed for quality and fragment size using an Agilent bioanalyser with the DNA High Sensitivity Kit (Agilent, Catalogue No. 5067‐4626). Fragment size and quantity measurements were used to calculate molarity for sequencing. Sequencing was performed using the NextSeq 2000 platform (Illumina, Catalogue No. SY‐415‐1002) using NextSeq 1000/2000 P1 Reagents (600 cycles) version 3 Kit (Catalogue No. 20075294). Loading concentration was 750 p m as recommended in the Quick‐16S Plus Kit User Guide. PhiX Control version 3 (Catalogue No. FC‐110‐3001) library was spiked in at 40% to help with cluster resolution and facilitate troubleshooting in case of any problems with the run.   Alignment and gene‐level counts RNA‐seq data were processed using the nf‐core ‘rnaseq’ pipeline version 3.8.1 . In brief, samples were aligned using STAR version 2.6.1d and gene‐based counts produced using Salmon version 1.5.2. Canine data were aligned to the ROS_Cfam_1.0 reference genome, annotated using the corresponding GTF file for build accession GCA_014441545.1. Human sequence data were aligned to hg38 and annotated with GENCODE version 41 . Unsupervised clustering Unsupervised consensus clustering of expression data was performed using the R/Bioconductor package cola . Prior to clustering, gene‐level count data were subjected to a variance‐stabilising transformation using DESeq2 , and the resulting matrix was used as input. Five different clustering algorithms were tested, with an evaluation of two to six clusters in each case. The cola algorithm resamples count data a fixed number of times, repeating the clustering process on each iteration. Upon completion, cola may present multiple stable solutions, with partitions varying in size and number. When determining which of these outputs to investigate, the best solution was selected based on the highest performance measured by the silhouette score, concordance, and proportion of ambiguous clustering (PAC) score. Differential expression analysis All downstream differential expression analysis was performed using DESeq2. Models were fitted treating the clusters identified by cola, sex, and, where appropriate, batch as factors. Log fold‐change (logFC) estimates were produced using the apeglm shrinkage method , which is intended to provide more robust estimates in the event of high within‐group variability. Shrunken logFC estimates were accompanied by s ‐values , an aggregate false sign rate, which are broadly analogous to q values. p values were also computed and adjusted using the independent hypothesis weighting (IHW) method . Results were annotated using biomaRt , and volcano plots were generated using EnhancedVolcano . Over‐representation analysis Over‐representation analysis (ORA) was performed using the R package ClusterProfiler . ORA was performed for up‐ and downregulated genes independently, taking the top differentially expressed genes with an s ‐value <0.01 in each case. These gene lists were then tested for over‐representation against the Biological Process (BP), Cellular Component (CC), and Molecular Function GO ontologies , in addition to known pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG) . Minimum and maximum gene set sizes were set to 10 and 1,000 respectively for all analyses.   RNA‐seq data were processed using the nf‐core ‘rnaseq’ pipeline version 3.8.1 . In brief, samples were aligned using STAR version 2.6.1d and gene‐based counts produced using Salmon version 1.5.2. Canine data were aligned to the ROS_Cfam_1.0 reference genome, annotated using the corresponding GTF file for build accession GCA_014441545.1. Human sequence data were aligned to hg38 and annotated with GENCODE version 41 .   Unsupervised consensus clustering of expression data was performed using the R/Bioconductor package cola . Prior to clustering, gene‐level count data were subjected to a variance‐stabilising transformation using DESeq2 , and the resulting matrix was used as input. Five different clustering algorithms were tested, with an evaluation of two to six clusters in each case. The cola algorithm resamples count data a fixed number of times, repeating the clustering process on each iteration. Upon completion, cola may present multiple stable solutions, with partitions varying in size and number. When determining which of these outputs to investigate, the best solution was selected based on the highest performance measured by the silhouette score, concordance, and proportion of ambiguous clustering (PAC) score.   All downstream differential expression analysis was performed using DESeq2. Models were fitted treating the clusters identified by cola, sex, and, where appropriate, batch as factors. Log fold‐change (logFC) estimates were produced using the apeglm shrinkage method , which is intended to provide more robust estimates in the event of high within‐group variability. Shrunken logFC estimates were accompanied by s ‐values , an aggregate false sign rate, which are broadly analogous to q values. p values were also computed and adjusted using the independent hypothesis weighting (IHW) method . Results were annotated using biomaRt , and volcano plots were generated using EnhancedVolcano .   Over‐representation analysis (ORA) was performed using the R package ClusterProfiler . ORA was performed for up‐ and downregulated genes independently, taking the top differentially expressed genes with an s ‐value <0.01 in each case. These gene lists were then tested for over‐representation against the Biological Process (BP), Cellular Component (CC), and Molecular Function GO ontologies , in addition to known pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG) . Minimum and maximum gene set sizes were set to 10 and 1,000 respectively for all analyses.   The computational framework CIBERSORTx was used to deconvolute the tumour‐infiltrating immune cells from bulk RNA‐sequencing data . Twenty‐two immune cell subtypes were parsed from the annotated gene signature matrix LM22 and 100 permutations of the CIBERSORTx web portal . After running, only samples with CIBERSORT p values <0.05 were included in subsequent analyses. B‐mode batch correction was applied, and quantile normalisation was disabled to generate absolute scores. CIBERSORT fractions for naïve B cells, memory B cells, and plasma cells were aggregated into B cells; monocytes, M0, M1, M2 macrophages into ‘monocytes and macrophages’; and CD8, CD4, T follicular helper cells, regulatory T cells, and gamma delta T cells into T cells; and resting and activated NK cells into NK cells. Wilcoxon tests were used to compare the immune cell fractions between species, metastatic status, and transcriptomic consensus clustering groups. A Cox regression model was used to describe the time‐to‐event outcomes, with the adjustedCurves package  used to show the effect of a continuous variable adjusted for age and gender/sex. Receiver operating characteristic (ROC) curves were plotted, and optimal cut‐off points were defined using Youden's index. Statistical significance was set at p < 0.05.   Using the HE slides to identify the tumour margins, the tumour areas were defined in QuPath . Positive cell detection was performed using the default parameters, except for sigma, which was set at 1.0 μm. Wilcoxon tests were used to compare the ‘number of positive cells detected per mm 2 ’ between transcriptomic consensus clustering groups. Optimal cut‐off points were defined as above. Statistical significance was set at p < 0.05.   analysis Cutadapt (version 4.5) was used to remove primers from paired‐end 16S rRNA gene reads . Both paired‐end reads were removed from further analyses if one or both of the pairs did not contain the 16S primers using the –discard‐untrimmed option. Mothur (version 1.48.0) was used for quality control, taxonomic assignment, and operational taxonomic unit (OTU) clustering , following an adjusted version of the mothur Miseq pipeline . Sequences were removed if they were >626 bp in length, <322 bp in length, contained ambiguous bases, had homopolymers of >9 bp in length, did not align to the correct variable region of the 16S rRNA gene, or did not originate from bacteria. Chimeras were identified via the ‘vsearch’ command, then removed. Sequence alignment and taxonomic assignment were conducted using the SILVA database  (version 138.1) formatted for use in mothur. OTUs were clustered using the ‘cluster.split’ command in mothur. The OTU file, phylogenetic tree file, taxonomy file, and sample metadata file were packaged into a phyloseq object . Data were decontaminated using decontam (version 1.22.0) , whereby contaminants were identified by comparing the prevalence (presence/absence across samples) of each sequence feature in the true positive samples with the prevalence in negative controls. The threshold of 0.5 was established to designate as contaminants exhibiting a higher prevalence in negative controls than positive controls. All contaminant sequence features were then removed from the entire dataset. Thirty‐four contaminants were removed from the dataset, resulting in 33,929 OTUs. Samples with <3,000 reads were discarded. Two human and one canine sample did not have >3,000 reads and were also discarded, resulting in 11 human OMM samples, 35 canine OMM samples, and four canine normal oral tissue samples. To adjust for sequencing depth across samples, raw counts for each OTU were converted into relative abundance values by dividing by the total amount of counts for each sample. Goods coverage was calculated using phyloseq‐coverage , and samples with a coverage value below 0.9 were considered to have insufficient sequencing depth and were discarded. To explore the diversity of microbial community composition within samples, we generated alpha‐diversity indices using both Shannon diversity and Chao1 richness. Differences in alpha diversity between cohorts were assessed by applying two‐sided Wilcoxon's rank‐sum tests. p values obtained from pairwise correlations were adjusted using the Benjamini–Hochberg procedure to control the false discovery rate. We presented significant findings as box‐and‐whisker plots. To evaluate whether samples clustered by group according to their overall microbiota compositions, we conducted a beta‐diversity analysis using permutational multivariate analysis of variance (PERMANOVA). This analysis was conducted with the ‘adonis2’ function and Bray–Curtis dissimilarity matrix, implemented in the vegan package in R (version 2.6‐4) . We presented significant findings as principal component analysis (PCA) plots. Statistical significance was defined as p < 0.05. To determine the OTUs most likely to explain differences between groups, we conducted linear discriminant analysis effect size (LefSE) analyses on variables exhibiting statistically significant alpha or beta diversity results. We used the R package microbiomeMarker (1.8.0)  to identify markers generated using all phylogenetic levels and presented results as tables and histographs.   OMM in humans and dogs share clinicopathological characteristics and poor prognosis We identified human (Edinburgh cohort; n = 17) and canine (Bowlt Blacklock cohort; n = 36) patients who presented for treatment of naturally occurring, treatment‐naive OMM and had archival FFPE tissue and clinical records available for analysis (Figure ). In both humans and dogs, OMM can appear as melanotic or amelanotic and is most frequently observed in the oronasal cavity in humans and in the oral cavity in dogs (Figure ). In both species, underlying bone lysis secondary to local tumour invasion is common (Figure ). Histopathological findings are also similar, consisting of proliferations of melanotic or amelanotic infiltrates of pleomorphic, round, and/or spindloid neoplastic cells that invade both the lamina propria and deep connective tissue (Figure ). RNA‐seq and 16S sequencing data were compared to clinical data. When considering survival, the time from OMM diagnosis to end of life is not directly comparable between species. Human patients die of their disease, whereas in dogs, 89.3% of recorded deaths are a result of euthanasia . We investigated the relationships between clinicopathological features, treatment, and survival. Clinicopathological data are presented in Figure . The best model to identify predictors of survival in humans included age, gender (male), presence of metastasis, and immunocyte infiltration (Figure ). Risk of death increased by 14% for every 1‐year increase in age, 16.5‐fold for male patients, and 40‐fold for patients with metastatic disease (Figure ). In dogs, the best model to identify predictors of survival included only age and presence of metastasis (Figure ). Human patients who had undergone only biopsy of their OMM without surgery had a significantly increased risk of death compared to those who received curative‐intent surgery (supplementary material, Figure ). We found no significant correlation between the relative hazard of death and provision of chemotherapy, radiation therapy, or immunotherapy in humans or dogs (supplementary material, Figure ). In our canine patients, reduced survival probability was identified in patients with WHO Stage 3 or 4 disease and ulcerated tumours (supplementary material, Figure ). Exophytic status, melanotic status, and melanoma histopathological subtype were not associated with survival in dogs (supplementary material, Figure ). Equivalent data were not available for human patients. In conclusion, for both human and canine patients with OMM, our clinical data reveal that prognosis is poor and metastasis common. In humans, our data show that curative surgery and being female are associated with improved survival, while radiation, immunotherapy, and chemotherapy do not significantly impact survival in either species. OMM can be stratified into two transcriptomic subgroups Despite no widespread actionable mutation targets in OMM, we hypothesised that OMM might have transcriptional subtypes, as found in CM . We conducted RNA‐seq transcriptomic analysis of OMM from human (Edinburgh) and canine (Bowlt Blacklock) patients and performed unsupervised hierarchical clustering from these cohorts as well as from a published canine OMM RNA dataset . Critically, we identified two stable consensus clusters (A and B) in both species using optimal consensus partitioning with median absolute deviation and spherical k‐means clustering (Figure ; see ). PCA showed minimal overlap between these transcriptional subgroups in all three datasets (supplementary material, Figure ). To understand the pathways that characterised the transcriptional subtypes, we applied ORA and showed that KEGG pathways and biological processes aligned with an immune response in consensus cluster group B (Figure ). In contrast, ORA of consensus cluster group A was less clear but highlighted pathways associated with the cell cycle and DNA/RNA processing and repair (supplementary material, Data ). To investigate whether the transcriptional subgroups shared between human and canine OMM cohorts were mediated through overlapping gene networks, we aimed to find the most important variables that drove the stratification across species and cohorts. We used the machine learning algorithm randomForest  and 812 human and canine homologous genes (supplementary material, Data ). Across five randomly selected seeds, this approach generated successful models (trained on 50% of the samples) that predicted the B subgroup versus the A subgroup in the test data (the remaining 50% of the samples), independently of species (accuracy 79–88%, binomial test, p < 0.001). Within these models, 41 genes were consistently ranked among the top 100 genes shared by at least four models, and 17 genes were shared by all five models (supplementary material, Figure and Data ). Using these 41 genes, optimal discrimination between transcriptomic subgroups was achieved by Euclidean distance clustering (Figure ). Using PCA clustering based on these 41 genes, we found that our human OMM and dog OMM RNA‐seq data closely clustered on PCA, indicating strong similarity (Figure ). The Prouteau dataset demonstrated a mildly distinct distribution on PC2 (9.19% variation) but not PC1, with the latter strikingly driven by the transcriptional subgroup variation (42.78%) (Figure ). Therefore, we concluded that the PC2 difference was mainly caused by technical variations during cohort establishment (Edinburgh and Bowlt Blacklock versus Prouteau), while biologically the two OMM transcriptional subtypes did not exhibit any species‐dependent variations. Notably, our signature highlighted CTLA‐4 gene expression in both human and dog OMM in subgroup B. This finding was validated by immunohistochemistry (IHC), showing a significantly greater number of CTLA4‐positive cells in the B transcriptomic subgroup (human: p = 0.01; canine: p = 0.03) (Figure ). We therefore termed this subgroup ‘CTLA4‐high’. Of particular interest in subgroup A is a small, enriched panel of genes ( GJC1, CLCN3, TYR, MET, PLEKHAB, FOSB, KCNA1 ) highlighted in the 41‐gene signature. Using g:Profiler, we mapped these genes to functional information sources to detect statistically significantly enriched terms, which included the hepatocyte growth factor receptor (HGFR, also known as MET) activity pathway and the MET complex. We performed MET IHC and found MET protein levels to be significantly increased in subgroup A in both humans and dogs (human: p = 0.04; canine: p = 0.02) (Figure ). Therefore, we termed the A subgroup ‘MET‐high’. Additionally, we noticed that the enriched panel of genes were microphthalmia‐associated transcription factor (MITF) target genes and suggest the gene expression signatures of this subgroup may be mediated by MITF  (supplementary material, Figure ). In conclusion, we show that human and dog OMM can be stratified into two distinctive transcriptional groups, termed CTLA4‐high and MET‐high, using the same species‐independent 41‐gene signature, providing actionable targets and highlighting the biological similarity of this disease between the two species. CTLA4 and macrophages are diagnostic biomarkers for OMM transcriptomic subtypes The presence of potential drug targets CTLA4 and MET within the subtypes indicated that the stratification of patients into each subtype could inform therapeutic decision‐making and clinical trial design. Because transcriptomic analysis is currently not practical in most clinical settings, we attempted to identify diagnostic IHC markers to discriminate between subtypes. We reasoned that the CTLA4‐high subtype could be identified by immunocyte infiltrate populations and utilised the computational framework CIBERSORTx to discern the composition of tumour‐infiltrating immune cells from bulk RNA‐seq data  (supplementary material, Data ). In both humans and dogs, total immune cell infiltration (the ‘Absolute Score’) was significantly larger in the CTLA4‐high transcriptomic subgroup ( p < 0.01 in both species), as well as the total B cells, monocyte/macrophages, and T cells (supplementary material, Figure ). Next, to extrapolate these data to a diagnostic test, we used IHC and digital pathology image analysis for immune infiltrate . Utilising species‐specific IHC for B cells, monocytes/macrophages, and T cells, we observed a higher density of B cells in dogs and higher densities of monocytes/macrophages and T cells in both species in tumour tissues of the CTLA4‐high transcriptomic cohort (Figure ). Then we utilised ROC analysis to determine how well these variables can discriminate between the two transcriptomic subgroups with clinically relevant sensitivity and specificity. Using the Youden index, we identified optimal cut‐off values for parameters with an area under the curve (AUC) of >0.75. ROC curve analysis demonstrated that monocytes and macrophages could differentiate between transcriptomic subgroups with high sensitivity in both humans and dogs and with high specificity in humans. This indicates that IHC using CD68 and IBA1 antibodies detecting monocyte and macrophages could be used to inform clinical decision‐making based on transcriptomic subtype. In addition, CTLA4 itself (but not MET) was associated with a clinically relevant sensitivity in both species and a clinically relevant specificity in humans (Figure and supplementary material, Figure ). Based on these findings and our RNA‐seq analysis, we recommend that the CTLA4‐high subtype OMMs be characterised as immune ‘hot’, with high levels of T‐cell infiltration and activation, and therefore provision of anti‐CTLA4 may show a favourable response in this subgroup . In contrast, our data suggest the MET‐high tumours are immune ‘cold’ and these patients may have a better response to targeted therapy. Tumour microbiota does not correlate with transcriptomic subgroup but may contribute to other pathological features of OMM The oral cavity hosts the largest and most diverse microbiota after the gut . Given the pivotal role of gastrointestinal tract microbiota in shaping both local and systemic immune system responses , we investigated whether the tumour microbiome correlated with our transcriptomic patient stratification, particularly the CTLA4‐high transcriptomic subgroup due to the high immune cell infiltration. To this end, we performed 16S sequencing of DNA extracted from FFPE cores harvested from our human and canine tumour cohorts as well as four normal canine oral tissue samples. Our analysis revealed a significant difference in alpha diversity between species, reflecting the microbiome's diversity within individual samples, with humans exhibiting significantly lower alpha diversity compared to dogs (Figure , p < 0.01). Beta diversity, a measure of the similarity or dissimilarity of groups of microbial communities, was significantly different between the two species ( p < 0.01; Figure and supplementary material, Data ). Due to the notable differences in alpha and beta diversity observed between species, we conducted separate analyses of the microbiomes associated with human and canine OMM. No significant difference was observed in the microbiome of human or canine OMM when stratified into transcriptomic subtype (Figure ). Consequently, the two transcriptomic subgroups of OMM were not attributed to an overall difference in the composition of the oral microbiota. Thus, the immunocyte infiltration that we identified above is likely attributed to inflammation with the primary tumour itself and not bacterial burden. Although tumour microbiota does not explain the OMM transcriptional subtypes, we did identify two clinically relevant associations that suggest that the tumour microbiome may have a function in OMM disease initiation and progression. First, in human OMM, the overall community diversity was significantly different between metastatic and non‐metastatic tumours (beta diversity p = 0.03; Figure ), with three OTU markers enriched in each group (Figure ). This finding suggests that the microbiota composition may influence the metastatic capacity of OMM . Second, we compared the microbiome profiles of canine OMM with those of normal canine oral tissue and found overall community diversity was significantly different between tumour and normal oral tissue (Figure ). These data indicate that OMM correlates with alterations in the oral microbiota composition, which warrants further investigation because of the emerging role of the microbiome in cancer development and response to therapy .   in humans and dogs share clinicopathological characteristics and poor prognosis We identified human (Edinburgh cohort; n = 17) and canine (Bowlt Blacklock cohort; n = 36) patients who presented for treatment of naturally occurring, treatment‐naive OMM and had archival FFPE tissue and clinical records available for analysis (Figure ). In both humans and dogs, OMM can appear as melanotic or amelanotic and is most frequently observed in the oronasal cavity in humans and in the oral cavity in dogs (Figure ). In both species, underlying bone lysis secondary to local tumour invasion is common (Figure ). Histopathological findings are also similar, consisting of proliferations of melanotic or amelanotic infiltrates of pleomorphic, round, and/or spindloid neoplastic cells that invade both the lamina propria and deep connective tissue (Figure ). RNA‐seq and 16S sequencing data were compared to clinical data. When considering survival, the time from OMM diagnosis to end of life is not directly comparable between species. Human patients die of their disease, whereas in dogs, 89.3% of recorded deaths are a result of euthanasia . We investigated the relationships between clinicopathological features, treatment, and survival. Clinicopathological data are presented in Figure . The best model to identify predictors of survival in humans included age, gender (male), presence of metastasis, and immunocyte infiltration (Figure ). Risk of death increased by 14% for every 1‐year increase in age, 16.5‐fold for male patients, and 40‐fold for patients with metastatic disease (Figure ). In dogs, the best model to identify predictors of survival included only age and presence of metastasis (Figure ). Human patients who had undergone only biopsy of their OMM without surgery had a significantly increased risk of death compared to those who received curative‐intent surgery (supplementary material, Figure ). We found no significant correlation between the relative hazard of death and provision of chemotherapy, radiation therapy, or immunotherapy in humans or dogs (supplementary material, Figure ). In our canine patients, reduced survival probability was identified in patients with WHO Stage 3 or 4 disease and ulcerated tumours (supplementary material, Figure ). Exophytic status, melanotic status, and melanoma histopathological subtype were not associated with survival in dogs (supplementary material, Figure ). Equivalent data were not available for human patients. In conclusion, for both human and canine patients with OMM, our clinical data reveal that prognosis is poor and metastasis common. In humans, our data show that curative surgery and being female are associated with improved survival, while radiation, immunotherapy, and chemotherapy do not significantly impact survival in either species.   can be stratified into two transcriptomic subgroups Despite no widespread actionable mutation targets in OMM, we hypothesised that OMM might have transcriptional subtypes, as found in CM . We conducted RNA‐seq transcriptomic analysis of OMM from human (Edinburgh) and canine (Bowlt Blacklock) patients and performed unsupervised hierarchical clustering from these cohorts as well as from a published canine OMM RNA dataset . Critically, we identified two stable consensus clusters (A and B) in both species using optimal consensus partitioning with median absolute deviation and spherical k‐means clustering (Figure ; see ). PCA showed minimal overlap between these transcriptional subgroups in all three datasets (supplementary material, Figure ). To understand the pathways that characterised the transcriptional subtypes, we applied ORA and showed that KEGG pathways and biological processes aligned with an immune response in consensus cluster group B (Figure ). In contrast, ORA of consensus cluster group A was less clear but highlighted pathways associated with the cell cycle and DNA/RNA processing and repair (supplementary material, Data ). To investigate whether the transcriptional subgroups shared between human and canine OMM cohorts were mediated through overlapping gene networks, we aimed to find the most important variables that drove the stratification across species and cohorts. We used the machine learning algorithm randomForest  and 812 human and canine homologous genes (supplementary material, Data ). Across five randomly selected seeds, this approach generated successful models (trained on 50% of the samples) that predicted the B subgroup versus the A subgroup in the test data (the remaining 50% of the samples), independently of species (accuracy 79–88%, binomial test, p < 0.001). Within these models, 41 genes were consistently ranked among the top 100 genes shared by at least four models, and 17 genes were shared by all five models (supplementary material, Figure and Data ). Using these 41 genes, optimal discrimination between transcriptomic subgroups was achieved by Euclidean distance clustering (Figure ). Using PCA clustering based on these 41 genes, we found that our human OMM and dog OMM RNA‐seq data closely clustered on PCA, indicating strong similarity (Figure ). The Prouteau dataset demonstrated a mildly distinct distribution on PC2 (9.19% variation) but not PC1, with the latter strikingly driven by the transcriptional subgroup variation (42.78%) (Figure ). Therefore, we concluded that the PC2 difference was mainly caused by technical variations during cohort establishment (Edinburgh and Bowlt Blacklock versus Prouteau), while biologically the two OMM transcriptional subtypes did not exhibit any species‐dependent variations. Notably, our signature highlighted CTLA‐4 gene expression in both human and dog OMM in subgroup B. This finding was validated by immunohistochemistry (IHC), showing a significantly greater number of CTLA4‐positive cells in the B transcriptomic subgroup (human: p = 0.01; canine: p = 0.03) (Figure ). We therefore termed this subgroup ‘CTLA4‐high’. Of particular interest in subgroup A is a small, enriched panel of genes ( GJC1, CLCN3, TYR, MET, PLEKHAB, FOSB, KCNA1 ) highlighted in the 41‐gene signature. Using g:Profiler, we mapped these genes to functional information sources to detect statistically significantly enriched terms, which included the hepatocyte growth factor receptor (HGFR, also known as MET) activity pathway and the MET complex. We performed MET IHC and found MET protein levels to be significantly increased in subgroup A in both humans and dogs (human: p = 0.04; canine: p = 0.02) (Figure ). Therefore, we termed the A subgroup ‘MET‐high’. Additionally, we noticed that the enriched panel of genes were microphthalmia‐associated transcription factor (MITF) target genes and suggest the gene expression signatures of this subgroup may be mediated by MITF  (supplementary material, Figure ). In conclusion, we show that human and dog OMM can be stratified into two distinctive transcriptional groups, termed CTLA4‐high and MET‐high, using the same species‐independent 41‐gene signature, providing actionable targets and highlighting the biological similarity of this disease between the two species.   and macrophages are diagnostic biomarkers for OMM transcriptomic subtypes The presence of potential drug targets CTLA4 and MET within the subtypes indicated that the stratification of patients into each subtype could inform therapeutic decision‐making and clinical trial design. Because transcriptomic analysis is currently not practical in most clinical settings, we attempted to identify diagnostic IHC markers to discriminate between subtypes. We reasoned that the CTLA4‐high subtype could be identified by immunocyte infiltrate populations and utilised the computational framework CIBERSORTx to discern the composition of tumour‐infiltrating immune cells from bulk RNA‐seq data  (supplementary material, Data ). In both humans and dogs, total immune cell infiltration (the ‘Absolute Score’) was significantly larger in the CTLA4‐high transcriptomic subgroup ( p < 0.01 in both species), as well as the total B cells, monocyte/macrophages, and T cells (supplementary material, Figure ). Next, to extrapolate these data to a diagnostic test, we used IHC and digital pathology image analysis for immune infiltrate . Utilising species‐specific IHC for B cells, monocytes/macrophages, and T cells, we observed a higher density of B cells in dogs and higher densities of monocytes/macrophages and T cells in both species in tumour tissues of the CTLA4‐high transcriptomic cohort (Figure ). Then we utilised ROC analysis to determine how well these variables can discriminate between the two transcriptomic subgroups with clinically relevant sensitivity and specificity. Using the Youden index, we identified optimal cut‐off values for parameters with an area under the curve (AUC) of >0.75. ROC curve analysis demonstrated that monocytes and macrophages could differentiate between transcriptomic subgroups with high sensitivity in both humans and dogs and with high specificity in humans. This indicates that IHC using CD68 and IBA1 antibodies detecting monocyte and macrophages could be used to inform clinical decision‐making based on transcriptomic subtype. In addition, CTLA4 itself (but not MET) was associated with a clinically relevant sensitivity in both species and a clinically relevant specificity in humans (Figure and supplementary material, Figure ). Based on these findings and our RNA‐seq analysis, we recommend that the CTLA4‐high subtype OMMs be characterised as immune ‘hot’, with high levels of T‐cell infiltration and activation, and therefore provision of anti‐CTLA4 may show a favourable response in this subgroup . In contrast, our data suggest the MET‐high tumours are immune ‘cold’ and these patients may have a better response to targeted therapy.   OMM The oral cavity hosts the largest and most diverse microbiota after the gut . Given the pivotal role of gastrointestinal tract microbiota in shaping both local and systemic immune system responses , we investigated whether the tumour microbiome correlated with our transcriptomic patient stratification, particularly the CTLA4‐high transcriptomic subgroup due to the high immune cell infiltration. To this end, we performed 16S sequencing of DNA extracted from FFPE cores harvested from our human and canine tumour cohorts as well as four normal canine oral tissue samples. Our analysis revealed a significant difference in alpha diversity between species, reflecting the microbiome's diversity within individual samples, with humans exhibiting significantly lower alpha diversity compared to dogs (Figure , p < 0.01). Beta diversity, a measure of the similarity or dissimilarity of groups of microbial communities, was significantly different between the two species ( p < 0.01; Figure and supplementary material, Data ). Due to the notable differences in alpha and beta diversity observed between species, we conducted separate analyses of the microbiomes associated with human and canine OMM. No significant difference was observed in the microbiome of human or canine OMM when stratified into transcriptomic subtype (Figure ). Consequently, the two transcriptomic subgroups of OMM were not attributed to an overall difference in the composition of the oral microbiota. Thus, the immunocyte infiltration that we identified above is likely attributed to inflammation with the primary tumour itself and not bacterial burden. Although tumour microbiota does not explain the OMM transcriptional subtypes, we did identify two clinically relevant associations that suggest that the tumour microbiome may have a function in OMM disease initiation and progression. First, in human OMM, the overall community diversity was significantly different between metastatic and non‐metastatic tumours (beta diversity p = 0.03; Figure ), with three OTU markers enriched in each group (Figure ). This finding suggests that the microbiota composition may influence the metastatic capacity of OMM . Second, we compared the microbiome profiles of canine OMM with those of normal canine oral tissue and found overall community diversity was significantly different between tumour and normal oral tissue (Figure ). These data indicate that OMM correlates with alterations in the oral microbiota composition, which warrants further investigation because of the emerging role of the microbiome in cancer development and response to therapy .   Our findings provide the conceptional advance that human and canine OMMs are the same disease that can be classified into two subtypes. These two distinct transcriptional subtypes have potential therapeutic implications, and we provide diagnostic tools to distinguish between these two subtypes. Our findings are important for personalised therapeutics for individual patients and for facilitating patient stratification that may inform clinical trial design and outcomes. The advent of T‐cell‐targeted immunomodulators, such as ipilimumab and nivolumab, and small‐molecule MAPK‐pathway inhibitors has revolutionised cancer care for patients with CM . However, the therapeutic efficacy of immune checkpoint inhibitors for MM remains unclear, and patients with MM are underrepresented in clinical trials . Therefore, outcome data for patients with MM are scarce and mainly based on retrospective studies with limited case numbers. Our two transcriptional subtypes suggest that outcomes and clinical trial design could be improved by patient stratification. The CTLA4‐high subtype exhibits elevated CTLA4 expression and increased immunocyte infiltration, indicative of a ‘hot’ immune tumour type. Patients with CTLA4‐high OMM may benefit from anti‐CTLA4 therapy. In contrast, the MET‐high subtype is characterised by high MET expression and, thus, is potentially responsive to small‐molecule tyrosine kinase inhibitors; the lack of immunocyte infiltration suggests this subtype is a ‘cold’ immune tumour type. Ipilimumab is currently used in human patients, but our review of the literature has not identified clinical data regarding ipilmumab or MET inhibitors as therapeutic approaches for MM. We propose a multicentre retrospective review of treatment efficacy utilising our diagnostic tool to stratify patients into subgroups and examining the impact on outcomes or treatment. We developed an inexpensive and readily available diagnostic tool which demonstrates that CTLA4, as well as monocytes and macrophages, distinguishes between transcriptomic subgroups with high sensitivity in both humans and dogs. This suggests that IHC utilising antibodies for CTLA4, CD68, and IBA1 offers practical utility for guiding clinical decision‐making based on transcriptomic profiles. The extent of immune cell infiltration into tumours has important prognostic value in several cancer types . We investigated whether transcriptomic subtype correlated with survival, but, unlike in CM , this did not reach significance in either human or canine patients. Transcriptional stratification alone may not predict survival because distinct mechanisms contribute to poor prognosis in each subgroup . However, notably, we did find that increasing age and metastatic disease was associated with poorer survival in both species. The reason for the survival disparity between genders in our human patients with OMM is unknown. While sex disparities in cancer mortality and survival have been reported, particularly for melanoma and cancers of the mouth and nasal cavity , the specific impact of sex on survival in OMM patients, as observed in our data, remains to be fully understood. Differential environmental exposures and/or physiological processes may explain these disparities . Interestingly, we found no correlation between microbiome diversity and transcriptomic subgroups, suggesting that tumour‐associated mechanisms, rather than bacterial infection, drive immune infiltration. However, the observed alterations in oral microbiota associated with OMM development in dogs and OMM metastasis in humans should prompt further investigation into the role of microbiota in OMM pathogenesis , although we appreciate that caution is advised when applying sequence‐based techniques to the study of microbiota present in low‐biomass environments . In conclusion, we provide valuable insights into the molecular and immune landscapes of OMM and their potential implications for therapeutic interventions. Our findings are potentially enhanced by the homogeneity of our human cohort, and therefore we underscore the importance of future studies prioritising inclusivity by encompassing diverse cohorts . Prospective data collection and tissue sampling efforts will augment the database size and mitigate inherent biases and limitations associated with retrospective data quality and completeness. Moreover, future investigations should explore MM from other anatomical sites to determine the generalisability of our findings across different contexts. The implications of our findings suggest that personalised therapeutic strategies targeting shared molecular pathways could hold promise for both human and canine patients with OMM.   KLBB and EEP conceptualised and designed the study. Data collection was performed by KLBB, JdP, GP, LS, J‐BT, DK, MP, JSM, IB, SZ, SMG, DJS, MST, DM, AM, KP and MC. Data analysis was conducted by KD, YL, DJS, AM, EEP and KLBB. Data interpretation was contributed by LG, MS, MEM, KD, YL, KLBB and EEP. KLBB and EEP conducted the literature search and generated the figures. All authors were involved in writing the paper and had final approval of the submitted and published versions. KLBB and EEP provided supervision and project administration.   Figure S1. Clinical data and Kaplan–Meier survival plots associated with human (Edinburgh) and canine (Bowlt Blacklock) OMM cohorts Figure S2. Two shared transcriptomic subgroups stratify OMM in human and canine patients Figure S3. Violin plots, ROC curves, and Kaplan–Meier survival plots associated with transcriptomic subgroup Data S1. Pathway analysis results associated with two transcriptomic subtypes identified in human and canine OMM Data S2. List of 812 human and canine homologous genes inputted into randomForest machine learning algorithm and genes utilised by the model Data S3. Immunocyte infiltration in human and canine OMM, parsed from annotated gene signature matrix LM22 and 100 permutations of CIBERSORTx web portal Data S4. Beta diversity metrics for microbiome data from human and canine OMM   