Epidemiology and Molecular Biology of HPV Variants in Cervical Cancer: The State of the Art in Mexico

Cervical cancer (CC) continues to be a major public health problem in Mexico, ranking second among cancers in women. A persistent infection with human papillomaviruses (HPV) is the main risk factor for CC development. In addition, a significant fraction of other cancers including those of the anus, oropharynx, and penis are also related to HPV infection. In CC, HPV-16 is the most prevalent high-risk HPV type, followed by HPV-18, both being responsible for 70% of cases. HPV intratype variant lineages differ in nucleotide sequences by 1–10%, while sublineages differ by 0.5–1%. Several studies have postulated that the nucleotide changes that occur between HPV intratype variants are reflected in functional differences and in pathogenicity. Moreover, it has been demonstrated that HPV-16 and -18 intratype variants differentially affect molecular processes in infected cells, changing their biological behavior that finally impacts in the clinical outcome of patients. Mexico has participated in providing knowledge on the geographical distribution of intratype variants of the most prevalent HPVs in premalignant lesions of the cervix and cervical cancer, as well as in other HPV-related tumors. In addition, functional studies have been carried out to assess the cellular effects of intratype variations in HPV proteins. This review addresses the state of the art on the epidemiology of HPV-16 and HPV-18 intratype variants in the Mexican population, as well as their association with persistence, precancer and cervical cancer, and functional aspects related to their biological behavior.


Introduction
Among the seven known human viruses associated with cancer, high-risk human papillomaviruses (HR-HPV) cause the greatest cancer burden, with more than 690,000 cancer cases attributed to HPV in 2018 [1]. Cervical cancer (CC) is the most common neoplasia related to HPV infections. A persistent infection with HR-HPV is the main risk factor for the development of most CCs (>85%) [2,3]. Additionally, other anogenital cancers, such as anal, vulvar, vaginal, penile, and a proportion of oropharyngeal and oral cancers, are also caused by HPV [4,5]. Therefore, HPV-related cancers affect both the female and the male populations.
Cervical cancer continues to be a public health problem, and it ranks fourth among cancers in women worldwide [6]. In Mexico, CC ranks second among cancers in women, where 9439 new cases and 4335 deaths from this neoplasia were estimated in 2020 [7].
Human papillomaviruses (HPV) are small non-enveloped particles which infect skin and mucous membranes. Virus particles contain approximately 8000 bp of double-stranded DNA, covered by a viral capsid composed of L1 and L2 proteins. According to their genomic sequence similarity, more than 200 types of HPV have been described [8], considering that each type of HPV differs from another by at least 10% in the highly conserved L1 gene [9]. Based on such similarity, HPVs are divided into five genera (α, β, γ, µ, and ν), where α-papillomaviruses predominantly contain HPV types that infect epithelial mucosa [10]. In addition, HPVs are classified into high-risk (HR) or low-risk (LR) types based on their oncogenic potential. Approximately 20 high-risk types can be found in HPV-related cancers, although HPV-16, is the most common type found in cervical cancer (50%) and premalignant lesions of the cervix as well as in 80% of oropharyngeal cancers [11][12][13][14]. The high prevalence of HPV-16 in cancer in relation to other oncogenic HPVs may be due to a greater oncogenicity, partly explained by the differential effects exerted by HPV-16 oncogenes on cellular targets. After HPV-16, the second most carcinogenic HPV type is HPV-18, which is found in approximately 12% of squamous cell carcinomas and 37% of cervical adenocarcinomas worldwide [15].
Sequence diversity is also found within each HPV type, where intratype variant lineages differ by 1-10%, while sublineages differ by 0.5-1% [16,17]. Being double-stranded DNA viruses, papillomaviruses use the efficient proofreading DNA polymerase of the host for their replication, which avoids high mutation rates; therefore, changes in HPV genomes have been acquired slowly, defining HPV types whose infection cycle has adapted to different host cells. Random mutations that eventually occur within viral types have generated the nucleotide polymorphisms contained in HPV intratype variants. When the intratype variants of HPV-16 and -18, began to be studied, they were named African, European, Asian, Asian-Amerindian, depending on the geographical origin and ethnicity of the infected population. Subsequently, several studies showed that HPV variants are spread throughout the world, and a more precise denomination was proposed, changing the nomenclature to lineages and sublineages [17].
Many HPV intratype variants have been isolated in various clinical studies and a careful classification of variants of the most frequent HPV types in cancer has been achieved [9,15,17]. It has been postulated that nucleotide variations that occur among HPV intratype variants are reflected in functional differences and in their pathogenicity [18][19][20][21]. Currently, there are many studies that report that HPV-16 and -18 intratype variants differentially alter host cells, which translates into different risks of persistence and progression to cancer [22][23][24][25][26][27][28][29].
Several epidemiological studies have been carried out on the geographical distribution of HPV types and intratype variants worldwide and Mexico is not the exception. In addition, functional studies have also been carried out regarding the participation of viral proteins and the transcriptional activity of the viral non-coding region with the goal of identifying differences in biological behavior towards the establishment of cancer, mainly focused on HPV-16 and -18 intratype variants.
Deciphering the impact that HPV intratype variations have on cancer development, as well as on the clinical outcome of patients with HPV-related cancers, could eventually influence patient prognosis and the development of therapeutic strategies. This review addresses the state of the art on the epidemiology of HPV-16 and HPV-18 intratype variants worldwide and particularly in the Mexican population, as well as their association with persistence, precancer and cervical cancer, and functional aspects related to their biological behavior.

