A Darwinian view of Behçet's disease

Introduction

The sine qua non of Behçet’s disease (BD) is recurrent oral aphthous ulceration, which must be accompanied by two or more of the following major criteria: genital ulceration, uveitis, and skin lesions such as erythema nodosum and a positive pathergy test. A multitude of other clinical features, although not essential for diagnosis, is now recognized to be associated with the disease, such as seronegative arthropathy, increased risk of venous thrombosis, and gastrointestinal inflammation. The cause of BD is unknown, but an autoinflammatory etiology, following infection by an as yet unidentified pathogen, has been suggested. In this scenario, a bacterial or viral infection of mucosa would induce an immune response, leading to cross-reactivity with self-proteins and producing the organ-specific damage, characteristic of BD.^[1] Although the precise etiology of BD remains mysterious, there is a strong genetic component to the disease, with a range of well-known associations with particular alleles.^[2]

HLA-B and chromosome 6

A heritable risk factor for BD was first identified in 1982 when an association with human leukocyte antigen-B5 (HLA-B5) in the human major histocompatibility complex (MHC) was first reported.^[3] The HLA-B*5 locus includes the family of HLA-B51 alleles and HLA-B52. In the majority of populations, the greatest heritable risk factor for disease, and indeed for ocular disease severity, is HLA B51. In contrast, HLA-B*52, which differs from HLA-B*51 by only two amino acids in the peptide-binding groove, is not associated with BD in any population.^[4] The geographical distribution of BD, concentrated mainly in countries from the Mediterranean basin to the Far East, led to the synonym of “Silk Road” disease. These regions lie between latitudes 30° and 45° north, and it is of interest that the reported prevalence of BD is associated with latitude and rarely encountered in more northern climes, the Americas, and Australasia. Furthermore, the geographical distribution of HLA-B*51 among healthy subjects roughly corresponds with global disease distribution.^[4] This association between geographical disease distribution and HLA-B*51 led to the hypothesis that genetic risk factors were propagated by migrant traders operating along the Silk Road some 2000 years ago.^[3] However, the new discipline of paleogenetics, involving recent advances in techniques for recovering and analyzing ancient DNA (aDNA), has produced evidence for much earlier acquisition of some HLA alleles, during Paleolithic migrations of modern humans (Homo sapiens) and admixture with archaic human species.

Following the draft sequence of the Neanderthal genome, it was estimated that the Neanderthal contribution represented between 1% and 4% of modern Eurasian genomes.^[5] Recent estimates of confidently inferred Neanderthal ancestry yield means of 1.38% and 1.15% in East Asian and European populations, respectively, but <0.1% of the genome in sub-Saharan African populations.^[6] Long stretches of modern Eurasian genomes lacking Neanderthal ancestry suggest that the genetic contribution was initially much larger, estimated at >3%, but that negative selection has operated to reduce the Neanderthal component now present in modern human genomes.^[6]

Evidence of another archaic population, genetically distinct from Neanderthals, comes from deoxyribonucleic acid (DNA) extracted from a finger bone found in the Denisova Cave in southern Siberia.^{[7, 8]} It is not clear how this population relates to other known hominin species, so the genetically identified population is currently known simply as the “Denisovans”. While present-day Melanesian genomes are estimated to contain a Denisovan contribution of 4–6%, the Denisovan contribution to European genomes appears to be negligible. The Denisovan contribution to the genomes of living East Asians and Native Americans is higher, at around 0.2%, but this may be due to admixture with the ancestors of the modern-day inhabitants of Oceania, long after interbreeding with the Denisovans themselves. There is now genetic evidence for three to five instances of interbreeding between hominin populations, including modern humans and Neanderthals, although this is accepted to provide a low-resolution reconstruction of demographic history.^[9]

Regardless of the imprecision surrounding the events that led to the admixture and incorporation of archaic DNA into modern human genomes, it is now clear that modern humans have received input from archaic genomes and these contributions have variously been both beneficial and deleterious within modern human genomes. Focusing on HLA genes, it seems likely that modern humans acquired the HLA-B*73 allele through admixture with Denisovans. HLA-B*51 was identified in Neanderthal samples and is found in high frequencies in modern human populations in Eurasia but not in Africa. Simulations of HLA gene introgression strongly supported a Neanderthal origin for HLA-B*51 in modern genomes, among others.^[10]

