In this study, we followed the evolutionary lineage of the contemporary human species since the organisms or taxonomic divisions under investigation represent important links of human evolution: Homo sapiens (human), Pan troglodytes (chimpanzee), Macaca mulatta (Rhesus monkey), Callithrix jacchus (marmoset), Mus musculus (mouse), Rattus norvegicus (rat), Canis lupus familiaris (dog), Equus caballus (horse), Sus scrofa (pig), Bos taurus (cattle), Tursiops truncatus (dolphin), Pteropus vampyrus (bat), Monodelphis domestica (opossum), Ornythorhynchus anatinus (platypus), Gallus gallus (chicken), Xenopus laevis (frog), Latimeria chalumnae (coelacanth), Danio rerio (zebrafish), Callorhinchus milii (shark), Petromyzon marinus (lamprey), Ciona intestinalis (vase tunicate), Strongylocentrotus purpuratus (sea urchin), Amphimedon queenslandica (sponge), Monosiga brevicollis (choanoflagellate), Saccharomyces cerevisiae (baker’s yeast), Schizosaccharomyces pombe (fission yeast), Actinobacteria, Chlamydiae, Cyanobacteria, Proteobacteria, Firmicutes, Bacteroidetes, Archaea.

The official names of the genes were initially retrieved from the HGNC database [32, 33] and then the accession numbers of the peptide sequences corresponding to the human genes were retrieved from the publicly available non-redundant NCBI Reference Sequence Database (RefSeq) [34]. The amino acid sequences were used subsequently as probes in an extensive series of BLASTp [35] reciprocal searches in order to obtain the canonical homologous amino acid sequences corresponding to the species included in this study. This process was reiterated until no novel sequences could be detected, ensuring that a full representation of each gene’s family is obtained.