HPV Infection and Life Cycle
The human papillomavirus genome consists of three functionally different regions: an early expressed region (E) that includes genes E1, E2, E4, E5, E6, and E7; a late expressed region (L) that contains the L1 and the L2 genes; and a long non-coding region (LCR) that contains elements that regulate viral replication and transcription [30].
The HPV replicative cycle is dependent on the differentiation process of the epithelium. During the initial HPV infection, E viral genes are expressed from the undifferentiated basal epithelium where E1 and E2 products regulate viral genome replication and E2 controls the expression of E6 and E7 oncoproteins. In addition, E6 and E7 oncoproteins interact with a variety of cellular proteins that favor the replication of the viral genome. L1 and L2 proteins are produced in the most differentiated layers of the epithelium and after being assembled, virions are released by cellular desquamation in the most superficial layers.
To establish a productive infection, HPV must infect the basal cells of a differentiating epithelium, where the viral genome is amplified and subsequently packaged, allowing for the generation of new viral particles. Through micro-wounds, HPVs reach and infect mitotically active basal cells that will eventually divide and enter the parabasal layer, undergoing differentiation. A specific receptor that mediates the entry of the virus into the host cell has not yet been described, although different receptors have been proposed [31]. HPV entry into the cell occurs by endocytosis independent of clathrin, caveolin, lipid raft, and dynamin [32]. It has been proposed that HPV binds to cell-free Heparan Sulfate (HS) or Heparan Sulfate Syndecan (Sdc) proteoglycan (HPSG), Sdc2, and Sdc4, bound to the cell membrane [33]. In addition, cell membrane receptors have been identified, including EGFR [34], α6-integrin [35], CD63 [36] and CD151 tetraspannin [37], and annexin A2/S100A10 heterotetramer (A2t), which are required for HPV uptake [38].
HPV interaction with HPSG exposes L2 on the capsid surface and mediates conformational changes in L1, allowing the kallikrein-8 (KLK8) protease to cleave L1, thereby exposing the RG-1 epitope within the N-terminus of the L2 protein [39]. Furin protease then cleaves L2 upstream of the RG-1 epitope [40], preparing the viral particle for entry and proper downstream trafficking of L2 that binds to the HPV episome. L2 then binds to γ-secretase protease and p120-catenin that act as chaperones of L2, allowing insertion of L2 into vesicular membranes [41,42]. Once in the endosome, L2 interacts with different proteins that ensure vesicular trafficking of the L2-HPV episome. For example, L2 has been shown to interact with Sortin Nexin 17 (SNX17), which prevents vesicle acidification and rapid lysosomal degradation of vesicle content, hence the L2-episome maintains its integrity [43]. Additionally, L2 binds to members of the retromer cargo complex, including Vps26, Vps29, and Vps35 [44], essential for retrograde trafficking of L2-HPV episomes to the Trans Golgi Network (TGN). Then, L2-HPV episomes exit the endosomal compartments and traffic to the TGN where they remain until the onset of mitosis [31]. During mitosis, the spindle microtubule is required to transport the L2-HPV episome, reaching the mitotic chromosomes and centrosomes, since L2 interacts with motor proteins Ran-binding protein 10 (RanBP10), Karyopherin alpha2 (KPNA2), and dynein light chain (DYNLT3) [45,46]. Subsequently, L2 binds to the cellular mitotic chromosomes, allowing L2-dependent chromosomal anchoring of the L2-HPV episome during open mitosis [47]. Finally, the viral episome is delivered into the nucleus of daughter cells in a mitosis-dependent manner, where HPV episomes are localized to highly transcriptionally active sites termed ND10/PML bodies [48].
Once inside the nucleus, HPV transcription is driven by the p97 and p105 promoters of HPV-16 and -18, respectively, allowing the expression of E6, E7, E1, and E2 genes, necessary for the early stages of viral genome amplification. E2 protein contains a DNA binding domain, and once it binds to its response elements located in the LCR, it negatively regulates the expression of the early expressed viral genes [49]. Meanwhile, E1 protein acts as an ATP-dependent DNA helicase, and it is considered the only enzyme encoded by papillomaviruses [50]. E2 interacts with E1, forming hexameric dimers that bind to the origin of replication and unwind the DNA ahead the replication fork, allowing recruitment of the replication machinery. The viral episomes are amplified in a low number of copies (50-100 episomes per nucleus) [51], which favors the evasion of immunological surveillance and, consequently, the establishment of persistent infections. Additionally, it has been shown that E1 proteins suppress the expression of genes related to immune response, such as IFNβ1 and IFNλ1 and Interferon-Stimulated Genes (ISG) [52].
Importantly, the expression of E6 and E7 is necessary for infected cells to enter mitosis, and the amplification of the HPV genome partially depends on the cell division program. Although the activity of E6 and E7 oncoproteins is crucial for the establishment of cancer, their expression is also necessary to carry out the replicative cycle of HPV. E7 is responsible for maintaining the continuity of the cell cycle, promoting the transition of the G1 to S phase, since E7 interacts with the tumor suppressor pRb and with the E3 ubiquitin ligase Cullin 2, which induces the degradation of pRb [53]. These events allow the release of the transcriptional factor E2F from its inhibitor, pRb, promoting the expression of target genes, such as cyclin E, cyclin A, and p16INK4A, among others [53].
Generally, E6 and E7 proteins are produced from a bicistronic mRNA and the activities of both viral oncogenes in HPV-infected cells ensure the generation of viral progeny. In this context, it may be that when E7 induces excessive proliferation, then E6 could contribute by mitigating such effect.
It has been described that E6 interacts with the tumor suppressor p53, forming a complex with E6AP (E6-associated protein), promoting the degradation of p53 through the proteasome. This event impedes p53 functions, such as induction of apoptosis and DNA damage repair [53].
Furthermore, viral E5 protein maintains cell proliferation in the upper layers of the epithelium. E5 prevents endosomal acidification by binding to the 16-kDa subunit of vacuolar ATPase (vATPase), which allows the continuous membrane recycling of the EGFR receptor, producing its activation and the activation of downstream effectors. In addition, E5 promotes immune system evasion by preventing antigen presentation by the MHCI, since it indirectly interacts with the heavy chain of MHCI, preventing its correct trafficking to the membrane resulting in its accumulation in the Golgi apparatus and endoplasmic reticulum [57]. At this point of the infection, HPV episomes increase up to thousands of copies per cell.
Meanwhile, the E4 protein has been shown to promote the collapse of the cytokeratin network, which makes cells fragile and more likely to release progeny virions [58]. Moreover, keratin has been shown to be hyperphosphorylated at residues K8 and S73, as well as ubiquitinated in the presence of E4, suggesting that its degradation could be promoted through the proteasome, allowing network disruption [59].
Finally, viral capsid proteins are synthesized from the late promoter p670 and p811 for HPV-16 and -18, respectively, and assembled in terminally differentiated cells. The HPV genome is encapsidated in an icosahedral outer shell made up of 360 L1 molecules organized into 72 pentameric capsomeres, where L2 is localized at the center [60].
During the infection process, reactive oxygen species, genomic instability, and DNA strand breaks are generated, which can favor the integration of the viral genome into the cellular genome. This integration process promotes the deregulated expression of E6 and E7 oncogenes, which is considered an initial step for cervical carcinogenesis. The interaction of E6 and E7 with tumor suppressors and with proteins involved in cellular signaling pathways alter normal cell function, leading to transformation [53,61].

