This case–control study included patients diagnosed with pemphigus vulgaris (PV) at the Bullous Diseases outpatient clinic, Department of Dermatology and Venereology, Karadeniz Technical University (KTU). Diagnosis was confirmed by clinical findings, histopathology, direct immunofluorescence, and/or ELISA. Only patients without other autoimmune diseases who provided informed consent were included. In total, 86 PV patients (aged 21–81 years) and 200 age- and sex-matched healthy controls participated. Peripheral blood samples (5–9 ml) were collected in EDTA tubes for genomic DNA extraction using the MiniFavorPrep™ kit (Favorgen, Australia). The study was approved by the KTU Faculty of Medicine Ethics Committee (protocol no: 2022/51).

HLA-DRB1 typing was performed by sequence-specific oligonucleotide probe (PCR-SSOP) using the Luminex platform (LABScan3D, One Lambda, USA) and analyzed with HLA Fusion software. The One Lambda SSP^® kit (Thermo Fisher Scientific) was used for allele identification (16, 17). Power analysis was conducted using the PGA v2.0 program to ensure 80% statistical power software (18). Hardy–Weinberg equilibrium in controls was tested with GENEPOP (19). Allele and genotype frequencies were compared using Pearson’s chi-square or Fisher’s exact test when expected counts were <5, implemented in RStudio (v2024.04.2). Significance was adjusted with the Benjamini–Hochberg FDR method; p_adj < 0.05 was considered significant.

