Advertisement

Genotype-Phenotype Study of the Middle Gangetic Plain in India Shows Association of rs2470102 with Skin Pigmentation

Open ArchivePublished:November 17, 2016DOI:https://doi.org/10.1016/j.jid.2016.10.043
      Our understanding of the genetics of skin pigmentation has been largely skewed towards populations of European ancestry, imparting less attention to South Asian populations, who behold huge pigmentation diversity. Here, we investigate skin pigmentation variation in a cohort of 1,167 individuals in the Middle Gangetic Plain of the Indian subcontinent. Our data confirm the association of rs1426654 with skin pigmentation among South Asians, consistent with previous studies, and also show association for rs2470102 single nucleotide polymorphism. Our haplotype analyses further help us delineate the haplotype distribution across social categories and skin color. Taken together, our findings suggest that the social structure defined by the caste system in India has a profound influence on the skin pigmentation patterns of the subcontinent. In particular, social category and associated single nucleotide polymorphisms explain about 32% and 6.4%, respectively, of the total phenotypic variance. Phylogeography of the associated single nucleotide polymorphisms studied across 52 diverse populations of the Indian subcontinent shows wide presence of the derived alleles, although their frequencies vary across populations. Our results show that both polymorphisms (rs1426654 and rs2470102) play an important role in the skin pigmentation diversity of South Asians.

      Abbreviations:

      LD (linkage disequilibrium), MGP (Middle Gangetic Plain), MI (melanin index), SNP (single nucleotide polymorphism)

      Introduction

      Human skin color varies remarkably within and among populations. This variation has been mainly attributed to the differences in the amount, type, and distribution of melanin in melanocytes and the ratio of eumelanin to phaeomelanin (
      • Rees J.L.
      Genetics of hair and skin color.
      ,
      • Thody A.J.
      • Higgins E.M.
      • Wakamatsu K.
      • Ito S.
      • Burchill S.A.
      • Marks J.M.
      Pheomelanin as well as eumelanin is present in human epidermis.
      ). Human skin color variation is a polygenic trait, and it has been characterized by a number of major genes, modifier genes, and environmental influences (
      • Miller C.T.
      • Beleza S.
      • Pollen A.A.
      • Schluter D.
      • Kittles R.A.
      • Shriver M.D.
      • et al.
      cis-Regulatory changes in Kit ligand expression and parallel evolution of pigmentation in sticklebacks and humans.
      ,
      • Sturm R.A.
      • Duffy D.L.
      Human pigmentation: painting by numbers or ancestry?.
      ). It has been suggested that skin pigmentation is largely shaped by different levels of UVR via natural selection (
      • Chaplin G.
      Geographic distribution of environmental factors influencing human skin coloration.
      ,
      • Jablonski N.G.
      • Chaplin G.
      The evolution of human skin coloration.
      ,
      • Jablonski N.G.
      • Chaplin G.
      Human skin pigmentation as an adaptation to UV radiation.
      ). Furthermore, there has been evidence of population-specific signatures of positive selection for different pigmentation genes (
      • Izagirre N.
      • García I.
      • Junquera C.
      • De La Rúa C.
      • Alonso S.
      A scan for signatures of positive selection in candidate loci for skin pigmentation in humans.
      ,
      • Lao O.
      • De Gruijter J.
      • Van Duijn K.
      • Navarro A.
      • Kayser M.
      Signatures of positive selection in genes associated with human skin pigmentation as revealed from analyses of single nucleotide polymorphisms.
      ,
      • Voight B.F.
      • Kudaravalli S.
      • Wen X.
      • Pritchard J.K.
      A map of recent positive selection in the human genome.
      ), reviewed in
      • Anno S.
      • Ohshima K.
      • Abe T.
      Approaches to understanding adaptations of skin color variation by detecting gene environment interactions.
      .
      SLC24A5 stands out in the list of 170 pigmentation genes, showing a major influence on pigmentation variation. The key single nucleotide polymorphism (SNP), rs1426654, in the third exon explains 25–38% of the pigmentation differences between Europeans and West Africans (
      • Lamason R.L.
      • Mohideen M.A.P.K.
      • Mest J.R.
      • Wong A.C.
      • Norton H.L.
      • Aros M.C.
      • et al.
      SLC24A5, a putative cation exchanger, affects pigmentation in zebrafish and humans.
      ). Functional assessment of the SNP suggests that it is involved in reduced melanogenesis, thereby making the skin lighter (
      • Cook A.L.
      • Chen W.
      • Thurber A.E.
      • Smit D.J.
      • Smith A.G.
      • Bladen T.G.
      • et al.
      Analysis of cultured human melanocytes based on polymorphisms within the SLC45A2/MATP, SLC24A5/NCKX5, and OCA2/P loci.
      ,
      • Ginger R.S.
      • Askew S.E.
      • Ogborne R.M.
      • Wilson S.
      • Ferdinando D.
      • Dadd T.
      • et al.
      SLC24A5 encodes a trans-Golgi network protein with potassium-dependent sodium-calcium exchange activity that regulates human epidermal melanogenesis.
      ,
      • Sturm R.A.
      Molecular genetics of human pigmentation diversity.
      ,
      • Tsetskhladze Z.R.
      • Canfield V.A.
      • Ang K.C.
      • Wentzel S.M.
      • Reid K.P.
      • Berg A.S.
      • et al.
      Functional assessment of human coding mutations affecting skin pigmentation using zebrafish.
      ).
      India has a wide variation in skin color ranging from fair to wheatish, brown, and dark brown skin tones (
      • Jaswal I.
      Pigmentary variation in Indian populations.
      ,
      • Hourblin V.
      • Nouveau S.
      • Roy N.
      • de Lacharrière O.
      Skin complexion and pigmentary disorders in facial skin of 1204 women in 4 Indian cities.
      ). Nevertheless, only a handful of studies featuring the skin pigmentation variation (
      • Jonnalagadda M.
      • Norton H.
      • Ozarkar S.
      • Kulkarni S.
      • Ashma R.
      Association of genetic variants with skin pigmentation phenotype among populations of west Maharashtra, India.
      ,
      • Mallick C.B.
      • Iliescu F.M.
      • Möls M.
      • Hill S.
      • Tamang R.
      • Chaubey G.
      • et al.
      The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
      ,
      • Stokowski R.P.
      • Pant P.
      • Dadd T.
      • Fereday A.
      • Hinds D.A.
      • Jarman C.
      • et al.
      A genomewide association study of skin pigmentation in a South Asian population.
      ) have been reported to date. The Indo-Gangetic Plain is one of the most densely populated regions in the world, accounting for 40% of the Indian population. It has been further divided into the Upper, Middle, and Lower Gangetic Plains. The Middle Gangetic Plain (MGP) extends over a stretch of approximately 600 km from east to west and 330 km from north to south, covering the states of Uttar Pradesh and Bihar and comprising the fertile banks of the Ganges and its tributaries (Ghaghra, Gandak, and Kosi) (Figure 1). MGP is a home to more than 100 million people following the complex, yet stratified and well-defined, caste system (
      • Chaubey G.
      • Metspalu M.
      • Kivisild T.
      • Villems R.
      Peopling of South Asia: investigating the caste–tribe continuum in India.
      ,
      • Schwartzberg J.
      Caste regions of the north Indian plains in structure and change in Indian society.
      ). The caste system, unique to South Asia, is characterized by multilayered endogamy (see Supplementary Materials online). Hence, this socio-culturally diverse region that includes various ethnic groups and castes provides us an excellent model to study the local skin pigmentation variation of the region and, in particular, investigate how the social hierarchy influences skin pigmentation variation.
      Figure 1
      Figure 1Map of India showing the region under study (Middle Gangetic Plain). Dots represent the sampling locations, which include seven districts of Uttar Pradesh and eight districts of Bihar.
      Our main objectives of the study were first, to characterize the skin pigmentation diversity of the region (cohort 1); second, to evaluate the association of the SLC24A5 variants with skin pigmentation (cohort 2); and third, to assess the frequency of distribution of SLC24A5 SNPs across the Indian subcontinent (cohort 3). The details of the three cohorts included in the study are specified in the Methods section (see Supplementary Tables S1–S3 online).

      Results

      Pigmentation diversity

      Analysis of cohort 1 (n = 1,167) shows wide variation in skin color, with melanin indices (MIs) ranging from 23.0−90.4 (see Supplementary Figure S1 online) and an average MI of 59.38 ± 10.78. Among the 27 ethnic groups assessed, the highest average MI was observed among Manjhis (70.83 ± 8.13), and the lowest average MI was among Brahmins of Uttar Pradesh (45.13 ± 5.91) (see Supplementary Table S4 online). The maximum variation in MI was seen among Bhagats (66.75 ± 10.26). However, there was no significant difference in average MI for males and females (males = 59.0, females = 60.1; P = 0.1036) included in cohort 1. One-way analysis of variance suggests that skin color varies significantly among ethnic groups (P < 2.2 × 10–16) and social categories (P < 2.2 × 10–16). This effect of ethnic group on skin color is in accordance with previous studies including six endogamous groups of Eastern Nepal (
      • Williams-Blangero S.
      • Blangero J.
      Skin color variation in eastern Nepal.
      ) and West Maharashtra (
      • Jonnalagadda M.
      • Ozarkar S.
      • Ashma R.
      • Kulkarni S.
      Skin pigmentation variation among populations of West Maharashtra, India.
      ). The comparison of skin color measurements among the four social categories (general, scheduled caste, other backward classes, and religious group) assessed in this study indicates that the general category (traditionally comprising of the upper and middle castes in the caste system) shows the lowest average MI (see Supplementary Figure S2 and Supplementary Materials online).

      Resequencing SLC24A5

      The human SLC24A5 gene includes nine exons spanning over 21.7 kilo base pairs. We sequenced all the exons (4.4 kilo base pairs) and the adjoining flanking regions. Of the 3,525 base pairs resequenced among 374 individuals (cohort 2), 10 variants were identified (see Supplementary Table S5 online). The variants included two common polymorphisms (Minor Allele Frequency > 5%), which were subsequently tested for association with skin pigmentation. One of them, rs1426654, is a nonsynonymous SNP in the third exon, whereas rs2470102 is located in intron 8 of SLC24A5 and 3′ untranslated region of MYEF2. Of these, the evidence for association of the rs1426654 SNP with skin color has been well reported (
      • Ang K.C.
      • Ngu M.S.
      • Reid K.P.
      • Teh M.S.
      • Aida Z.S.
      • Koh D.X.
      • et al.
      Skin color variation in Orang Asli tribes of peninsular Malaysia.
      ,
      • Beleza S.
      • Johnson N.A.
      • Candille S.I.
      • Absher D.M.
      • Coram M.A.
      • Lopes J.
      • et al.
      Genetic architecture of skin and eye color in an African-European admixed population.
      ,
      • Lamason R.L.
      • Mohideen M.A.P.K.
      • Mest J.R.
      • Wong A.C.
      • Norton H.L.
      • Aros M.C.
      • et al.
      SLC24A5, a putative cation exchanger, affects pigmentation in zebrafish and humans.
      ,
      • Mallick C.B.
      • Iliescu F.M.
      • Möls M.
      • Hill S.
      • Tamang R.
      • Chaubey G.
      • et al.
      The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
      ,
      • Stokowski R.P.
      • Pant P.
      • Dadd T.
      • Fereday A.
      • Hinds D.A.
      • Jarman C.
      • et al.
      A genomewide association study of skin pigmentation in a South Asian population.
      ), with the derived A allele making the skin lighter. However, to our knowledge, the role of rs2470102 in skin pigmentation has not been elucidated to date. The remaining eight variants were observed only in the heterozygous state. Graphical representation of the SNPs in relation to the exon-intron structure of the gene and the region resequenced are shown in Figure 2. The electropherograms for rs1426654 and rs2470102 are provided in Supplementary Figure S3 online.
      Figure 2
      Figure 2The structure of the human SLC24A5 gene (chromosome 15, 48120372–48143272: ENST00000341459). Exons of the gene are shown in pink and introns in purple. The SNPs rs1426654 (A>G) in the third exon and rs2470102 (A>G) in the eight intron and 3′ untranslated region variant of MYEF2 (chromosome 15: ENST00000324324) has been highlighted. bp, base pairs; UTR, untranslated region.

      Genotype-phenotype association

      The beeswarm plot illustrates the distribution of MI for the rs1426654 and rs2470102 genotypes (Figure 3a and b). We found that the mean MI for individuals with AA genotypes is lower than for those with AG and GG genotypes (Figure 3a and b) for both SNPs. However, to further our understanding, we used a linear model (general genotype model) to assess the effect of genotype on skin pigmentation tested individually for each of the two SNPs (Table 1, Models 1 and 2). After controlling for sex and population, we found both SNPs, rs1426654 and rs2470102, to be significantly associated with skin color (P = 4.493 × 10–7 and P = 5.79 × 10–7, respectively). These P-values remained significant even after using the conservative Bonferroni corrections for multiple testing (the corrected P-values being 8.9 × 10–7 and 1.158 × 10–6, respectively). rs1426654 and rs2470102 SNPs lie 7,010 base pairs apart and were found to be in high linkage disequilibrium (LD) (D′ = 1, r2 = 0.87) when examined using Gujarati Indians included in HapMap (GIH; Gujarati Indians in Houston, HapMap release 27, Feb 2009, available at http://www.hapmap.org). A relevant question to ask, then, was whether the observed association signal for rs2470102 is entirely due to the known causative SNP rs1426654. For this, we used another linear model (model 3), which was adjusted for rs1426654. We found that the effect of rs2470102 was still statistically significant (P = 0.03), thereby suggesting that rs2470102 has its independent effect on skin pigmentation variation among Indian populations. The details of the results for the full model and their summary are given in Supplementary Table S6 online.
      Figure 3
      Figure 3Association of rs1426654 and rs2470102 genotypes with melanin index. A beeswarm dotplot showing (a) association of the rs1426654 genotype with melanin index in the studied cohort and (b) association of the rs2470102 genotype with melanin index in the studied cohort. The black line indicates the average melanin index for each genotype: AA, AG, and GG.
      Table 1Linear models tested (Models 1–3) and their descriptions
      ModelModel Description
      Model 1MI ∼ rs1426654 (df = 2) + Ethnic.group (df = 22) + Sex (df = 1)
      Model 2MI ∼ rs2470102 (df = 2) + Ethnic.group (df = 22) + Sex (df = 1)
      Model 3MI ∼ rs2470102 (df = 2) + rs1426654 (df = 2) + Ethnic.group (df = 22) + Sex (df = 1) + rs1426654 × Ethnic.group (df = 36)
      In line with the evidence of the association of SNP rs2470102 with observed skin pigmentation, the effect of rs2470102 was further assessed by comparing the difference in the estimated means of skin measurements for individuals with the genotypes (AA, AG, and GG). It appeared that presence of two copies of the derived allele (AA) of rs2470102 made an individual, on average, 7.10 melanin units lighter compared with GG homozygotes (P = 0.0011) (see Supplementary Table S7 online). After adjusting for rs1426654 (model 3), individuals with AA genotypes were still, on average, 4.52 units lighter than those with GG genotypes (P = 0.03) (see Supplementary Table S7). Hence, taken together, we conclude that rs2470102 has an independent effect on skin pigmentation variation.
      rs2470102 has been earlier reported to create novel MRESS (microRNA recognition element seed site) (
      • Richardson K.
      • Lai C.-Q.
      • Parnell L.D.
      • Lee Y.-C.
      • Ordovas J.M.
      A genome-wide survey for SNPs altering microRNA seed sites identifies functional candidates in GWAS.
      ). To address this, we examined five web-based algorithms to identify the possible microRNAs targeting the SNP (see Supplementary Materials, Supplementary Table S8 and Supplementary Figure S4 online). Our in silico analyses are suggestive of mir-1180 interacting with rs2470102; however, it requires further confirmation using functional or gene reporter assays to have a complete understanding of the role of the SNP.

      Haplotype analysis

      We performed haplotype analyses to determine the relationship of haplotypes with skin color and social categories. We identified nine haplotypes (H1–H9) in the resequenced region. A median-joining network was constructed to study the relationship of the haplotypes (Figure 4). Most of the chromosomes in this study (97.6%) contribute to three major haplotypes (H1–H3) that differ in two polymorphisms (rs1426654 and rs2470102), and the others (2.4%) belong to the remaining haplotypes (H4–H9). Earlier studies have shown that rs2470102 is phylogenetically ancestral to rs1426654 (
      • Canfield V.A.
      • Berg A.
      • Peckins S.
      • Wentzel S.M.
      • Ang K.C.
      • Oppenheimer S.
      • et al.
      Molecular phylogeography of a human autosomal skin color locus under natural selection.
      ).
      Figure 4
      Figure 4Network analysis showing the relationship of nine haplotypes and distribution of the major haplotypes among four social categories. The upper panel shows the median-joining network of the nine haplotypes (H1–H9) defined by 10 variants (made using NETWORK software, available at www.fluxus-engineering.com). The colors represent the haplotype distribution among the four social categories included in the study. Each circle represents a haplotype. Circles are proportional to the number of chromosomes. The numbers in the figure correspond to the position of the variants as specified in . The lower panel shows the distribution of the three major haplotypes (H1–H3) among four social categories having zero, one, or two chromosomes in each haplotype. OBC, other backward classes; RG, religious group; SC, scheduled caste.

      Haplotype versus skin color

      When comparing the frequencies of haplotypes across the skin color distribution, we found that the light-skinned individuals had higher frequency of the H1 haplotype. Additionally, we observed a clear trend of gradual decrease in frequency of H1 haplotypes with increase in MIs (see Supplementary Figure S5 online). The frequency of H1 haplotypes for individuals with MI (30–40) was 96% and it declined to 37% for individuals with MI (70–80) (see Supplementary Figure S5).

      Haplotype versus social category

      When comparing the haplotypes across the four social categories studied (Figure 4), we noted the presence of all the four social categories in three major haplotypes (H1–H3); therefore, social category is not exclusive to any particular haplotype. A chi-square test further showed that the distribution of the haplotypes was significantly different across the social categories for H1 (P = 3.5 × 10–10) and H3 (P = 1.449 × 10–7) but not for H2 (P = 0.3475) (Figure 4). This further attests to the effect of social categories on skin pigmentation, because the predominance of individuals from general category (including traditionally mostly the upper and middle castes in the caste system) in haplotype H1 was observed, which is correlated with lighter skin color (Figure 4).

      Phylogeography of rs1426654 and rs2470102

      Because our association analyses ascertained that both rs1426654 and rs2470102 SNPs are important determinants of skin pigmentation variation among Indian populations, we were interested to study the phylogeography of the SNPs. For this, we genotyped 1,825 individuals belonging to 52 diverse populations, representing the Indian subcontinent (see Supplementary Table S3). The frequencies across the populations are shown in Figure 5.
      Figure 5
      Figure 5Geographical distribution of allele frequencies for rs1426654 and rs2470102 SNPs. (a) Pie charts showing the allelic frequency distribution of the rs1426654 polymorphism among caste populations. (b) Pie charts showing the allelic frequency distribution of the rs1426654 polymorphism among tribe populations. (c) Pie charts showing the allelic frequency distribution of the rs2470102 polymorphism among caste populations. (d) Pie charts showing the allelic frequency distribution of rs2470102 polymorphism among tribe populations. These caste and tribe populations are in total, represented by 1,825 individuals from 52 diverse populations of India (see ). SNP, single nucleotide polymorphism.
      The allele frequencies of the SNPs suggest that both are widely spread and highly polymorphic across the subcontinent (Figure 5, and see Supplementary Table S3). Both the SNPs exhibit similar allele frequencies among Indian populations (when grouped by language or social status) and among other world populations (see Supplementary Table S9 online). This can be further explained by the fact that these SNPs are in high LD in most of the populations studied (D′ = 0.64–1), except for Warli (D′ = 0.53) and Adi-Dravida (D′ = 0.23). The derived allele (A) frequencies ranged from 0.04–1.00 and 0.13–0.98 for rs1426654 and rs2470102, respectively (see Supplementary Table S3). The range of allele frequencies observed in this study (0.0–1.0) for rs1426654 is concordant with those observed in our previous study (0.03–1.00) (
      • Mallick C.B.
      • Iliescu F.M.
      • Möls M.
      • Hill S.
      • Tamang R.
      • Chaubey G.
      • et al.
      The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
      ). When populations were grouped by social status, caste populations showed higher A allele frequencies compared with tribes (0.72 vs. 0.39 for rs1426654 and 0.76 vs. 0.45 for rs2470102). This finding is consistent with another study, where researchers found that the caste populations were significantly lighter than tribal populations at a localized level, among the populations of West Maharashtra (
      • Jonnalagadda M.
      • Ozarkar S.
      • Ashma R.
      • Kulkarni S.
      Skin pigmentation variation among populations of West Maharashtra, India.
      ), on the basis of their skin color measurements. A study focusing on 11 endogamous Indian populations also found higher rs1426654-A allele frequencies among castes than tribes (
      • Mukherjee M.
      • Mukerjee S.
      • Sarkar-Roy N.
      • Ghosh T.
      • Kalpana D.
      • Sharma A.K.
      Polymorphisms of four pigmentation genes (SLC45A2, SLC24A5, MC1R and TYRP1) among eleven endogamous populations of India.
      ). However, some tribes show an exceptionally high frequency of the rs1426654-A allele, for example, Gujjar (Jammu and Kashmir) having 1.00 and Meena (Rajasthan) as 0.91 (see Supplementary Table S3). These tribes have been also known to be fair skinned (Joshua

      Joshua Project, https://joshuaproject.net/; 2016 (accessed 14 April 2016).

      ). Brahmins traditionally belonging to upper castes in the social hierarchy of the caste system (see Supplementary Materials), irrespective of their geographical locations (north: Kashmiri Pandits of Jammu and Kashmir, Pandits of Haryana, Brahmins of Uttar Pradesh; south: Havik of Karnataka), show similar frequencies of the rs1426654-A variant (range = 0.83–1.00) (see Supplementary Table S3). When populations were grouped by their linguistic affiliations, the frequency of the rs1426654- and rs2470102-derived alleles (A) were found to be very low among Austroasiatic and Tibeto-Burman speakers (see Supplementary Table S9).

      Discussion

      This study, undertaken in the MGP, helps us understand microdifferentiation patterns in skin pigmentation measures of the region. The apparent phenotypic variation (MI = 23.0–90.4) observed in this small geographical region covering an area of 144,409 km2 can be mainly attributed to the socio-cultural boundaries, superimposed with a high level of endogamy (
      • Bhasin M.
      • Walter H.
      Genetics of castes and tribes of India (Peoples of Indian region).
      ,
      • Chaubey G.
      • Metspalu M.
      • Kivisild T.
      • Villems R.
      Peopling of South Asia: investigating the caste–tribe continuum in India.
      ,
      • Karve I.
      Kinship organization in India.
      ). Our data show that skin color varies significantly across different social categories (see Supplementary Figure S2) and ethnic groups (see Supplementary Table S4). On the basis of skin color measurements, we found that the mean MIs of individuals of the general category were significantly different compared with other castes (P = 0.00019). On the other hand, the difference in skin measurements of individuals under other backward classes and scheduled caste was not significant (P = 0.43). A similar picture was also evident in the haplotype analyses. The number of individuals having both the chromosomes in the H1 haplotype was significantly higher in the general category (64%) than other social groups (27%) (P = 1 × 10–10) (Figure 4). Taken together, our findings provide evidence for the primacy of socio-cultural factors in skin pigmentation patterns observed in the subcontinent. This is also concordant with previous studies, which suggested that a UVR-based selection model alone cannot explain the entire skin pigmentation variation in Indian populations but is rather an interplay between selection and the demographic history of the populations (
      • Mallick C.B.
      • Iliescu F.M.
      • Möls M.
      • Hill S.
      • Tamang R.
      • Chaubey G.
      • et al.
      The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
      ,
      • Mukherjee M.
      • Mukerjee S.
      • Sarkar-Roy N.
      • Ghosh T.
      • Kalpana D.
      • Sharma A.K.
      Polymorphisms of four pigmentation genes (SLC45A2, SLC24A5, MC1R and TYRP1) among eleven endogamous populations of India.
      ).
      Our association tests suggest that both rs1426654 and rs2470102 have independent effects on skin color; that of rs1426654 was shown to be convincingly replicated earlier in other Indian populations (
      • Jonnalagadda M.
      • Norton H.
      • Ozarkar S.
      • Kulkarni S.
      • Ashma R.
      Association of genetic variants with skin pigmentation phenotype among populations of west Maharashtra, India.
      ,
      • Mallick C.B.
      • Iliescu F.M.
      • Möls M.
      • Hill S.
      • Tamang R.
      • Chaubey G.
      • et al.
      The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
      ,
      • Stokowski R.P.
      • Pant P.
      • Dadd T.
      • Fereday A.
      • Hinds D.A.
      • Jarman C.
      • et al.
      A genomewide association study of skin pigmentation in a South Asian population.
      ) and world populations (
      • Adhikari K.
      • Fontanil T.
      • Cal S.
      • Mendoza-Revilla J.
      • Fuentes-Guajardo M.
      • Chacón-Duque J.-C.
      • et al.
      A genome-wide association scan in admixed Latin Americans identifies loci influencing facial and scalp hair features.
      ,
      • Ang K.C.
      • Ngu M.S.
      • Reid K.P.
      • Teh M.S.
      • Aida Z.S.
      • Koh D.X.
      • et al.
      Skin color variation in Orang Asli tribes of peninsular Malaysia.
      ,
      • Beleza S.
      • Johnson N.A.
      • Candille S.I.
      • Absher D.M.
      • Coram M.A.
      • Lopes J.
      • et al.
      Genetic architecture of skin and eye color in an African-European admixed population.
      ,
      • Lamason R.L.
      • Mohideen M.A.P.K.
      • Mest J.R.
      • Wong A.C.
      • Norton H.L.
      • Aros M.C.
      • et al.
      SLC24A5, a putative cation exchanger, affects pigmentation in zebrafish and humans.
      ). rs2470102 was found to be associated with brown eye color in the African-European admixed population of Cape Verde (
      • Beleza S.
      • Johnson N.A.
      • Candille S.I.
      • Absher D.M.
      • Coram M.A.
      • Lopes J.
      • et al.
      Genetic architecture of skin and eye color in an African-European admixed population.
      ) and to be associated with melanoma, in a versatile gene-based test called VEGAS (
      • Liu J.Z.
      • Mcrae A.F.
      • Nyholt D.R.
      • Medland S.E.
      • Wray N.R.
      • Brown K.M.
      • et al.
      A versatile gene-based test for genome-wide association studies.
      ) but not with skin pigmentation. Therefore, our association and haplotype analyses taken together allow us to infer that rs1426654 and rs2470102 as a “two SNP model” can better explain the variation in skin color among Indian populations than each SNP individually. We found that the social category and associated SNPs explain 32% and 6.4%, respectively, accounting for a total of 38.4% of the variability in skin pigmentation. Of the 32% of phenotypic variance explained by the social category, 37.4% is due to variation in pigmentation among the social categories (akin to 11.97% of the total variability) (see Supplementary Figure S6 online). We found that rs1426654 (5.37%) had a slightly larger effect than rs2470102 (4.99%), and both the SNPs taken together (6.41%) can better explain the variation in skin pigmentation than each of the SNPs individually.
      In summary, this study involving the populations of the MGP refines our existing knowledge of skin pigmentation genetics in South Asia. We report the association of the rs2470102 SNP with skin pigmentation, and spatial distribution of the SNP marks its ubiquitous presence over the Indian subcontinent. Our haplotype analyses show the presence of three major haplotypes (H1–H3), among which the frequency of haplotype H1 is higher in individuals with lower MIs. Our data reflects effect of socio-cultural factors in shaping the skin pigmentation variation of the subcontinent. This aspect observed in our study among Indian populations helps us understand the possible paradigms contributing to the global spectrum of skin color.

      Materials and Methods

      Study area

      This study was conducted in the MGP across North Bihar and the eastern part of Uttar Pradesh (Figure 1). These samples were collected through multiple field visits at various districts of the MGP during October 2010 through October 2011.

      Study subjects and cohorts

      This study comprised three cohorts. Cohort 1 consisted of 1,167 unrelated subjects (27 ethnic groups, 759 men and 408 women), who were used for studying the phenotypic variation of the MGP (see Supplementary Table S1). Out of these, a subset of individuals was included in cohort 2 (n = 448, 25 ethnic groups), which was recruited for phenotype-genotype study (see Supplementary Table S2). About a 2- to 4-ml blood sample was taken from each individual. Cohort 3 comprised 1,825 individuals, including 52 populations (24 caste and 27 tribe populations) across India (see Supplementary Table S3). This cohort was used for studying the phylogeography of the associated SNPs observed in the study. These samples are part of the DNA bank of CSIR-Centre for Cellular and Molecular Biology (Hyderabad). The Institutional Ethical Committee of Centre for Cellular and Molecular Biology (Hyderabad) approved the study, and prior permission was also obtained from local authorities. Informed written consent was obtained from each subject.

      Phenotype data collection

      The skin reflectance readings were measured using Dermaspectrometer (Cortex Technology, Hadsund, Denmark). Measurements were performed in duplicate at the inner upper arm and forehead. However, only readings from the former were used for genetic analyses. The quantitative assessment of melanin content, as obtained with the Dermaspectrometer, is expressed as the MI. Subjects were excluded if they reported the use of skin ointments on the measurement areas. Also, individuals reported to have any pigmentation disorders or skin diseases were excluded. Skin measurements were taken for cohorts 1 and 2.

      Sequencing

      Genomic DNA was extracted from whole blood using standard methods (
      • Thangaraj K.
      • Joshi M.B.
      • Reddy A.G.
      • Gupta N.J.
      • Chakravarty B.
      • Singh L.
      CAG repeat expansion in the androgen receptor gene is not associated with male infertility in Indian populations.
      ). Primers for all the nine exons of SLC24A5 were designed using the National Center for Biotechnology Information’s Primer-BLAST (available at http://www.ncbi.nlm.nih.gov/tools/primer-blast), MacVector (MacVector, Apex, NC), and the Amplify 3X Software (available at http://engels.genetics.wisc.edu/amplify) (see Supplementary Table S10 online). All nine exons and its flanking regions were amplified through PCR using these primer pairs and Emerald AMP GT Master Mix (DSS Takara Bio India (Pvt) Ltd, Delhi, India) among individuals of cohort 2. PCR amplification was performed using the following conditions: initial denaturation at 95 °C for 5 minutes, followed by 35 cycles of 95 °C for 30 seconds, 62 °C for 25 seconds, and 72 °C for 3 minutes, and final extension at 72 °C for 3 minutes. PCR product was cleaned using Exo-SAP (Exonuclease-Shrimp Alkaline Phosphatase, USB Corporation, Cleveland, OH), according to the manufacturer’s protocol. Exo-SAP–treated amplicons were sequenced directly using a BigDye terminator cycle sequencing kit, version 3.1 (Applied Biosystems, Carlsbad, CA) on an ABI 3730xl DNA analyzer (Applied Biosystems, Foster city, CA). Sequences were assembled with the reference sequence and analyzed using AutoAssembler software (Applied Biosystems). Observed variants were validated and confirmed by visual confirmation of electropherograms. The genotype frequencies were in Hardy-Weinberg equilibrium (P > 0.05) among the variants studied.

      Statistical analysis

      The difference in average MI of males and females was calculated by Welch two-sample t test. One-way analysis of variance was used to evaluate whether skin pigmentation varied across ethnic groups and social categories. The associations between polymorphisms of SLC24A5 and skin pigmentation measures were analyzed using three linear models (models 1–3) in R version 3.2.4 (available at https://www.r-project.org/) using F test (Table 1). The summary and the full description of the models is given in Supplementary Table S6. For multiple comparisons, we used post hoc Tukey like-comparison tests, available in the multicomp package in R. All statistical analyses were done using R version 3.1.2 (R Development Core Team). For all statistical tests, significance was defined as P < 0.05 unless specified. Sequence data (3,525 base pairs) were phased using PHASE v 2.1.1 (
      • Stephens M.
      • Donnelly P.
      A comparison of Bayesian methods for haplotype reconstruction from population genotype data.
      ) implemented in DnaSP v 5.1 (
      • Librado P.
      • Rozas J.
      DnaSP v5: a software for comprehensive analysis of DNA polymorphism data.
      ) and haplotypes were inferred. A haplotype network was drawn using the median-joining algorithm by NETWORK software (Fluxus Engineering, Suffolk, UK; available at fluxus-engineering.com) to study the relationship between haplotypes. LD block structure was examined using Haploview software, version 4.2 (
      • Barrett J.C.
      • Fry B.
      • Maller J.
      • Daly M.J.
      Haploview: analysis and visualization of LD and haplotype maps.
      ) using default parameters. For calculating LD estimates for populations in cohort 3, only those populations that had at least 10 individuals were included. Pairwise LD between polymorphisms is expressed as D′ (a normalized measure for assessing LD).

      Conflict of Interest

      The authors state no conflict of interest.

      Acknowledgments

      We thank B. B. Thakur, K. C. Sinha, Gyan Bhusan, Ajit, Jitendra Kumar, Ram Lakhan Jha, Sushil Upadhyay, Arun Pandey, Anil Upadhyay, Anand Rai, Shakti P. Singh, and Manoj Rai who assisted in the fieldwork. We kindly acknowledge Priya Moorjani and Urmos Vosa for their suggestions. We are also grateful to Richard Villems and Toomas Kivisild for their helpful comments on the manuscript. This work was supported by the CSIR Network Project (EpiHeD-BSC0118 to KT), DBT-RA program (to AM), 2016 Post-Doc Development Program-Pusan National University, S.Korea (fellowship to AM), P. R. Foundation (Uttar Pradesh), ICMR-SRF program (to SN), and Estonian Personal Grant (PUT-766 to GC). CBM and GC acknowledge financial support from European Union Regional Development fund through the Centre of Excellence in Genomics to Estonian Biocentre and University of Tartu and Estonian Institutional Research grant (IUT24-1).

      Supplementary Material

      References

        • Adhikari K.
        • Fontanil T.
        • Cal S.
        • Mendoza-Revilla J.
        • Fuentes-Guajardo M.
        • Chacón-Duque J.-C.
        • et al.
        A genome-wide association scan in admixed Latin Americans identifies loci influencing facial and scalp hair features.
        Nat Commun. 2016; 7: 10815
        • Ang K.C.
        • Ngu M.S.
        • Reid K.P.
        • Teh M.S.
        • Aida Z.S.
        • Koh D.X.
        • et al.
        Skin color variation in Orang Asli tribes of peninsular Malaysia.
        PloS One. 2012; 7: e42752
        • Anno S.
        • Ohshima K.
        • Abe T.
        Approaches to understanding adaptations of skin color variation by detecting gene environment interactions.
        Expert Rev Mol Diagn. 2010; 10: 987-991
        • Barrett J.C.
        • Fry B.
        • Maller J.
        • Daly M.J.
        Haploview: analysis and visualization of LD and haplotype maps.
        Bioinformatics. 2005; 21: 263-265
        • Beleza S.
        • Johnson N.A.
        • Candille S.I.
        • Absher D.M.
        • Coram M.A.
        • Lopes J.
        • et al.
        Genetic architecture of skin and eye color in an African-European admixed population.
        PLoS Genet. 2013; 9: e1003372
        • Bhasin M.
        • Walter H.
        Genetics of castes and tribes of India (Peoples of Indian region).
        India:Kamla-Raj Enterprises, New Delhi2001
        • Canfield V.A.
        • Berg A.
        • Peckins S.
        • Wentzel S.M.
        • Ang K.C.
        • Oppenheimer S.
        • et al.
        Molecular phylogeography of a human autosomal skin color locus under natural selection.
        G3 Genes Genomes Genet. 2013; 3: 2059-2067
        • Chaplin G.
        Geographic distribution of environmental factors influencing human skin coloration.
        Am J Phys Anthropol. 2004; 125: 292-302
        • Chaubey G.
        • Metspalu M.
        • Kivisild T.
        • Villems R.
        Peopling of South Asia: investigating the caste–tribe continuum in India.
        Bioessays. 2007; 29: 91-100
        • Cook A.L.
        • Chen W.
        • Thurber A.E.
        • Smit D.J.
        • Smith A.G.
        • Bladen T.G.
        • et al.
        Analysis of cultured human melanocytes based on polymorphisms within the SLC45A2/MATP, SLC24A5/NCKX5, and OCA2/P loci.
        J Invest Dermatol. 2008; 129: 392-405
        • Ginger R.S.
        • Askew S.E.
        • Ogborne R.M.
        • Wilson S.
        • Ferdinando D.
        • Dadd T.
        • et al.
        SLC24A5 encodes a trans-Golgi network protein with potassium-dependent sodium-calcium exchange activity that regulates human epidermal melanogenesis.
        J Biol Chem. 2008; 283: 5486-5495
        • Hourblin V.
        • Nouveau S.
        • Roy N.
        • de Lacharrière O.
        Skin complexion and pigmentary disorders in facial skin of 1204 women in 4 Indian cities.
        Indian J Dermatol Venereol Leprol. 2014; 80: 395
        • Izagirre N.
        • García I.
        • Junquera C.
        • De La Rúa C.
        • Alonso S.
        A scan for signatures of positive selection in candidate loci for skin pigmentation in humans.
        Mol Biol Evo. 2006; 23: 1697
        • Jablonski N.G.
        • Chaplin G.
        The evolution of human skin coloration.
        J Hum Evol. 2000; 39: 57-106
        • Jablonski N.G.
        • Chaplin G.
        Human skin pigmentation as an adaptation to UV radiation.
        Proc Natl Acad Sci USA. 2010; 107: 8962-8968
        • Jaswal I.
        Pigmentary variation in Indian populations.
        Acta Anthropogenet. 1983; 7: 75
      1. Joshua Project, https://joshuaproject.net/; 2016 (accessed 14 April 2016).

        • Jonnalagadda M.
        • Norton H.
        • Ozarkar S.
        • Kulkarni S.
        • Ashma R.
        Association of genetic variants with skin pigmentation phenotype among populations of west Maharashtra, India.
        Am J Hum Biol. 2016; 28: 610-618
        • Jonnalagadda M.
        • Ozarkar S.
        • Ashma R.
        • Kulkarni S.
        Skin pigmentation variation among populations of West Maharashtra, India.
        Am J Hum Biol. 2015; 28: 36-43
        • Karve I.
        Kinship organization in India.
        India:Asia Publishing House Bombay, Bombay1968
        • Lamason R.L.
        • Mohideen M.A.P.K.
        • Mest J.R.
        • Wong A.C.
        • Norton H.L.
        • Aros M.C.
        • et al.
        SLC24A5, a putative cation exchanger, affects pigmentation in zebrafish and humans.
        Science. 2005; 310: 1782-1786
        • Lao O.
        • De Gruijter J.
        • Van Duijn K.
        • Navarro A.
        • Kayser M.
        Signatures of positive selection in genes associated with human skin pigmentation as revealed from analyses of single nucleotide polymorphisms.
        Ann Hum Genet. 2007; 71: 354-369
        • Librado P.
        • Rozas J.
        DnaSP v5: a software for comprehensive analysis of DNA polymorphism data.
        Bioinformatics. 2009; 25: 1451-1452
        • Liu J.Z.
        • Mcrae A.F.
        • Nyholt D.R.
        • Medland S.E.
        • Wray N.R.
        • Brown K.M.
        • et al.
        A versatile gene-based test for genome-wide association studies.
        Am J Hum Genet. 2010; 87: 139-145
        • Mallick C.B.
        • Iliescu F.M.
        • Möls M.
        • Hill S.
        • Tamang R.
        • Chaubey G.
        • et al.
        The light skin allele of SLC24A5 in South Asians and Europeans shares identity by descent.
        PLoS Genet. 2013; 9: e1003912
        • Miller C.T.
        • Beleza S.
        • Pollen A.A.
        • Schluter D.
        • Kittles R.A.
        • Shriver M.D.
        • et al.
        cis-Regulatory changes in Kit ligand expression and parallel evolution of pigmentation in sticklebacks and humans.
        Cell. 2007; 131: 1179-1189
        • Mukherjee M.
        • Mukerjee S.
        • Sarkar-Roy N.
        • Ghosh T.
        • Kalpana D.
        • Sharma A.K.
        Polymorphisms of four pigmentation genes (SLC45A2, SLC24A5, MC1R and TYRP1) among eleven endogamous populations of India.
        J Genet. 2013; 92: 135
        • Rees J.L.
        Genetics of hair and skin color.
        Annu Rev Genet. 2003; 37: 67-90
        • Richardson K.
        • Lai C.-Q.
        • Parnell L.D.
        • Lee Y.-C.
        • Ordovas J.M.
        A genome-wide survey for SNPs altering microRNA seed sites identifies functional candidates in GWAS.
        BMC Genomics. 2011; 12: 504
        • Schwartzberg J.
        Caste regions of the north Indian plains in structure and change in Indian society.
        in: Singer M. Cohn B. Aldine Publishing Co., Chicago1965: 81-113
        • Stephens M.
        • Donnelly P.
        A comparison of Bayesian methods for haplotype reconstruction from population genotype data.
        American Journal of Human Genetics. 2003; 73: 1162-1169
        • Stokowski R.P.
        • Pant P.
        • Dadd T.
        • Fereday A.
        • Hinds D.A.
        • Jarman C.
        • et al.
        A genomewide association study of skin pigmentation in a South Asian population.
        Am J Hum Genet. 2007; 81: 1119-1132
        • Sturm R.A.
        Molecular genetics of human pigmentation diversity.
        Hum Mol Genet. 2009; 18: R9
        • Sturm R.A.
        • Duffy D.L.
        Human pigmentation: painting by numbers or ancestry?.
        Pigment Cell Melanoma Res. 2013; 26: 605-606
        • Thangaraj K.
        • Joshi M.B.
        • Reddy A.G.
        • Gupta N.J.
        • Chakravarty B.
        • Singh L.
        CAG repeat expansion in the androgen receptor gene is not associated with male infertility in Indian populations.
        J Androl. 2002; 23: 815-818
        • Thody A.J.
        • Higgins E.M.
        • Wakamatsu K.
        • Ito S.
        • Burchill S.A.
        • Marks J.M.
        Pheomelanin as well as eumelanin is present in human epidermis.
        J Invest Dermatol. 1991; 97: 340-344
        • Tsetskhladze Z.R.
        • Canfield V.A.
        • Ang K.C.
        • Wentzel S.M.
        • Reid K.P.
        • Berg A.S.
        • et al.
        Functional assessment of human coding mutations affecting skin pigmentation using zebrafish.
        PloS One. 2012; 7: e47398
        • Voight B.F.
        • Kudaravalli S.
        • Wen X.
        • Pritchard J.K.
        A map of recent positive selection in the human genome.
        PLoS Biol. 2006; 4: e72
        • Williams-Blangero S.
        • Blangero J.
        Skin color variation in eastern Nepal.
        Am J Phys Anthropol. 1991; 85: 281-291