HPV Pathogenesis
Most women will be infected by HPV at least once in their lives [62]; however, it has been reported that about 90% of these HPV infections are transient and will disappear within the first two years due to the host immune system [63]. Although progression to cancer is relatively rare, a persistent HPV infection is a key event determining progression to cervical cancer. Therefore, cancer caused by HPV is a sporadic event which is not beneficial for the viral replicative cycle or for the generation of viral progeny.
A minority of infections become persistent (10-20%) and are at increased risk of progression to cancer [64]. It has been demonstrated that in women who were persistently infected with one HR-HPV and tested over 2 years, the risk of cervical cancer was 12.4% compared to 0.14% in women with a repeatedly negative HPV test. Moreover, the risk of developing cancer increases with age, being 5.5%, 14.4%, and 18.1% for women aged 30 to 44 years, 45 to 54 years, and 55 years or older, respectively [65].
Among the factors proposed to determine persistence that are inherent to the virus are viral load, viral type, and intratype variations, while those inherent to the host include genetic predisposition and lifestyle factors [66].
Prior to the establishment of cancer, a series of clinical characteristics inherent to HPV infection occur, where the normal epithelium evolves into epithelial precursor lesions, causing cervical disease. Such lesions have been histologically classified into three groups, according to the progression of the lesion: cervical intraepithelial neoplasia (CIN) 1, 2, and 3 (mild, moderate, and severe epithelial dysplasia, respectively). Additionally, another classification based in cytological evaluation, divides the lesions into two groups known as low-grade squamous intraepithelial lesions (mild dysplasia) (LGSIL) and high-grade squamous intraepithelial lesions (moderate to severe dysplasia) (HGSIL) [67].
CIN1 or LGSIL lesions have a low probability of progressing to cervical cancer and are caused by transient HPV infections and 90% of them revert to healthy epithelium. Conversely, CIN2 and 3 or HGSIL are caused by persistent and unproductive HPV infections, yet 60% of such lesions clear spontaneously in immunologically competent individuals [68]. In this scenario, 0.6% of HPV infections are known to progress to cancer. A study with 16 years follow-up showed that women with persistent infections with any type of carcinogenic HPV have a 75.4-fold increased risk of developing cancer compared to HPV-negative women [65]. Interestingly, some HPV-16 and -18 variants have been associated with a high risk of persistence, which may be influencing the clinical outcome of the infection.