The evidence suggests that Neanderthal admixture contributed HLA-B*07, B*51, C*0702, and C*1602 to modern human genomes. Of these alleles, HLA-B*51 and HLA-C*1602 are associated with BD. Importantly, the frequency of these alleles in modern humans matches the prevalence of BD very closely^[10] [Figure 1]. The existence of these HLA alleles in archaic genomes and their preservation in modern genomes, when other Neanderthal genes were subject to negative selection, strongly suggest that these alleles were beneficial in the context of local pathogens encountered by both Neanderthal populations and modern human populations living in similar environments.

Figure 1

BD, HLA-B*51, and archaic DNA. The frequency of HLA-B*51 is highest in the Middle East (red—high and green—low) and decreases across the Silk Road (dotted line) to China. The prevalence of BD follows a similar pattern (number/100,000 in boxes). Circles show the sites of bones used in typing archaic genomes. BD, Behçet’s disease; HLA-B51, human leukocyte antigen-B51.

Indeed, Neanderthals and their ancestors had been living in Europe and Western Asia for hundreds of thousands of years before modern humans, expanding out of Africa or the Middle East, arrived in that territory. It seems likely that the immune system of Neanderthals had adapted to the local pathogen pool and that this genetic advantage was passed on to modern human populations through interbreeding. The links between admixture, retention of archaic DNA, and disease in human populations, archaic and modern, are likely to be related to genes that were and are protective against infectious diseases.^[10] Piga and Mathieu^[11] had proposed such a link with Neanderthal admixture and negative selection of HLA-B*5101 and Plasmodium falciparum infection and a positive selection in response to Yersinia pestis to describe the geoepidemiology distribution of BD.

The identification of single nucleotide polymorphisms (SNPs), of Neanderthal origin, that are the risk factors for a range of diseases including systemic lupus erythematosus, biliary cirrhosis, Crohn’s disease, smoking behavior, and type 2 diabetes (T2D) further supports a model of the introgression of Neanderthal genes, exerting a positive influence on the survival of early modern humans but leading to increased susceptibility to disease in modern individuals.^[6] The current Sars-CoC2 pandemic has supported a role for Neanderthal-derived gene mutation and disease severity. Increased levels of a Neanderthal haplotype of 2′5′ oligoadenylate synthase (OAS1) in Europeans have a protective effect in Coronavirus-19 (COVID) infection.^[12] By comparison, a Neanderthal core haplotype of a gene cluster on Chr 3 is positively selected in South Asians, particularly Bangladeshis, is associated with a more severe response to COVID.^[13]

The link between archaic alleles and disease in modern humans is paradoxical; if such alleles were retained due to their protective nature, why are they linked with the disease in current populations? Several possibilities exist: first, until very recently, life expectancy was very short, and in particular, infant mortality very high. Genetic polymorphisms that helped to control pathogens in childhood and adolescence, allowing adulthood to be reached, would be retained—regardless of whether those polymorphisms caused “damage” later in life, a concept known as antagonistic pleiotropy.^[14] A second possibility is that genes that were protective in the Paleolithic have a very different effect today, where diet and lifestyle are much changed from those of our ancestors. For instance, modern diseases such as Crohn’s and T2D have been linked to changes in diet including increased sugar intake, which leads to changes in the gut microbiome that play an important role in changing profiles of disease.^[15] Finally, combinations of gene variants maintained by the mechanisms described above are responsible for many modern diseases. This Darwinian concept of many small, mainly positive, changes accumulating over time is supported by the large number of SNPs identified in many complex diseases.

HLA-B*51 and Behcet’s disease