To identify candidate cadherin-like proteins with potential for molecular mimicry, we performed BLASTP searches against the VectorBase database (https://vectorbase.org/vectorbase/app) using the EC1 domain of human DSG3 as the query sequence. We employed PSI-BLAST to increase sensitivity for detecting structurally conserved and evolutionarily related sequences. PSI-BLAST utilize position-specific scoring matrices to capture deeper pattern-level homology not evident from linear sequence identity alone. The top-scoring result was a cadherin-like protein from Aedes albopictus, annotated based on homology to Aedes aegypti (protein: AAEL000597; Gene symbol: LOC5564848/neural-cadherin) (Supplementary File). This protein, with a predicted length of 1569 amino acids, was selected for further analysis based on its high similarity score and cadherin domain architecture. From the BLAST alignment, we extracted the homologous region to human DSG3 ectodomain 1 (hDSG3 EC1) and converted it to FASTA format. This sequence was subjected to MHC class II peptide-binding predictions using the IEDB MHC-II binding tool (www.iedb.org), applying the same parameters used for hDSG3. Based on our and other genetic association findings, HLA-DRB1*04:02 and HLA-DRB1*14:01 were selected as the MHC class II alleles of interest. Candidate antigenic peptides were prioritized by (1): strong predicted binding affinity (2), high immunogenicity scores, and (3) maximal alignment similarity with hDSG3 EC1 sequences (20, 21).

Structural modeling of EC1 peptides from hDSG3 and Aedes albopictus cadherin-like protein was performed using AlphaFold3 and AlphaFoldMultimer tools (22, 23) on the Colab platform (24). Molecular models were visualized using UCSF ChimeraX (v1.10) (25). All-atom molecular dynamics (MD) simulations were carried out using GROMACS (v2022.5) with the AMBER99SB-ILDN force field. Systems were solvated with TIP3P water, neutralized, energy-minimized, equilibrated, and subjected to a 100 ns production run at 298 K and 1 atm using the V-rescale thermostat and Parrinello–Rahman barostat (26). MD trajectories were analyzed to characterize structural stability and residue-level flexibility of the peptide–HLA complexes. Peptide backbone root mean square deviation (RMSD) was calculated following alignment of each trajectory frame to the HLA backbone atoms. Root mean square fluctuation (RMSF) was computed for peptide residues using the equilibrated portion of the trajectories. For all analyses, the final 47.5 ns of the 100 ns production simulations were considered. Detailed trajectory analysis procedures, including atom selection criteria and alignment settings, are provided in the Supplementary Methods.

Binding free energies (ΔG_binding) for peptide–protein complexes were estimated using gmx_MMPBSA (v1.5.7) with the MM/GBSA approach (GB-OBC II, igb=5). Representative frames from the equilibrated trajectories were used for binding energy estimation and per-residue energy decomposition analyses (27, 28). All additional parameters, scripts, and analysis settings are detailed in the Supplementary File.

In the control group, HLA-DRB1 alleles were in Hardy–Weinberg equilibrium. A total of 33 distinct HLA-DRB1 alleles were detected across all participants (Table 1). The ten most frequent alleles in PV patients and controls are listed in Table 1. The most common alleles observed in PV patients were HLA-DRB1*04:02 (40%) and HLA-DRB1*14:01 (18%). In controls, the predominant alleles were HLA-DRB1*11:01, 11:04, and 07:01 (each approximately 9%). Allele frequencies and corresponding statistical analyses are summarized in Table 2. Carriage of HLA-DRB1*04:02 was significantly associated with PV in this cohort (OR = 12.95; 95% CI: 7.40–24.57). Similarly, HLA-DRB1*14:01 was more frequent among PV patients compared with controls (p = 1.23×10^-5; OR = 4.55; 95% CI: 2.35–8.74). Genotype analysis was subsequently performed to evaluate the co-occurrence of both alleles (HLA-DRB1*04:02/14:01). Genotype frequencies and statistical results are shown in Table 3. The heterozygous HLA-DRB1*04:02/14:01 genotype was observed in PV patients but not detected among controls. This distribution reflects an enrichment of this genotype in the patient group within the studied cohort; however, estimates based on low genotype counts should be interpreted with caution.

MHC class II peptide-binding predictions performed using the IEDB tool identified 19 unique predicted peptides within the hDSG3 EC1 domain and 167 within the A. albopictus EC1-like domain (Aa_P) (Supplementary Table S1). From these predicted peptides, candidates were prioritized based on a combination of predicted binding affinity, immunogenicity score, and sequence similarity to the hDSG3 EC1 domain (Figure 1).

Structural comparison of the modeled hDSG3 and A. albopictus EC1 domains revealed a high degree of similarity in overall domain architecture and secondary structure (Figures 2A-C; Supplementary File). Following docking of the predicted peptides to HLA-DRB1*04:02, stable receptor–ligand complexes were obtained for both peptide–HLA pairs (Figure 3). Subsequent molecular dynamics simulations indicated comparable binding poses for the two peptides throughout the analyzed trajectories. Despite these structural similarities, calculated binding free energies differed between the two complexes (Table 4).

The dynamic stability of the peptide–HLA complexes was evaluated by peptide backbone root mean square deviation (RMSD) analysis following alignment to the HLA backbone (Figure 4A). RMSD traces calculated over the equilibrated portion of the trajectories (last 47.5 ns) indicated that both the hDSG3 EC1-derived peptide and the Aa_P peptide fluctuated within comparable RMSD ranges, without evidence of sustained structural drift. Occasional transient RMSD excursions were observed for the hDSG3 EC1-derived peptide during the analyzed window; however, these fluctuations were short-lived and did not result in a sustained increase in RMSD, supporting comparable backbone stability of both peptides under the simulated conditions. Residue-level flexibility of the bound peptides was assessed by root mean square fluctuation (RMSF) analysis after alignment to the HLA backbone (Figure 4B). Both the Aa_P and hDSG3 EC1-derived peptides exhibited low RMSF values across the central region of the peptide sequence, indicating constrained fluctuations within the HLA binding groove. Increased flexibility was primarily observed at the N- and C-terminal residues for both peptides, a pattern commonly associated with reduced anchoring at peptide termini. Overall, the RMSF profiles of the two peptides were broadly similar, with only modest residue-specific differences. No pronounced or regionally distinct increases in flexibility were observed for either peptide, suggesting comparable dynamic accommodation within the HLA binding site under the simulated conditions.

To assess the plausibility of the VCM hypothesis, we constructed HLA-DRA/DRB1*04:02–peptide complexes and performed peptide-binding predictions, structural modeling, and molecular dynamics simulations. Binding free energy (ΔG_binding) calculations indicated that the Aa_P peptide derived from A. albopictus exhibited a more favorable binding energy than the hDSG3 EC1-derived peptide (Table 4), suggesting a high affinity for the risk-associated HLA-DRB1*04:02 molecule. Structural similarity and comparable binding modes support the potential for molecular mimicry (Figures 2, 3). These computational results provide preliminary support for the VCM hypothesis and guided the development of a mechanistic model (Figure 5).

In genetically susceptible individuals (e.g., carriers of HLA-DRB1*04:02 or *14:01), repeated exposure to cadherin-like peptide candidates potentially present in the saliva of Aedes albopictus or other hematophagous insects (e.g., sandflies) may be associated with a subclinical or priming immune response against these foreign antigens. One plausible scenario involves immune priming during vector bites, in which tissue injury or microbial co-inoculation generates damage-associated molecular patterns (DAMPs) or pathogen-associated molecular patterns (PAMPs). These signals may promote dendritic cell activation and favor the differentiation of naïve T cells toward immunogenic rather than tolerogenic phenotypes. In the presence of sequence and structural similarity to hDSG3 and hDSG1, certain B- and T-cell clones may therefore exhibit features consistent with molecular mimicry and potential cross-reactivity. Over time, somatic hypermutation, affinity maturation, and class-switch recombination within germinal centers (44) may enhance the specificity and binding strength of these B-cell receptors, potentially increasing the likelihood of recognition of self-antigens (DSG3/DSG1) (45). Initial tissue damage may further expose cryptic desmoglein epitopes, which are processed and presented by antigen-presenting cells, thereby promoting epitope spreading. This process could drive clonal expansion of autoreactive B cells and facilitate a transition from vector-derived immune recognition to sustained autoimmunity. Ultimately, persistent autoantibody production against DSG3 and DSG1 could contribute to disruption of keratinocyte adhesion, consistent with mechanisms underlying blister formation in pemphigus vulgaris (Figure 5).

To test our VCM hypothesis, at least in silico, we employed a combination of generating HLA-DRA/DRB1*04:02 complexes, peptide-binding predictions, structural modeling, and molecular dynamics analyses (Table 4). Our analysis revealed that the Aa_P peptide exhibited lower estimated binding free energy than the hDSG3 EC1-derived peptide in our in silico analyses. These data provide preliminary support for the VCM hypothesis, although experimental validation is required. Additional contextual support for the plausibility of the hypothesis comes from transcriptomic data reported by Sim and colleagues, indicating regulation of cadherin-related transcripts in the salivary glands of Aedes aegypti (46) (Supplementary Table S3). For experimental validation, an appropriate animal model is essential. Humanized mice expressing the HLA-DRB1*04:02 allele can be utilized, as previously demonstrated by Eming et al. in a preclinical mouse model of PV (47). The VCM hypothesis could be conceptually evaluated in future studies using these mice with both recombinant hDSG3 ectodomains and A. albopictus–derived cadherin-like (CAD-like) peptides (Figures 1A, 2B).

In regions such as Brazil and Tunisia, endemic pemphigus foliaceus (EPF)—including its localized form Fogo Selvagem (FS)—has long been suspected to arise from environmental exposures, particularly to insect vectors. Multiple studies have identified associations between anti-DSG1 autoantibodies and exposure to vector-borne antigens. For example, in Brazil, Lutzomyia longipalpis (the sandfly vector of visceral leishmaniasis) was shown to elicit anti-Dsg1 IgG4 antibodies through salivary proteins such as LJM11 and LJM17, which cross-react with epitopes on human Dsg1. In murine models, immunization with LJM11 induced anti-Dsg1 antibodies, lending support to the hypothesis of vector-induced molecular mimicry in genetically susceptible individuals (e.g., HLA-DRB1*04:04 carriers) (12–15). Similarly, in Tunisia, anti-Dsg1 antibodies have been detected in sera of patients with leishmaniasis and hydatidosis (48). Salivary proteins such as PpSP32 from Phlebotomus papatasi were also shown to bind both Dsg1 and Dsg3 in vitro. Interestingly, although vector antigens triggered desmoglein-recognizing antibodies in mice, definitive cross-reactivity and pathogenic potential remain unconfirmed. These observations suggest a possible mimicry mechanism but fall short of providing direct evidence for pathogenic antibody induction and acantholysis (13, 14, 48, 49). While these studies contribute valuable epidemiological and immunological insights, they do not fully explain the molecular mechanisms of epitope recognition or autoimmune activation. Notably, sandfly-derived proteins studied thus far have not been conclusively shown to generate immunodominant or acantholytic responses, and structural mimicry has not been conclusively demonstrated through high-resolution computational or crystallographic methods (48, 50). Against this backdrop, our current study advances the field by focusing not on salivary proteins per se, but rather on cadherin-like proteins encoded by the vectors themselves, particularly those with structural homology to the EC1 domain of human Dsg3—the primary autoantigen in PV. This novel approach addresses a major gap in the literature by proposing that vector-derived cadherin homologs may act as molecular mimics of hDsg3. Unlike earlier research, our study integrates peptide–HLA binding predictions, homology modeling, and molecular dynamics simulations to explore the plausibility of cross-reactive interactions. By centering on HLA-DRB1*04:02, a well-established susceptibility allele for PV, we propose a mechanistically and structurally coherent framework for how vector-derived proteins may break immune tolerance and initiate autoimmunity in genetically predisposed individuals.

Several limitations of the present study should be acknowledged. First, although statistically significant associations were observed for HLA-DRB1*04:02 and HLA-DRB1*14:01 alleles, genotype-level analyses—particularly for the heterozygous HLA-DRB1*04:02/14:01 genotype—were based on low-frequency observations. As a consequence, estimates of effect size for rare genotypes are subject to statistical uncertainty and imprecision, and confidence intervals should be interpreted with caution. Under conditions of sparse data, odds ratio estimates may be inflated and do not necessarily reflect stable population-level effects.

Second, the study cohort size, while sufficient to detect allele-level associations, limits the robustness of genotype-level inference. Confirmation of the observed genotype-level enrichment will therefore require replication in larger, independent cohorts and in populations with different genetic backgrounds.

Third, the molecular dynamics simulations and binding free energy analyses were conducted over a finite simulation timescale and were designed for comparative, hypothesis-generating purposes. These in silico analyses do not provide direct evidence of pathogenicity, immunogenicity, or causal involvement in disease initiation, and should be interpreted as supportive structural and energetic observations rather than definitive mechanistic proof.

Finally, the proposed Vector-derived Cadherin Mimicry (VCM) framework is conceptual in nature. While it integrates immunogenetic associations with structural modeling, experimental validation—using appropriate cellular or animal models—will be necessary to assess the biological relevance of this hypothesis and its potential contribution to pemphigus vulgaris pathogenesis.

In conclusion, while vector-associated molecular mimicry has been most extensively explored in the context of endemic pemphigus foliaceus, the present study proposes a computationally informed and structurally plausible extension of this conceptual framework to PV, with a shift in focus from DSG1 to DSG3. Rather than establishing a definitive mechanistic link, this hypothesis expands existing models of vector-associated autoimmunity by integrating immunogenetic context with structural bioinformatics. To the best of our knowledge, this study also reports a significant enrichment of the heterozygous HLA-DRB1*04:02/14:01 genotype in a PV cohort, although confirmation of genotype-level effects will require validation in larger, independent populations. In parallel, our analyses suggest a previously underexplored connection between genetic susceptibility and repeated environmental exposure to vector-derived cadherin-like proteins. Through the integration of immunogenetic data, structural modeling, and comparative energetic analyses, we provide preliminary in silico observations consistent with the proposed VCM framework. If supported by future experimental studies, this model could contribute to a broader understanding of autoimmune pathogenesis, in which hematophagous insect vectors may represent context-dependent environmental factors rather than direct causal agents. Such insights may ultimately inform HLA-based risk stratification approaches and guide the development of targeted preventive strategies in genetically susceptible populations.