Classification of HPV-16 and -18 Lineages and Sublineages
Among the first studies that described the existence of intratype variants of HPV-16, that performed in 1993 by Ho et al. [69] stands out. Based on the LCR sequence, they classified intratype variants into five main branches, named according to the geographical origin of the populations where the variants were found with a higher prevalence: European (E), Asian (As), Asian-American (AA), and two African lineages (Af1 and Af2). In 1997, Yamada et al. [70,71]  HPV-18 variants were initially classified according to nucleotide changes in the L1 gene as European (E), Asian-Amerindian (AsAi), and African (Af) [76]. Subsequently, through phylogenetic analysis of the complete LCR and the E6 gene, the HPV-18 variants were grouped into 3 lineages (A, B and C) and 9 sublineages, where A1/A2 were previously AsAi variants; A3/A4, were the European variants; and the B sublineages were formerly the African variants [15].
According to the Papillomavirus Episteme (PaVE) portal [74], the E6 variants of HPV-16 contain at least 17 nucleotide changes in their sequence with respect to the reference variant (A1), while in E7 at least 7 mutations have been observed. Meanwhile, for HPV-18, at least 14 and 9 mutations have been reported in the E6 and E7 variants, respectively, relative to the reference variant.
Regarding HPV-16, the greatest number of variations are detected in the B4 sublineage with 7 nucleotide substitutions in the E6 gene, impacting in 4 amino acid changes. Meanwhile for the E7 gene, sublineages B2, C1, C4, D2 and D3 present 3 mutations each. In the E7 protein, the A4, B2, C1 and C4 sublineages were affected with one amino acid change. In relation to HPV-18, the B1 sublineage is the most different in relation to the A1 reference sublineage, with 9 nucleotide changes in E6 gene and 2 amino acid substitutions; while in the E7 gene, B1, B2 and B3 sublineages exhibit 4 sequence mutations. Interestingly, the B2 sublineage contains 3 amino acid substitutions in the E7 protein, being the most different variant. The nucleotide and amino acid alignments of HPV-16 and -18 sublineages available in PaVE are described in detail in the Supplementary Material.
In regard to HPV-18, several studies, although controversial, support an association of certain sublineages with infection persistence or CC risk. A Brazilian study shows that HPV-18 European variants are more likely to develop persistent infections [26]; in contrast, a study with a Portuguese population indicates that HPV-18 European variants are associated with a lower risk for CC development [84]. Xi et al. (2007) [85] stated that the risk of developing > CIN3 was higher for European and Asian-American variants, compared with HPV-18 Af variants. Another study, which analyzed cervical samples from the International Agency for Research on Cancer (IARC) biobank collected from different regions worldwide, affirms that HPV-18 variant sublineages cannot discriminate cancer risk [15].
Furthermore, it has been determined that some HPV-16 intratype variants are preferentially associated with specific histological types of cancer [22,86]. However, the information is controversial for HPV-18 variants, since some authors support an association of certain variants with histological types of CC [87][88][89], while others refute it [15].
Knowing the distribution of high-risk HPV intratype variants in different populations, particularly of HPV-16 and HPV-18, will help to identify their risk of developing cervical cancer. A study by Mirabello et al. (2016) included cervical smears from 3200 HPV-16 positive women who routinely attended a clinical center in the US [22]. In this study, the association of HPV variant lineages with precancer and cancer was analyzed. The authors found that 87.2% of the women were infected with an HPV-16 A variant, including European (A1-A3) and Asian (A4) variant sublineages. The presence of A4 variants was significatively associated with an increased risk of CC, compared to A1/A2 (OR = 3.16, 95% CI = 1.05 to 9.54, p = 0.04); while women with non-A sublineages (B, C, D) had a higher risk of CIN 3 and CC (OR = 6.87, 95% CI = 4.42-10.68, p < 0.0001), which is in agreement with data that non-European variants are associated with a higher risk of developing cancer.
One of the studies with the largest number of samples carried out by Clifford G. et al. (2019) [86], and it analyzed 7116 cervical smears positive for HPV-16 from women from 52 countries. The authors found that the A lineage was the most prevalent, present in 78.7% of the samples, followed by lineages D (9.2%), C (6.4%) and B (5.8%). Regarding the distribution of HPV-16 sublineages, A1 was the most widespread, in Europe, South/Central America, North America, South Asia and Oceania, while sublineages A3 and A4 were mainly prevalent in East Asia. Meanwhile, B and C sublineages were only found in African samples. Increased risk for developing CC was observed for A3, A4 and D sublineages, compared with A1, in regions where those sublineages were more prevalent. When analyzing associations with histological types, sublineage D was more frequent in adenocarcinoma than in squamous cell carcinoma.
All this information supports that even though the distribution of the different variant sublineages of HPV-16 varies widely throughout the world, non-European variants have a higher risk of developing high-grade premalignant lesions and CC.

Distribution of HPV Variants in Mexico
In Mexico, some studies have been performed to determine HPV intratype variant distribution. An initial work, based on E6 and L1 genes stated that nearly 40% of HPV-16 positive cervical cancer cases in Mexican women belong to Asian-American lineages D2 (AA-c) and D3 (AA-a), which also presented a higher risk of developing cervical cancer than European lineages A1/2 [90], with an OR of 27 and 3.4, respectively. Interestingly, D2 was associated with a younger age of CC development (younger than 50 years old) than D3 [91].
Similar results were found when the HPV-16 E6 gene was sequenced in 40 cases of CC and cervical premalignant lesions from women from the southeast of Mexico (Yucatan) [92]. In this study, the most frequent variants found overall were European-prototype (E-P) (sublineage A1), followed by E-350G (sublineage A2), with no differences among the groups. Furthermore, the most frequent variant in CC was AA (D2), with a 44% rate of positivity in this group, which in turn was not found in premalignant lesions.
Another analysis carried out in 277 samples of premalignant lesions and 133 samples of CC of women who attended primary health care centers and a cancer institute in Mexico City [93] showed that in HPV-16 positive cases, the most common E6 variants were E, followed by AA-a, although no differences were found in their distribution among the groups. When analyzing AA variants together (AA-a plus AA-c), their proportion in CC was up to 41%, similar to that of E variants in this group. Interestingly, although in a low proportion, AA-c (9%) was found exclusively in CC, and not in premalignant lesions, suggesting a more aggressive behavior for this variant.
Subsequently, Escobar Escamilla et al. [94] analyzed the LCR and E6 sequences to study the distribution of HPV-16 lineages in samples with normal cytology, squamous intraepithelial lesions, and CC from 94 patients in Mexico City, the state of Mexico, and the state of Jalisco. Sublineage A was also the most prevalent, representing 82.93% of the samples analyzed, including combined sublineages A1 and A2. The second place in prevalence was occupied by D2 variants (9.75%). Regarding LCR analysis, the most common mutations found were G7191T and G7518A within sublineages A1 and D2, respectively. On the other hand, variant E6-T350G (E-350G) was observed in 80.48% of the cases, where 100% of the CC samples contained this variant. Interestingly, in asymptomatic women in the state of Sinaloa (northwestern Mexico), the E-350G variant has been reported as the most prevalent among HPV-16 variants [95]. This variant has also been found frequently in women with other conditions, such as those analyzed in a study that included samples from the oral cavity and cervix of HIV+ women in Mexico City [96], where an HPV prevalence as high as 96.6% in the cervix and 92.5% in the oral cavity, mainly with high-risk HPV types was found. In 52 of the 174 women from this group that were HPV-16 positive, the most frequent E6 intratype variant was E-350G, being present in both anatomical sites in 80.77% of the cases.
In contrast, in a population of northeastern Mexico (state of Nuevo Leon) where cervical smears of women attending primary health care centers were analyzed, the most common HPV-16 variants found had an AA origin (87%), while European variants were present in 13% of the samples [97]. This confirms that HPV-16 variant distribution varies within the different geographical areas of Mexico.
Little information is available regarding the distribution of HPV-18 variants in Mexico and their association with the risk of developing CC. In a study including 33 cases of premalignant lesions and CC that were positive for HPV-18, significant data was obtained in this regard. The HPV-18 Af variant predominated in normal cytologies and premalignant lesions, while the E variant was exclusively found in CC, being present in 62% of CC cases, which suggests a higher oncogenicity for this E variant [93].
Moreover, cervical swab analysis from a cohort of women attending primary health care centers in Monterrey, Nuevo Leon, showed that in 15 cases positive to HPV-18 the most frequent variants were European, with 13 cases, while 2 cases corresponded to the African variant; however, no information is available regarding HPV-18 variants and cancer risk in this population [97].
It's worth mentioning that intratype variants of HPV types 31, 35, 62, and 58 have been also characterized from cervical samples of Mexican women from different states of Mexico [97,100,101]. Additionally, variations in HPV-51 have been analyzed in HIV positive Mexican patients, this viral type being highly prevalent in this group [102]. Although no significant associations with the severity of the lesions have yet been found, these studies set precedents to continue analyzing intratype variations of different types of HPV in relation to their distribution and pathogenicity.
The most prevalent variants of HPV-16 and -18 found in the Mexican population are summarized in Table 1.