While an obvious interpretation of the role of HLA-B*51 in BD is the presentation of antigen to CD8^+ve T cells, the association with HLA-B*51 could be due to its interaction with killer cell immunoglobulin-like receptor (KIR) genes, which are also highly polymorphic and diversely expressed on natural killer (NK) cells, which for inhibitory receptors is under genetic control, whereas activatory receptors are under environmental control. NK cell function is therefore strictly controlled by the expression of the inhibitory receptors. Alterations in either allele expression of inhibitory KIR or environmental-induced changes in activatory KIR could influence disease prevalence.^[16] While HLA-B loci have been identified in multiple genetic studies as the strongest association, recent deep sequencing of the HLA loci in BD patients with imputation studies based on the 1000 Genomes Project was performed in two different independent BD cohorts. The robust association with HLA-B^*51 was explained by a variant between B*51 and MHC class 1-related A (MICA; rs116799036), which was independent of B*51. Additional independent associations in psoriasis susceptibility locus (PSORS1), HLA-f-AS1, and HLA-Cw*1602 were also identified and replicated.^[17] However, the independence of rs116799036 was not supported in the analysis of a larger cohort of patients that strongly backed HLA-B as the primary gene association on chromosome 6. Importantly, the amino acids relevant to the HLA-B*51 association were identified and included HLA-B*67, which is one of the two amino acid differences between B*51 and B*52. The latter is protective in Middle Eastern patients with BD.^[18] In support of the concept that HLA-B has been the major risk for BD, imputation studies identified SNP in endoplasmic reticulum-expressed amino-peptidase (ERAP1) (rs17482078). There was evidence of epistasis between ERAP1 and HLA-B*51. ERAP1 is a protein involved in trimming peptides to fit MHC class I molecules; mutations conferred risk for BD in HLA-B*51^+ve individuals, with homozygosity of ERAP1 rs17482070 giving an OR of 3.78.^[19] A haplotype, based on ERAP-1 sequencing, designated Hap10 with five non-ancestral amino acids was reported to be associated with BD.^[20] Individual variants encoded by this haplotype rs2287987, rs10050860, and rs17482078 were previously reported by studies in various ethnic groups to be associated with BD.^{[21, 22]} Such variants influence the binding of low-affinity peptidomes to HLA-B*5101, which may be relevant in BD. Specifically, 2 subpeptidomes were identified, whereby only peptides in the subpeptidome with Ala at position 2 are extensively destroyed, except when their position 1 residues are ERAP-1 resistant.^{[23, 24]} Recently, KIR3DL1 allele-level analysis revealed the KIR3DL1^LOW/KIR3DS1 functional genotype to be implicated in disease pathogenesis and KIR3DL1^HIGH/KIR3DL1^NULL to be protective.^[25] These data suggest a complex mechanism by which HLA-B*5101 is involved in BD and more functional studies are required to elucidate the process.

Other loci on chromosome 6

Tumor necrosis factor (TNF) is a pro-inflammatory cytokine whose serum levels are raised in patients with BD and inhibition of TNF by various biologic drugs is effective in controlling the disease. The gene for TNF is on human chromosome 6. Many SNPs in TNF have been analyzed since the −1031 promoter allele was associated with BD in UK patients.^[26] Subsequently, −1031C has been linked with the disease in Turkish and Tunisian with BD.^{[27, 28]} Similarly, −308A/G in TNF has been associated with disease in Iranian Azeri Turks but not in Palestinian patients.^{[4, 29]} In Korean patients, none of the tested TNF SNP were associated with the disease.^[30] Therefore, SNPs in this important proinflammatory cytokine are associated with BD in many populations, but which polymorphism would appear to differ in different geographical cohorts.

Non-chromosome 6 associations

(a) Lymphocyte Signaling

Genome-wide analysis studies (GWAS) in a cohort of 379 Korean patients with Behçet’s and 800 controls have found a significant association with the GTPase of immunity associated protein family (GIMAP) gene cluster on chromosome 7q36.1. Genome-wide significance was found using a dominant model of inheritance for 5 SNPs within the GIMAP cluster (rs1608157, rs1522596, rs10266069, rs10256482, and rs2286900) with subsequent fine-mapping efforts finding an association with 10 additional markers although these were distributed across 3 blocks of linkage disequilibrium (LD) within this ~200 kbp region. Replication of the signal was achieved for SNPs mapping to the GIMAP2 and GIMAP4 regions, but not GIMAP1, in a Japanese cohort. Additionally, functional data of 31 BD patients showed significantly reduced expression of GIMAP1, GIMAP4 (both in CD4 T cells), and GIMAP2 (in CD8 T cells) when compared to controls (P < 0.05). Further to this, using a luciferase reporter assay they found that a construct having the risk allele (C) within rs1608157 locus (the most significant GWAS SNP) had reduced luciferase activity compared to the G allele implying that the polymorphism negatively affected GIMAP4 promoter activity. However, GIMAP SNPs were not validated in the UK or European patients with BD.^{[31, 32]}