To investigate the evolutionary history of each gene, we conducted comprehensive phylogenetic analysis using the entire length protein sequences of the species under study. The full-length amino acid sequences were aligned with MAFFT v.7 (https://mafft.cbrc.jp/alignment/server/) [36]. The alignments were subsequently used to reconstruct phylogenetic trees by employing two separate methods, Maximum Likelihood, a method based on a heuristic approach for finding the optimal tree that fits the observed data, and Neighbor Joining, a method based on a hierarchical clustering algorithm [37], as implemented in the software package MEGA (https://www.megasoftware.net/), version 10 [38]. Furthermore, MEGA v.10 [38] was used to estimate the best-fit substitution model, which best describes the number of observed amino acid substitutions per position. For both methods of phylogenetic tree reconstruction, bootstrap analyses (150 pseudo-replicates) were conducted in order to evaluate the robustness and the statistical significance of the inferred reconstructed trees. Finally, the phylogenetic trees were visualized with MEGA v.10 [38].

The pathways and biological processes in which each gene is involved were retrieved from the Reactome Pathway database [39], the KEGG Pathways database [40, 41] and the biomedical literature.

RNA sequencing data for tumor and corresponding normal tissue samples from the TCGA (The Cancer Genome Atlas) and GTEx (Genotype-Tissue Expression) databases, respectively, were retrieved from the GEPIA2 website (http://gepia2.cancer-pku.cn/). The atavistic genes that are differentially expressed (DE) between tumor and normal samples were identified using ANOVA; threshold value for absolute log fold change $\mid$log${}_2$FC$\mid$$\ge$2 and FDR-adjusted p-value (or q-value) $\le$0.05.

The twelve hub genes were provided as input to the STRING (version 11.0) (https://string-db.org/) database [42] of investigating and visualizing both known and predicted gene/protein associations.

By conducting a thorough and comprehensive review of studies published in PubMed (https://pubmed.ncbi.nlm.nih.gov/) and The Human Protein Atlas [43], we found those cancers each of the 12 genes is involved in; even different types of cancers that can affect the same organ (Table 1, Ref. [43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 182, 183, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195, 196, 197, 198, 199, 200, 201, 202, 203, 204, 205, 206, 207, 208, 209, 210, 211, 212, 213, 214, 215, 216, 217, 218, 219, 220, 221, 222, 223, 224, 225, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 236, 237]). In this context, certain genes have been shown to exhibit either oncogenic or tumor suppressive effects. In fact, both qualitative and quantitative modifications were identified in the genes under study which are associated with the development and progression of different cancers in different tissues. In particular, the types of alterations observed in those genes and/or their protein products, include aberrant, increased or decreased, expression (Table 2), epigenetic modifications, mutations, gene inactivation or amplification, copy number alterations, altered protein interactions and post-translational modifications, characteristic protein subcellular localization. Furthermore, these genes constitute diagnostic and prognostic biomarkers in cancer, regarding malignancy, disease stage, clinical subtypes, disease progression, metastasis, clinical outcome and survival and prediction of patients’ response to therapeutic treatments. In addition, each gene and/or its products constitute drug targets or are considered as novel potential drug targets for at least one type of cancer.

The critical role of these genes and their corresponding products in carcinogenesis is highlighted by the fact that they are all involved in cancer-relevant pathways and processes, usually more than one (Table 3). More specifically, the identified cancer signaling pathways include the Wnt pathway (RCC2, DSP, CTNND1), the MAPK pathway (TLN1, CTTN, PKP3), the VEGF pathway (TLN1, CTNNA1, CTNND1), the PI3K-Akt pathway (VASP, ILK, PKN2), the Hippo pathway (ACTG1, CTNNA1), the PPAR pathway (ILK) and the Rap1 pathway (VASP, CTNND1, ACTG1, TLN1). Moreover, particularly vital processes of the cell, which also play a role in carcinogenesis, are eukaryotic cell cycle regulation and mitosis (RCC2, ILK, CTNNA1, CTNND1), apoptosis (ACTG1, PLEC, DSP, CTTN, CTNND1), immune system processes (TLN1, VASP, ACTG1, DSP, PKP3) and interaction with proteoglycans in cancer (ACTG1, CTTN). Finally, the proteins that constitute the adherens junctions (ACTG1, CTNNA1, CTNND1) and contribute to focal adhesions (VASP, TLN1, ACTG1) play a critical role especially in metastasis.

Of note, based on the STRING database analysis, the most significantly enriched pathway of the 12 gene/proteins is ‘cell junction assembly’ (GO: 0034329; FDR = 4.10 $\times$ 10${}^{-15}$). In this network (Fig. 1), ten hub genes/gene products appear to interact either physically or functionally in the process of cell juction organization, having an integral role in EMT which is essential in cancer progression and metastasis [238].

Therefore, the importance of the examined genes is enhanced by their involvement in pathways related to carcinogenesis and the progression of the disease in various types of cancer.

In terms of overall topology, there is congruence between the phylogenetic trees that were generated with both methods. The trees constructed with the maximum likelihood method are considered more accurate and reliable and they display higher bootstrap values (Figs. 2,3,4, 5,6,7,8,9,10,11,12). These high bootstrap values indicate that the tree nodes are statistically significant and the inferred topology is biologically significant. The trees generated with the neighbor-joining method can be found at the Supplementary material.

The ILK, CTNNA1, CTNND1 and PKP3 genes are detected exclusively in multicellular Animals as shown in the phylogenetic trees in Figs. 2,3,4, respectively. Therefore, those genes most likely first appeared in a eukaryotic multicellular organism which was the common ancestor of Metazoa.

Catenins A (CTNNA) and catenins D (CTNND) are members of the catenin superfamily (Fig. 3). CTNNA1, CTNNA2, CTNNA3, CTNND1, CTNND2 genes are paralogs which probably occurred through a series of gene duplication events. The corresponding orthologous proteins of each gene form distinct clades, and as a result the phylogenetic tree is divided into two major subtrees, comprised of protein sequences encoded by CTNNA and CTNND genes, respectively. The subtree of CTNNA includes the protein sequences encoded by CTNNA1, CTNNA2 and CTNNA3. CTNNA2 appears to be the primordial gene as it was first detected in Amphimedon queenslandica, thus it probably first emerged in an ancestor of Porifera after their divergence from Choanoflagellates, since an ortholog was not found in Monosiga brevicollis. On the other hand, CTNNA1 and CTNNA3 were detected for the first time in Callorhinchus milli, and thus, they appeared later in evolution and probably occurred due to duplication events of the CTNNA2 gene, as they most likely first appeared in an ancestor of Chondrichthyes after their divergence from Tunicates, since orthologs were not detected in Ciona intestinalis and in fact they demonstrate high similarity to each other. The subtree of catenins D includes the protein sequences encoded by CTNND1 and CTNND2. CTNND2 is apparently the primordial gene as it was first detected in Amphimedon queenslandica, and thus it probably first arose in an ancestor of Porifera after their divergence from Choanoflagellates, since an ortholog was not detected in Monosiga brevicollis. On the other hand, CTNND1 appeared later in evolution and probably occurred due to CTNND2 gene duplication, as it was detected for the first time in Ciona intestinalis. Therefore, CTNND1 might have arisen in an ancestor of Tunicates after their divergence from Echinodermata, given that CTNND1 orthologs were not detected in Strongylocentrotus purpuratus.

The PKP1, PKP2, PKP3, PKP4 genes are likely paralogs (Fig. 4) which probably occurred through a series of gene duplication events of an ancestral PKP gene. The corresponding orthologs of each gene form distinct clades. All PKP genes were detected for the first time in Callorhinchus milli, suggesting that they probably appeared in an ancestor of Chondrichthyes after their divergence from Tunicates, since PKP orthologs were not detected in Ciona intestinalis. Based on the inferred phylogenetic tree, however, PKP1 and PKP2 exhibit the highest similarity, followed by PKP3. We could speculate that PKP4 is the primordial gene of this family, since it is basal to PKP1, PKP2 and PKP3. PKP3 appeared later in evolution, probably due to PKP4 gene duplication, while PKP1 and PKP2 have emerged even more recently through another series of duplication events in a chondrichthyan PKP1/2 gene and have not accumulated a large number of mutations.

The VASP, CTTN and DSP genes are detected during the relatively late transition phase from unicellularity to multicellularity, based on the inferred phylogeny shown in Figs. 5,6,7, respectively.

VASP and Enah/Vasp-Like are likely paralogs (Fig. 5). In Monosiga brevicollis, Strongylocentrotus purpuratus, Ciona intestinalis and Callorhinchus milli the corresponding detected genes were annotated as VASP-like in RefSeq, but display high similarity to the VASP genes of the rest of Metazoa.

The PLEC and TLN1 genes were detected in unicellular and multicellular organisms, as shown in Figs. 8,9, respectively, and they probably first appeared in a eukaryotic unicellular organism that was the common ancestor of Animalia and Fungi.

PLEC homologs were found in Animalia as well as in Fungi (SAC6 in Saccharomyces cerevisiae and FM1 in Schizosaccharomyces pombe) (Fig. 8).

The TLN1 and TLN2 genes are likely paralogs (Fig. 9), and their corresponding orthologous protein sequences form distinct clades. TLN1/2 is probably the primordial gene of this family, since it was first detected in Amphimedon queenslandica, Monosiga brevicollis, Ciona intestinalis and Strongylocentrotus purpuratus which gave rise to TLN1 and TLN2 after successive gene duplications. Talin homologs were also identified in Fungi, that is, SLA2P in Saccharomyces cerevisiae and END4 in Schizosaccharomyces pombe, as also confirmed in literature [239, 240]. Therefore, although members of the TLN1 gene family were found only in Metazoa, it has deep evolutionary roots in some ancient unicellular eukaryotic organism, which is the common ancestor of Animalia and Fungi.

The characteristically long clade corresponding to Monosiga brevicollis, as shown in Figs. 5,6,7,8,9, is due to the evolutionary course of the organism, as it evolved for many millions of years separately from the rest of the Metazoa.

The evolutionary origin of ACTG1, RCC2 and PKN2 genes can be traced to a universal common ancestor, namely an ancient prokaryotic unicellular organism, that is the common ancestor of Animalia, Bacteria and Archaea as illustrated in Figs. 10,11,12, respectively.

The ACTA1, ACTA2, ACTB, ACTC1, ACTG1, ACTG2 genes are likely paralogs (Fig. 10), the orthologs of which cluster together in distinct clades. The phylogenetic tree is divided into two major subtrees, one comprised of protein sequences encoded by genes of different types of cytoplasmic actins and the other includes sequences encoded by genes of different types of muscle actins.

The subtree of cytoplasmic actins harbors the protein sequences encoded by ACTG1 (Actin cytoplasmic 2), ACTB (Actin cytoplasmic 1) as well as various other actin genes. ACTG1 is apparently the ancestral gene of this family as it was detected in Animalia, Bacteria (Microbacterium arborescens (Actinobacteria), Chlamydia trachomatis (Chlamydiae), Leptolyngbya sp. PCC 7376 (Cyanobacteria), Kangiella spongicola (Proteobacteria), Staphylococcus aureus (Firmicutes)), as well as Archaea (Candidatus Heimdallarchaeota (Archaea), Candidatus Lokiarchaeota (Archaea)). The sequences detected in bacteria do not cluster together, but are rather scattered among those corresponding to genes that encode the two types of cytoplasmic actins in Metazoa. This is indicative of high similarity among the aforementioned genes, since the encoded protein sequences appear to be highly conserved. Additionally, Actin genes found in Fungi (Saccharomyces cerevisiae and Schizosaccharomyces pombe) exhibit high similarity to Actin genes found in Archaea, which encode cytoplasmic actin 2.

On the other hand, the ACTB gene likely appeared later in evolution and probably occurred due to a duplication event of the ACTG1 gene, as it was detected for the first time in Petromyzon marinus. Also, in Ciona intestinalis, a cytoplasmic actin gene and a muscle actin gene were identified, and in fact the gene encoding the cytoplasmic actin shows high similarity to the only ACTIN gene detected in Monosiga brevicollis. Finally, Actin-85C-Like gene found in Amphimedon queenslandica also demonstrates high similarity to types of cytoplasmic actin.

The subtree of muscle actins, includes the protein sequences encoded by ACTA1 (Actin alpha 1, skeletal muscle), ACTA2 (Actin alpha 2, smooth muscle), ACTC1 (Actin alpha cardiac muscle 1) and ACTG2 (Actin gamma 2, enteric smooth muscle). Muscle actin genes appear to be the primordial genes of all muscle actin encoding genes as they were first identified in Strongylocentrotus purpuratus and Ciona intestinalis; thus, the gene probably first appeared in a common ancestor of Echinodermata and Chordata after their divergence from Porifera, since an ortholog was not detected in Amphimedon queenslandica. ACTA2 and ACTC1 genes were first detected in Callorhinchus milli, ACTA1 in Danio rerio, and ACTG2 in Latimeria chalumnae, respectively. We would expect to find the ACTA1 and ACTG2 genes in Callorhinchus milli, as it is the first organism in this study to have differentiated organs, but the sequences were not detected using BLAST. In any case, we believe that these genes occurred through a series of gene duplication events yielding four paralogs that probably first appeared in an ancestor of the Chondrichthyes after their divergence from Echinodermata and Tunicates.

In summary, those genes encoding different types of cytoplasmic actins are evolutionarily older, since those genes that encode different types of muscle actins appeared later in evolution. Specifically, ACTG1 is the ancestral gene of the entire family, since it is the evolutionarily older out of all paralogs found in Metazoa.

The RCC1 and RCC2 genes are likely paralogs, and their corresponding orthologous protein sequences form distinct clades (Fig. 11). RCC1 is probably the primordial gene of this family as it was detected in Metazoa, Bacteria Bifidobacterium asteroides (Actinobacteria), Parachlamydia sp. C2 (Chlamydiae), Synechococcus sp. WH 8103 (Cyanobacteria), Myxococcus xanthus (Proteobacteria), Cohnella sp. CAU 1483 (Firmicutes), Hymenobacter chitinivorans (Bacteroidetes) and Methanocella conradii (Archaea). RCC1 and RCC2 homologs have also been identified in Fungi, that is, SRM1 in Saccharomyces cerevisiae and PIM1 in Schizosaccharomyces pombe, as confirmed in literature [241, 242]. Furthermore, although the corresponding gene in Monosiga brevicollis and Amphimedon queenslandica is annotated as RCC in RefSeq, it displays high similarity to the RCC1 gene. Thus, RCC1 probably first appeared in a common ancestor of Eukaryotes, Bacteria and Archaea. RCC2 on the other hand, might have appeared later in evolution, as it was detected for the first time in Strongylocentrotus purpuratus, and probably occurred due to duplication of the RCC1 gene in an ancestor of Echinodermata after their divergence from the phylum of Porifera, since an RCC2 ortholog was not detected in Amphimedon queenslandica. Therefore, although RCC2 is found only in Metazoa, its ancestor is traced to some ancient prokaryotic unicellular organism, that is, the common ancestor of Eukaryotes, Bacteria and Archaea.

The PKN1, PKN2, PKN3, PKN, PKC genes are likely paralogs, the orthologs of which cluster together in distinct clades (Fig. 12). PKN was detected in Amphimedon queenslandica, in the unicellular Choanoflagellate Monosiga brevicollis, in Bacteria Streptomyces seoulensis (Actinobacteria), Chlamydia ibidis (Chlamydiae), Chloracidobacterium thermophilum (Cyanobacteria), Escherichia coli (Proteobacteria), Paenibacillus donghaensis (Firmicutes), Spirosoma panaciterrae (Bacteroidetes) and in Archaea (Methanoregula boonei and Thermococcus thioreducens) and is likely the primordial gene of this family and precursor of the PKN1, PKN2 and PKN3 genes which probably occurred through a series of PKN gene duplication events. The PKC gene was also detected in Petromyzon marinus, as well as in Fungi (Saccharomyces cerevisiae and Schizosaccharomyces pombe). PKN2 is apparently the oldest one as it was first detected in Ciona intestinalis and Strongylocentrotus purpuratus, so it probably first appeared in a common ancestor of Tunicates and Echinodermata after their divergence from Choanoflagellates, since an PKN2 ortholog was not detected in Monosiga brevicollis. PKN1 most likely emerged later as it was first detected in Callorhinchus milli; it probably first appeared in an ancestor of Chondrichthyes after their divergence from Tunicates, since an PKN1 ortholog was not detected in Ciona intestinalis. PKN3 was detected for the first time in Latimeria chalumnae, so it probably first emerged in an ancestor of Actinists after their divergence from Osteichthyes, since an PKN3 ortholog was not detected in Danio rerio. Therefore, although PKN2 was found exclusively in Metazoa, it has very deep evolutionary roots found in some ancient prokaryotic unicellular organism, that is the common ancestor of Animalia, Bacteria and Archaea.