Functional Analysis of HPV Variants: From Non-Coding Region Activity to Viral Protein Products, an Interesting Approach in HPV Research in Mexico
Until now, many research groups have focused on functional studies of intratype variants, mainly on the most prevalent types found in the Mexican population. These studies focus on elucidating the mechanisms by which these variants participate in cancer progression and the establishment of persistent infections. Assays that assess transcriptional and replication activity in the LCR, as well as early proteins, such as E6, E7, E1, and E2, are essential to understand the outcome of HPV-related cancers. Notably, most of the performed functional studies focus on deciphering the involvement of E6 variants in processes related to cancer, probably because among HPV oncogenes, E6 harbors more nucleotide variations compared to E7 (Supplementary File).

Long Control Region
Epidemiological studies reveal that AA HPV-16 variants are more oncogenic than E variants, and such oncogenicity could be due to E6 and E7 activity, viral replication, or transcription enhancement. It was previously suggested that E2 may be involved in the increased oncogenicity of the AA variants through the regulation of E6/E7 transcription. The presence of intact E1/E2 in cervical tumors was shown to correlate with a higher viral copy number than those with a disruption in the E1/E2 region, where a low viral copy number was observed. Moreover, in tumors positive for HPV-16, it was shown that those containing E2 AA (170 ± 38.9) had a higher number of viral copies per cell than those positive for E2 E (67.6 ± 13.5) [103] (Figure 1a). Furthermore, it was suggested that the E6/E7 overexpression is a late event in cases harboring the E variant, and that E6/E7 up-regulation depends on HPV episome integration into the host genome causing the E2 gene to be disrupted [104]. In contrast, carcinomas harboring AA variants retain the E2 gene, suggesting an alternative mechanism for the regulation of E6/E7 transcription. Interestingly, it was observed that in tumors containing AA variants, the levels of E6/E7 transcripts are 2.1 times higher than those with E variants, and that no differences in the status of the viral genome were found between them (Figure 1b). Moreover, the E2 E variant highly repressed the LCR transcriptional activity compared to E2 AA-a and AA-c, and such effect was unrelated to the capacity of E2 to induce apoptosis [104]. it was demonstrated that E2 binding in the LCR did not differ, and E2 repressor activity was unaffected [106].

E6
One way to characterize the specific involvement of high-risk HPV E6 intratype variants is through bioinformatic analyses of gene expression profiles that may ultimately impact in the modulation of cellular processes and cell signaling pathways associated with carcinogenesis.
Expression patterns were assessed in C-33 A cells stably transfected with HPV-16 E6 In addition to containing E2 binding sites that regulate E6/E7 expression, the LCR also harbors an origin of replication (Ori site), important for the amplification of the viral genome. Another study revealed that variations in HPV-18 E1 and in the LCR differentially affect viral replication. The Ori site of the AsAi variant was found to exhibit a higher replication rate than the Ori site of the E and Af variants, when the E1 and the E2 of AsAi were present (Figure 1c). Furthermore, variations in E2 were shown not to affect Ori function. When evaluating variant homologous elements, that is, when the E1 and Ori sites came from the same variant, E1 Af induced the highest Ori viral replication compared to E1 E and AsAi variants, and such activity was associated with an increased nuclear localization of E1 Af (Figure 1d) [105].
Previously, it was reported that the LCR of the HPV-18 Af variant contains the highest number of nucleotide changes compared to E and AsAi variants, especially in the tissue enhancer region, impacting low transcriptional activity, probably affecting the binding of transcription factors, including KRF-1/OCT1 and OCT-1/TEF1 (Figure 1e). Interestingly, no differences in transcriptional activity were found between LCR AsAi and LCR E variants. Even though several LCR nucleotide changes are found in the HPV-18 variants, it was demonstrated that E2 binding in the LCR did not differ, and E2 repressor activity was unaffected [106].