Several other genes linked to the severity of particular manifestations of BD, associated with T cell signaling, have been reported and many of these have been validated in multiple studies. Certain gene polymorphisms associated with other autoimmune diseases, including protein tyrosine phosphatase non-receptor type 22 (PTPN22) and cytotoxic T-lymphocyte antigen 4 (CTLA-4), were shown to have a negative association—results that support the concept that BD is an autoinflammatory rather than an autoimmune disease.^[33] PTPN22 R620W polymorphism is most prevalent in northern Europe and protective against Mycobacteria tuberculosis but induces susceptibility to Mycobacteria leprae, which was more common in the Middle East, where PTPN22 R620W is very rare. PTPN22 R620W was inversely associated with BD in patients particularly in a UK cohort.^[34] This association was seen neither in patients from the Middle East and Southern Europe, where the prevalence of the polymorphism is much lower, nor in Chinese patients with BD, where the polymorphism does not exist.^{[35, 36, 37]} PTPN22 encodes a protein Lyp that binds to Csk, an inhibitory kinase, and together disrupt T cell signaling through dephosphorylation of Lck. The R620W SNP disrupts this initial interaction.^[38] In myeloid cells, PTPN22 positively increases type 1 interferon production, and dependent viral infection and gut homeostasis.^[39] In mice, carriage of 620W led to increased autoreactive T cell clones in new bone-marrow emigrants.^[40] Inhibition of PTPN22 activity restored central B cell tolerance.^[41] Superresolution microscopy showed that while wild type PTPN22 forms large clusters in unstimulated T cells and dissociates on stimulation by LFA-1, on ICAM-1 monolayers, enabling association with Lck and ZAP70 at the leading edge. By comparison, PTPN22 620W was not retained at the leading edge resulting in increased LFA-1 clustering and integrin-mediated adhesion.^[42] PTPN22 is also highly expressed in neutrophils. R620W enhanced neutrophil migration and increased Ca+ release, in humans, and is a critical regulator of FcR activation and interaction with immune complexes by neutrophils in mice.^{[43, 44, 45]}

For CTLA4 polymorphisms, the situation is even more complex. Early studies on Turkish and Tunisian patients with BD found associations with certain manifestations of the disease.^{[46, 47]} However, studies on the UK, Middle Eastern, and Chinese patients showed no association.^{[48, 49, 50]} Meta-analysis of eight studies including those above showed no association with BD.^[51] However, whether certain features of BD are linked to CTLA4 polymorphisms or they are associated with other combinations of genes will need larger studies.

A signal transducer and activator of transcription 4 (STAT4) association are supported by a study on Han Chinese patients with BD, with SNP rs7574070 and additional SNPs rs7572482 and rs897200, all associated although in strong LD. In addition, Hou et al. showed that BD patients’ homozygous for the rs897200 risk allele had increased expression of the STAT4 mRNA and protein, as well as IL-17, which is regulated by the gene.^[52] Furthermore, ENCODE data reveal that the rs897200 polymorphism determines the formation of a transcription factor binding site, and the functional studies by the authors indicate that this SNP is likely to be an expression of quantitative trait locus (eQTL) of STAT4. Association of STAT4 was further reported in patients with BD from Turkey and Iran.^{[24, 53]} Although a STAT4 gene association (rs7574865) has been reported in a number of autoimmune diseases including RA and systemic lupus erythematosus, variants associated with BD are in low LD, with this particular marker implying that such associations are independent of each other and, therefore, may confer different mechanisms of pathogenesis.^[24]

Therefore, several polymorphisms in genes involved with T cell signaling are reported in BD but not in all populations. However, that many genes signaling through different, though interconnected, pathways in different geographical populations influence the same outcome that supports a Darwinian concept of BD [Figure 2].

Figure 2

GIMAP, and PTPN22, identified with BD in different geographic populations are linked to many similar molecules that are involved in lymphocyte activation such as STAT3 (signal transducer and activator of transcription 3). Image in Biorender.

Conclusion

The strongest association with BD is HLA-B*51 on chromo-some 6, a link that has been confirmed in several studies in different geographical groups. HLA-B*51 is derived from Neanderthal admixture and is likely to have been maintained in modern humans as part of a protective mechanism against a disease that affected both archaic and modern humans. Beyond chromosome 6, SNPs in genes involved in lymphocyte signaling and function, along with ubiquitination and pathogen recognition, which may be protective in evolutionary time but when brought together represent further evidence supporting a Darwinian view of BD. The identification of different polymorphisms in different geographical cohorts of patients with BD supports the accumulation of small mutations that, when combined with HLA-B*51, lead to the disease. However, it is now becoming clear that these different, positively selected SNPs are, in many cases, affecting the same biological pathways. An excellent, comprehensive genetic and epigenetic analysis identified a heritability of 60% in BD incorporating mutations in genes associated with the immune response, including interferon-gamma and IL-12 production, regulation of lymphocyte-mediated immunity cytokine–cytokine receptor interaction, and regulation of innate immune response. Genetic cumulative risk scores supported the known prevalence of BD globally.^[74] These new insights cast light on the biochemical pathway involved in BD, beyond links with individual genes. This emerging understanding of the evolutionary, genetic, and biochemical basis of BD may explain much of the variation in BD prevalence with geography and may even suggest new directions for therapeutic research. In conclusion, protective genes inherited from archaic humans in the distant past, interacting with other polymorphisms in modern humans, may explain the strong link that we observe between HLA-B and BD.