E6
One way to characterize the specific involvement of high-risk HPV E6 intratype variants is through bioinformatic analyses of gene expression profiles that may ultimately impact in the modulation of cellular processes and cell signaling pathways associated with carcinogenesis.
To determine the genome wide expression modulated by E6 variants of HPV-18 found in the Mexican population, an analysis of expression patterns in HaCaT immortalized cells transiently transfected with HPV-18 E6 AsAi and Af variants was carried out. It was found that E6 Af regulates the expression of 2236 genes, while the E6 AsAi variant affected the expression of 1945 genes, and both variants regulated 414 genes in common, suggesting that the E6 HPV protein regulates cellular transcriptional programs regardless of the variant in question. As expected, both variants were able to alter genes related with cancer pathways, E6 AsAi preferentially altered elements of PI3K, Notch, and ErbB signaling, while E6 Af altered processes, such as cell migration, cell cycle, DNA replication, and mTOR signaling and Wnt signaling pathways. In cell proliferation assays, it was observed that E6 AsAi had a greater increase in proliferation than E6 Af. Conversely, migration was higher in cells harboring E6 Af than those with AsAi [108].
Another study showed that HPV-16 E6 variants E-P, AA-c, E-G350, and E-C188/G350 differentially regulate the expression of long non-coding RNAs, including MALAT1, FLVCR1-AS1, CASC15, ZNF667-AS1, SNHG8, CASC19, and MINCR. Interestingly, it was found that MINCR was overexpressed (8-fold) in C-33 A cells harboring E6 AA-c in comparison with the control group. Through bioinformatic analysis, it was found that MINCR is a competitive endogenous RNA (ceRNA) of miR-28-5p, and that RAP1B transcript is a target of miR-28-5p. Furthermore, using data from TCGA of HPV-16 positive cervical cancer cases, it was demonstrated that RAP1B expression was higher, compared to HPV negative cases ( Figure 2b); however, no differences were found when evaluating miR-28-5p expression. Moreover, C-33 A cells expressing the E6 AA-c variant exhibited a 0.3-fold decrease expression of miR-28-5p compared to the control group, while the expression and the protein levels of RAP1B increased 21.8-and 2.3-fold, respectively. These data demonstrated an inverse correlation of miR-28-5p and RAP1B in E6 AA-c expressing cells. Additionally, when MINCR was knocked-down in these cells, miR-18-5p expression recovered by 90%, while RAP1B expression was reduced by 60%. Finally, when MINCR was ablated, a reduction of 50% in migration was observed in C-33 A cells harboring HPV-16 E6 AA-c variant. These data are consistent with those previously reported that Asian-American HPV-16 variants exhibit a greater oncogenic potential compared to European variants [109].
Using proteomic analysis, it was demonstrated that E-P, E-G350, E-A176/G350, E-C188/G350, AA-a, and AA-c HPV-16 E6 variants differentially regulate proteins mainly involved in metabolism, including KCNH1, PPP2R1A, CCT8, ALDH1A1, ZN554, ENOA, G6PD, ARMX1, GRHPR, CPT1A, CENP, and ZNF93. In addition, qPCR assays showed that E-A176/G350 and AA-a variants repressed ALDHA1 compared to C-33 A cells transfected with the empty vector, while CCT8 expression was diminished in cells containing E6 E-A176/G350, E-P, and AA-a [110]. Currently, the structure of E6 proteins and their variants has not been characterized under native experimental conditions, which is a limitation to study the interactions with cellular proteins involved in the establishment of cancer. A recent study revealed that variations in E6 genes affect the thermodynamics and structural stability of the 3D protein structure. Through Molecular Dynamics (MD) simulations and docking analyses, it was found that E6 variants, including E-P (reference), E-G350 (L83V), E-C188/G350 (E29Q/L83V), E-A176/G350 (D25N/L83V), AA-a (D25N/L83V), and AA-c (Q14H/I27R/H78Y/L83V) interact with a high affinity with the PDZ-1 domain within MAGI-1. Interestingly, it was observed that E6 E-C188/G350, AA-a, and AA-c proteins showed a greater difference in the C-terminal region in comparison to the E-P variant, while E-G350 and E-A176/G350 maintain similar 3D structures to the E-P variant. Moreover, it was found that residues (141-151) in the C-terminal region of E-P are less flexible than the other tested variants, probably affecting the capacity of the PBM motif (ETQV) to interact with PDZ-containing proteins such MAGI-1. It was found that the E6 non-E-P variants required less energy to interact with MAGI-1, having more affinity for MAGI-1 than the E-P variant, particularly E-A176/G350. The residues responsible for the E6-MAGI-1 interaction were found to be close to the changed nucleotides in the variants, in addition to the fact that the interacting residues differ among the variants analyzed. The detailed study of the interactions that occur between E6 and cellular proteins will allow us to understand the differential processes modulated by E6 variants during the replicative cycle and the establishment of cancer [111]. MCF-7 and C-33 A cells stably expressing AsAi, E, and Af variants were used to investigate the specific contribution of HPV-18 E6 variants. E6 splicing patterns were shown to be different among the tested variants. Cells expressing E6 AsAi and E6-E exhibited high levels of full-length E6 transcript with low levels of E6*I transcript. In contrast, E6 Af produced high levels of E6*I with low amounts of full-length E6 (Figure 2c). Interestingly, these expression patterns were similar in tumors containing HPV-18 variants [18]. The E6 mRNA profile of the tested variants was in concordance with E6 protein patterns. The reduction of functional full-length E6 may partly explain the biological behavior observed in HPV-18 premalignant lesions and cancer. Interestingly, p53 levels decreased in the presence of E6 AsAi and E, while in cells harboring E6 Af, p53 levels were unchanged (Figure 2c). An in vitro p53 degradation assay was performed to evaluate the specific effect of variations in full-length E6, where no splicing should occur; surprisingly, no differences were found between E6 AsAi and Af variants in their ability to induce p53 degradation, confirming that the reduced proportion of full-length E6 in the African variant impacts the reduction of p53 degradation. Additionally, biological assays revealed that cells with AsAi and E variants were able to induce colony formation and tumor growth more quickly than the E6 Af variant. Interestingly, when the E6 Af variant was mutated (C491A), a change in E6 expression pattern was obtained that was similar to that observed for the E6 AsAi variant, and even the capacity to induce p53 degradation was restored. Taken together, these data suggest that differences in oncogenic capacity displayed by HPV-18 E6 variants are highly dependent on their splicing patterns [18].
In agreement with these data, using MCF-7 cells stably transfected with HPV-18 E6 AsAi, E, and Af variants, differences in E6/E6*I ratio among variants were confirmed, showing that E6 Af produced a higher proportion of the E6*I transcript. Subsequently, p53 protein levels and its downstream effector p14ARF were evaluated in the presence of E6 variants. As expected, p53 protein was ablated in E6 AsAi and E expressing cells, whereas a low reduction of p53 levels was observed in cells with E6 Af. In contrast, p14ARF, a gene negatively regulated by p53, was augmented in the presence of E6 AsAi and E variants and decreased in cells with E6 Af (Figure 2c). Those effects were attributed to a higher proportion of E6*I and low full-length E6 levels in cells with E6 Af, where E6 is limited in its ability to degrade p53 [112].
HPV-18 E6 variants were shown to differentially regulate the activation of cell signaling pathways involved in proliferation, such as Wnt/β-catenin and PI3K/Akt. In C-33 A cells transiently transfected with the E6 AsAi variant, a 1.5-fold increase in TCF-4 transcriptional activity was observed compared to those cells transfected with the empty vector. Interestingly, the E6 Af variant induced a significant increase (2.5-fold) in TCF-4 transcriptional activity above that seen with the E6 AsAi variant (Figure 2d). Such activity was enhanced in E6 AsAi and Af expressing cells when β-catenin was exogenously overexpressed (1.6-and 2.25-fold above control, respectively). Cyclin D1 and Axin2 genes were similarly up-regulated in the presence of E6 AsAi and Af variants. E6 variants were found to interact with β-catenin and TCF-4. Additionally, TCF-4 stability and nuclear levels were increased in the presence of E6 AsAi; and, surprisingly, E6 AsAi was also found to bind to the SP5 promoter, which is transcriptionally regulated by TCF-4. Finally, the contribution of E6 AsAi in the regulation of the Wnt/β-catenin pathway was to increase cell proliferation [113].
The effect of HPV-18 E6 variations on the PI3K/Akt pathway was evaluated. Protein levels of PTEN, an inhibitor of PI3K kinase, were found to remain unchanged in MCF-7 and C-33 A cells expressing E6 AsAi, E6 E variants, or without E6. However, when phospho-PTEN was assessed, virtually no phospho-PTEN was found in cells with E6 AsAi or E6 E, while reduced phospho-PTEN levels were found in cells expressing E6 Af compared to cells with no E6 (Figure 2e). Additionally, hDlg, a negative PTEN regulator, was reduced in the presence of E6 variants, although E6 Af produced a smaller reduction of hDlg levels relative to E and AsAi variants. In addition, the PTEN/hDlg interaction was decreased in cells expressing the Af variant. A similar effect occurred in the levels of phospho-Akt/PKB, phospho-PI3K, and phospho-ERK1/2, since the presence of E6 AsAi and E variants increased those levels, suggesting a strong activation of the pathway, while the Af variant induced a modest activation of Akt/PKB and no activation of ERK1/2 or PI3K. Consistent with these results, cells expressing E6 AsAi and E6-E increased cell proliferation, while the lowest proliferation rate was observed for cells expressing E6 Af. As described above, the splicing patterns in E6 transcripts differ among E6 HPV-18 variants. Therefore, those differences in phospho-PTEN levels observed among AsAi and E variants compared to Af, could be attributable to the fact that E6 AsAi and E retained the full-length E6 isoform, while E6 Af produced practically only the E6*I transcript [114]. It is likely that the effect observed on cell proliferation by each of the variants is partially due to an effect exerted by changes in the proportions of E6 isoforms. In this sense, different proportions of E6 and E6* isoforms are found depending on various cellular stimuli. It has been reported that both E6 and E6* exhibit controversial functions, for example, E6 induces p53 degradation, while E6* interferes with E6-dependent p53 degradation. Furthermore, unlike E6, E6* does not induce keratinocyte immortalization, but instead reduces cell proliferation and the tumor size of HPV-16 positive cancer cells. Additionally, E6 and E6*I have been shown to increase the levels and the activity of caspase 8, in addition to inducing its nuclear translocation, however, they do not promote apoptosis [115].  (c) HPV-18 E6 AsAi and E variants exhibit different E6 splicing patterns, being E6 full-length the most prevalent isoform, in contrast E6 Af showed a higher proportion of the small-spliced isoform E6*I, impacting p53 and p14ARF protein levels. (d) HPV-18 E6 Af variant induces higher TCF-4 transcriptional activity than E6 AsAi. Both variants bind to TCF-4/β-catenin and induce Cyclin D1 and Axin2 expression. (e) HPV-18 E6 AsAi and E variants promote overactivation of the PI3K/AKT pathway compared to E6 Af, differentially impacting proliferation. When indicated, up arrows mean augmented levels, where enlarged arrows refer to higher levels than the short arrows. Down arrows indicate decreased levels. Figure created in BioRender.

E7
To elucidate the specific participation of the intratype variants of E7 in cancer processes, one in silico study has been carried out to establish differences in the three-dimensional structure of HPV-16 E7 N29S and H51N variants in comparison with the E-prototype. Those variants have been most frequently associated with cervical cancer. The amino acid changes within variants are located in CR2 and CR3 domains for N29S and H51N, respectively. Interestingly, it was reported that the change of N29S is susceptible for phosphorylation that may reinforce nearby sites (S31 and S32) also susceptible to phosphorylation which could affect the structure of E7. Molecular dynamics (MD) analysis revealed that the N29S variant exhibits a more compact structure than the E7 E-P and H51N variants. Furthermore, 3D structure was evaluated comparing E7 variants, revealing major changes in the CR2 and CR3 domains that may affect the interaction with different cellular proteins, including pRb. Further, when analyzing the 2D structure by MD, it was demonstrated that E7 E-P lost its α-helix and β-sheets, while the H51N variant conserved the α-helix but lost the β-sheets, while most of the secondary structures of N29S were preserved. These results indicate that N29S variant is more structurally stable than the reference E7 and H51N. Finally, using analysis of electrostatic potential, differences in the charge on the LXCXE motif were also observed. The E7 E-P variant essentially presents a negative charge in this region, while the N29S variant harbors a combination of negative, electrostatically positive, and neutral charges in this region. Meanwhile, the H51N variant contains negative, neutral, and positive regions in the LXCXE motif. These interesting findings suggest that the affinity of cellular proteins for E7, which is dependent on the LXCXE motif may differ between variants, which may affect their oncogenic potential [116].

E2
In relation to the participation of E2 variants in cancer, hitherto only one functional study has been carried out. In order to have a more comprehensive idea of the participation of E2 AsAi and-Af variants of HPV-18 in processes associated with cancer, the transcriptional profiles altered by E2 and its variants were assessed in different cell lines, including MCF-7, C-33 A, and HEK293. The analyses included the results of the three cell lines in order to evaluate the effect of the E2 variants regardless of the cellular context, that is, the aim was to show the genes that were exclusively deregulated by the E2 variants regardless of the cell type. It was demonstrated that E2 variants differentially regulated gene transcription patterns where expression of 1361 genes were affected by E2 AsAi, while E2 Af altered the expression of 235 genes, and only 34 genes were shared between the E2 tested variants. Interestingly, gene ontology analysis demonstrated that both AsAi and Af variants regulated cellular processes, including metabolism (523 vs. 73 genes), proliferation (78 vs. 15 genes), cellular communication (252 vs. 59 genes), and death (93 vs. 18 genes), respectively. It was demonstrated that the effect of E2 variants on gene regulation was maintained in different cell contexts. Moreover, it was previously demonstrated that E2 induced cell death; however, when comparing the effects of E2 variants it was found that E2 AsAi favors gene expression related to the extrinsic apoptotic pathway, while E2 Af favors the intrinsic apoptotic pathway. Through viability assays, it was shown that both variants reduced cell viability in comparison to control cells (77.63% and 46.27% for E2 AsAi and E2 Af, respectively); however, the strongest effect was observed in those cells with E2 Af. Similar results were obtained when apoptosis was evaluated, it was found that both variants differentially induced apoptosis in relation to the control cells, being 29.4% attributable to E2 AsAi, while E2 Af induced 53.3%. These results suggest that E2 variants differentially regulate cellular viability and apoptosis, probably due to the different gene expression profiles induced by the variants [117].

Concluding Remarks: Where We Are and Where We Are Going?
Cervical cancer continues to be a serious public health problem in Mexico, since mortality related to this type of cancer has had a significant increase from 3357 cases reported in 2012 to 4335 cases in 2020 [118]. The vaccination program in Mexico was introduced in 2012 for 11-year-old girls, who have not started their sexual life, and it includes two doses of the quadrivalent vaccine, where the second dose is applied at six months up to one year after the first dose. According to official reports provided by the World Health Organization portal, a coverage of more than 80% of vaccination against papillomavirus has been achieved in Mexico in this population, during the period 2012-2019. Notably, these numbers decreased drastically in 2020 to 11.22% and, unfortunately, 0.45% was reported for 2021, probably due to the health policy in our country being mainly focused on dealing with the COVID pandemic [119]. To our knowledge, there are so far no studies available that evaluate vaccination coverage in addition to highlighting the efficacy of papillomavirus vaccines in Mexico. Unfortunately, the incidence and the mortality of cervical cancer and other HPV-related cancers continue to rise. In addition, the effect of immunization against papillomavirus granted by vaccination programs in the prevention of cervical cancer could be reflected in the coming years (more than 20 years), when the first vaccinated population reaches the age at which the greatest number of cases of cervical cancer currently occur. Likewise, contemplating the beginning of the vaccination scheme in young men will allow for controlling the spread of papillomavirus infection, in addition to encompassing the prevention of the development of other HPV-related tumors, including penile, anus, and oropharynx cancers [120][121][122]. The participation of various research groups throughout Mexico, has provided cutting-edge knowledge in the study of papillomavirus and cancer, either through epidemiological or functional studies of HR-HPV and its variants, which allow for defining the geographical distribution of HPVs and understanding the molecular and the cellular events affected by these viruses. Interestingly, as highlighted in the previous sections, HPV-16 is the most prevalent viral type cancer in the Mexican population, and even more interesting is that some intratype variants are commonly found in cancer samples that confer risk of progression to develop cancer. This suggests that variants of the same viral type differentially and finely regulate cellular processes that ultimately lead to cancer. Furthermore, it is of interest to observe the effect of vaccines on the protection against tumors that contain the same HPV type but have a different variant. Understanding the biological behavior of intratype variants of papillomaviruses associated with cancer will eventually allow for the development of targeted strategies aimed not only at preventing cervical cancer but also at predicting the outcome of patients currently suffering from this type of cancer.

Conflicts of Interest:
The authors declare no conflict of interest.