If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Division of Molecular and Forensic Genetics, Department of Forensic Medicine, Nicolaus Copernicus University in Toruń, Collegium Medicum in Bydgoszcz, Skłodowskiej-Curie St. 9, 85-094 Bydgoszcz, Poland
Ephelides are hyperpigmented spots on skin observed in Europeans and Asians.
Pigmentation genes explain significant portion of freckles heritability.
Binomial and multinomial models for prediction of freckles are developed.
Non- medium- and heavy-freckled people predicted with AUC = 0.75, 0.66 and 0.79.
Several epistatic effects and sex important for freckles prediction.
Freckles or ephelides are hyperpigmented spots observed on skin surface mainly in European and Asian populations. Easy recognition and external visibility make prediction of ephelides, the potentially useful target in the field of forensic DNA phenotyping. Prediction of freckles would be a step forward in sketching the physical appearance of unknown perpetrators or decomposed cadavers for the forensic DNA intelligence purposes. Freckles are especially common in people with pale skin and red hair and therefore it is expected that predisposition to freckles may partially share the genetic background with other pigmentation traits. The first proposed freckle prediction model was developed based on investigation that involved variation of MC1R and 8 SNPs from 7 genes in a Spanish cohort [
]. In this study we examined 113 DNA variants from 46 genes previously associated with human pigmentation traits and assessed their impact on freckles presence in a group of 960 individuals from Poland. Nineteen DNA variants revealed associations with the freckle phenotype and the study also revealed that females have ∼1.8 higher odds of freckles presence comparing to males (p-value = 9.5 × 10−5). Two alternative prediction models were developed using regression methods. A simplified binomial 12-variable model predicts the presence of ephelides with cross-validated AUC = 0.752. A multinomial 14-variable model predicts one of three categories – non-freckled, medium freckled and heavily freckled. The two extreme categories, non-freckled and heavily freckled were predicted with moderately high accuracy of cross-validated AUC = 0.754 and 0.792, respectively. Prediction accuracy of the intermediate category was lower, AUC = 0.657. The study presents novel DNA models for prediction of freckles that can be used in forensic investigations and emphasizes significance of pigmentation genes and sex in predictive DNA analysis of freckles.
The presence of ephelides is a very distinctive pigmentation feature, which is observed in Asian and European populations and is particularly common among individuals with lighter pigmentation. Freckles show inter-individual variation in terms of the number, from few to hundreds as well as in colour, from red to light or dark brown [
]. Typically, ephelides appear in childhood, suggesting mainly genetic determination of their appearance. Although the number of freckles can increase in adolescence, they may also disappear with age [
Ephelides represent important target for the emerging field of forensic DNA phenotyping (FDP). As with the other pigmentation traits and externally visible characteristics in general, prediction of freckling can be valuable for forensic DNA intelligence. FDP is a relatively new area of forensic genetics that explores DNA variation in terms of its potential to predict biogeographic ancestry, physical appearance and age [
]. Prediction of freckles may complement the description of pigmentation in unknown individuals and this might be useful in both criminal investigations and identification of human remains. Until now, only one prediction model for freckles presence has been proposed in forensic genetics based on investigation of MC1R variation and additional 8 SNPs from 7 pigmentation genes in a Spanish cohort. The discovery cohort involved people with freckles reported in the childhood who could lose ephelides during the sample collection in the adulthood. The developed binomial regression model includes 5 genetic predictors in 4 pigmentation genes (MC1R R alleles, MC1R r alleles, rs12203592 in IRF4, rs2153271 in BNC2, and rs4911442 in ASIP), and information about sex, which was found relevant for prediction of ephelides in women. Validation of that model on an independent group of Spanish adults showed the model accuracy at AUC = 0.756 and the proportion of correctly predicted individuals equalled 69.6%. The accuracy values were even better when information on the presence of freckles was taken from childhood (AUC = 0.809), showing that the model could better detect the presence of the childhood phenotype [
]. In order to take a step forward in predicting the presence of freckles, we conducted a study covering a wider spectrum of candidate pigmentation loci, and because freckles can have varying degrees of intensity, we have sought to differentiate between the low and significant abundance of ephelides. In this work we report results of association testing and prediction modelling that involved 113 candidate markers analysed in a group of 960 adult individuals of Polish origin. Two models for freckle prediction were developed based on categorisation of freckle phenotype as a binomial trait or with assumption of intermediate state with moderate freckle presence. These models may find direct application in forensic genetics and anthropology.
2. Materials and methods
2.1 The study cohort and phenotyping regime
The study cohort consisted of 960 samples: 290 females (30.2%) and 670 males (69.8%), collected in the project NEXT [grant number DOB-BIO7/17/01/2015]. The study was approved by the Ethics Committee of the Jagiellonian University in Krakow (decision no. KBET/122/6120/11/2016) and all donors gave informed consent. The investigated group included 617 freckled individuals and 343 non-freckled ones. Facial freckles presence was assessed by a combination of self-declared status in the questionnaire and by evaluation of high-quality photographs of upper part of the face in the eye area, collected during the phenotype investigation. Two categories of freckles: 1) few (2–3) to several freckles (n = 501); and 2) numerous (>50–60) freckles (n = 116) were identified, based on the thorough analysis of the photographs (Supplementary Figure S1). Dependence tests were performed for the freckling prevalence between males and females as well as light and dark-skinned individuals using Pearson’s χ2 test. Higher prevalence of freckles in females comparing to males was observed (73.5% vs. 60.3%; p-value = 9.5 × 10−5). Skin colour was assigned according to the Fitzpatrick scale by an experienced dermatologist (AB). Light skin colour was predominant in the studied cohort with 65.0% of the samples assigned to category I or II of the Fitzpatrick scale. The difference in freckles presence between lighter and darker skin colour individuals was significant between extreme phenotypes (non-freckled vs. highly freckled: p-value = 0.001), but not for the entire cohort (p-value = 0.262). Additional information from the questionnaire included in this study involved age and sex. Characteristics of the studied sample set are summarized in Table 1.
], for an overview). Altogether, 113 DNA candidates previously associated with pigmentation traits (eye, hair or skin colour) were subjected to analysis. They were tested for agreement with Hardy-Weinberg equilibrium and linkage disequilibrium using PLINK 1.07 [
]. Deviation from Hardy-Weinberg equilibrium was noted for rs4936890 located in the intergenic region between OR10G7 and OR10DP5 genes on chromosome 11 (Supplementary Table S1). Strong linkage disequilibrium (r2>0.8) concerned 29 SNP pairs (Supplementary Table S2), among which 12 pairs referred to SNPs in HERC2 gene and 9 pairs - SNPs in OCA2 gene. One pair (namely, rs12931267 and rs1805007) represent variants from different, but neighbouring genes (FANCA and MC1R, respectively). For further analyses polymorphic DNA variants within the MC1R gene were divided into two subgroups and assigned according to their penetrance as high-penetrance MC1R R (rs1805006, rs11547464, rs1805007, rs1805008, and rs1805009) and low-penetrance MC1R r (rs1805005, rs2228479, rs1110400 and rs885479), as applied in many previous papers [
]. Nine rare SNP markers were excluded at the level of statistical analyses (monomorphic in our dataset: rs1800414 in OCA2 and rs3212355, rs312262906, rs201326893 in MC1R gene or minor allele frequency < 1%: rs74653330, rs121918166, rs1448484 in OCA2, rs6119471 in ASIP and rs1426654 in SLC24A5 genes). Ninety-seven markers analysed in the study are listed in Supplementary Table S3.
2.3 Collection of genetic data
Whole blood was collected from the volunteers and subjected to DNA extraction using PrepFiler Express™ Forensic DNA Extraction Kit (Thermo Fisher Scientific) according to the manufacturer's protocol. Genotypic data were collected using targeted next-generation sequencing with Ion AmpliSeq™ technology and Ion S5™ Platform. Ion AmpliSeq™ custom, single pool panel comprising 294 SNPs covered by 287 amplicons was designed using Ion AmpliSeq™ Designer tool (https://www.ampliseq.com) with support of Thermo Fisher Scientific team. DNA libraries were prepared using the Ion AmpliSeq™ Library Kit 2.0 according to the manufacturer’s protocol but a modified PCR reaction volume of 10 μl and DNA input of 5–10 ng. DNA libraries were quantified using the Agilent High Sensitivity DNA Kit (Agilent Technologies) or Qubit dsDNA High-Sensitivity Assay Kit (Thermo Fisher Scientific) and normalized to 40 pM. Then, 48 libraries were combined in equal ratios and subjected to the template preparation reaction with the Ion 520™ & 530™ Kit-Chef and the Ion Chef System. Sequencing was conducted with the Ion S5™ platform using Ion 530™ Chips. Raw data was analysed with Torrent Server default settings and variants were called with Variant Caller (v126.96.36.199) and HID_SNP_Genotyper (v5.2.2) plugins.
2.4 Single-marker association testing and analysis of epistasis
DNA variants as well as age, and sex were tested for association with freckles on the entire dataset of 960 samples, using univariate binomial logistic regression analysis. Allelic odds ratios (ORs) with 95% confidence intervals (CIs) and respective p-values were determined for minor alleles categorized in an additive manner. The proportion of total variance in freckles presence explained by the tested variable was estimated using the Nagelkerke pseudo-R2 statistic. P-value < 0.05 was considered as nominally significant. Association analyses were carried out with PS IMAGO 4 software (IBM SPSS Statistics 24). SNP-SNP interactions were also tested based on the entire dataset using Multifactor Dimensionality Reduction which is the data mining strategy used for detection and characterisation of non-linear relationships between variables [MDR 3.0.2; http://multifactordimensionalityreduction.org; [
]]. Furthermore, we used multivariate logistic regression for simultaneous testing all possible pairwise interactions between SNPs and identification of those that expressed independent effects on freckles prediction. This analysis was carried out with PS IMAGO 4 software (IBM SPSS Statistics 24).
2.5 Prediction modelling
Two different models were developed for predicting the presence of freckles in a total of 960 individuals. The binomial logistic regression model assumed simplistic classification of the subjects as freckled or non-freckled. The multinomial logistic regression model was developed for the study participants classified as heavily freckled, medium freckled or non-freckled. Multivariate logistic regression was applied for simultaneous analysis of all the SNPs including all possible pairwise interactions. The minimal sets of DNA predictors for optimizing the prediction accuracy were chosen using forward step-wise regression analysis and selection of markers based on the likelihood-ratio test. This approach provides optimal selection of a set of markers keeping balance between the accuracy of prediction and the number of predictors and protects against the model overfitting. In order to assess the contribution of markers to the final models’ prediction accuracy, they were ranked according to their importance using -2log likelihood of reduced model statistic. The cumulative contribution of each variable to the final prediction accuracy of the models was evaluated on the whole dataset by the area under ROC curve (AUC).
Since the use of the whole dataset may produce overestimated prediction accuracy, cross-validation approach was used to evaluate performance of predictors included in the models. The data was partitioned into 10 equal-sized parts and in each of 10 rounds of the process, 9/10 of the data (n = 864) served as a training set while the remaining 1/10 of the data (n = 96) as a testing set. AUC was calculated based on the combined classifier scoring outputs determined for each testing set. Sensitivities, specificities, positive predictive values (PPV) and negative predictive values (NPV) were calculated based on confusion matrices with freckle phenotype prediction provided by the cross-validated probability of >0.5 when binomial categorisation was applied (freckled, non-freckled) or the highest probability approach in case of multinomial categorisation (heavily freckled, medium freckled and non-freckled). Additionally, alternative thresholds of probability for freckles prediction were tested ranging from p > 0.5 to p > 0.9 with 0.05 interval. Prediction modelling was carried out using PS IMAGO 4 software (IBM SPSS Statistics 24).
The association testing revealed 19 DNA variants nominally associated with the freckle phenotype and in the studied cohort females had ˜1.8 higher odds of freckles presence comparing to males (p-value = 9.5 × 10−5). MDR analysis indicated a synergistic interaction between rs1042602 in TYR and rs9894429 in NPLOC4 important for freckles prediction, which was coherently confirmed by regression analysis. The details of association testing and gene-gene interaction analyses are presented in supplementary materials (Supplementary Table S3, Table S4 and Figure S2).
Two alternative prediction models were developed using regression methods. Binomial logistic regression, assuming individuals classified as freckled and non-freckled identified 12 independently contributing predictors, including 6 main effects of SNPs in IRF4, DEF8, NPLOC4, NCOA6, HERC2 and OCA2 genes, pooled MC1R R variants, 4 SNP-SNP interactions and sex. Fig. 1 and Table 2 shows that among the three strongest predictors of freckles was the interaction between rs1042602 in TYR and rs9894429 in NPLOC4, consistently identified by MDR and regression, as well as main effects of rs9894429 in NPLOC4 and rs12203592 in IRF4. These predictors altogether provide prediction accuracy of AUC = 0.660 and explain 12.0% of the total variation in freckling. Sex ranked at 4th position was found to explain additional 2.1% variation and improved prediction accuracy by ΔAUC = 0.023, while rs4911442 in NCOA6 gene, ranked 5th, explained additional 1.4% of variation and increased prediction accuracy of the model to AUC = 0.694. Substantial improvement of prediction accuracy was noted after including variables located on chromosome 16, ranked at 6th – 9th positions: rs8049897 in DEF8 and MC1R represented by MC1R R polled variants, interaction between MC1R R and r variants and interaction of MC1R r with rs10874518 in OLFM3 (chr1), which altogether explained further 9.3% of total variation in freckling and rose up the prediction accuracy by ΔAUC = 0.061. Other predictors, including rs7183877 in HERC2, rs1448485 in OCA2 and interaction between rs7183877 in HERC2 and rs4911442 in NCOA6 taken together improved prediction accuracy by additional ΔAUC = 0.013 and explained 2.6% of variation.
Table 2Parameters of the BLR and MLR prediction models.
rs9894429 × rs1042602
NPLOC4 × TYR
T × A
MC1R r × rs10874518
MC1R × OLFM3
r × C
MC1R R × MC1R r
R × r
rs7183877 × rs4911442
HERC2 × NCOA6
A × G
rs12350739 × rs1885120
BNC2 × MYH7B
G × C
Sex × rs6059655
Sex / RALY
male × A
BLR – binomial logistic regression, MLR – multinomial logistic regression, b1, b2 – betas for freckled and heavily freckled phenotypes, respectively, non-freckled category as a reference; ΔR2 - Nagelkerke pseudo R2 gain while the current variant is included in the model. Sex was coded as follows: male - 1, female - 0.
The overall prediction accuracy of the binomial model as expressed by AUC was 0.768 while after 10-fold cross-validation this value was slightly lower – 0.752. The total variation in freckles presence explained by this model accounted for 27.4%. The cumulative contribution of each variable to the final prediction accuracy and the explained phenotypic variance is presented in Fig. 1 and Table 2. Sensitivity of the model reached 84.0% and specificity of 46.4%, while PPV and NPV accounted for 73.8% and 61.6%, respectively. Parameters describing prediction accuracy of the model are shown in Table 3. Alternative thresholds of probability for freckles prediction tests revealed that increasing threshold to p > 0.6 caused raising of all the prediction model parameters, i.e. sensitivity up to 86.2%, specificity to 54.6%, PPV to 80.4% and NPV to 64.8%, however further increase in the threshold (p > 0.65 and larger) negatively affected the specificity value. Using the threshold p > 0.6, 24.7% of cases remained inconclusive, while the error rate decreased from 29.5% to 23.8%. Results obtained for the all tested thresholds are presented in Supplementary Table S5.
Table 3Accuracy estimations calculated for the studied cohort using the developed prediction models and the Hernando’ predictors [
Multinomial logistic regression assuming classification of individuals as heavily freckled, medium freckled and non-freckled identified 14 significantly contributing predictors. Ten of them were overlapping with binomial regression model. The four additional predictors included rs6059655 in RALY, MC1R r variants and two interactions: rs12350739 in BNC2 with rs1885120 in MYH7B and sex with rs6059655 in RALY. Similar to the binomial model, interaction between rs1042602 in TYR and rs9894429 in NPLOC4, rs12203592 in IRF4 and rs9894429 in NPLOC4 were most important for freckles prediction, explaining 12.7% of the total variance in freckles presence. The interaction had the highest impact on non-freckled category (ΔAUCnon-freckled = 0.072), while two other predictors affected extreme freckle categories (rs12203592 ΔAUCheavily freckled = 0.089, ΔAUCnon-freckled = 0.054 and rs9894429 ΔAUCheavily freckled = 0.029, ΔAUCnon-freckled = 0.033). Variation in the MC1R gene, included in the model as pooled MC1R R variants (4th), MC1R r variants (12th) and two interactions: MC1R R with MC1R r (11th) and MC1R r with rs10874518 (OLFM3) (8th), altogether explained additional 8.1% of variance and the impact of MC1R R on prediction accuracy was particularly strong (ΔR2 = 4.1%) in the category of heavily freckled individuals (ΔAUCheavily freckled = 0.051). Interaction between rs12350739 in BNC2 and rs1885120 in MYH7B was ranked 5th and explained 1.1% of total variation. The multinomial model also included sex, which increased the total variation explained by 2.0% boosting mostly prediction accuracy of medium freckled individuals (ΔAUCmedium freckled = 0.022). The category of non-freckled individuals was highly influenced by rs8049897 in DEF8, which explained further 3.6% of the phenotype variation (ΔAUCnon-freckled = 0.026). Finally, rs6059655 in RALY and interaction between this DNA variant and sex, were ranked at 9th and 13th positions. rs6059655 had strong impact on the category of heavy freckled individuals (ΔAUCheavily freckled = 0.023) while in interaction with sex slightly increased prediction accuracy of all categories (ΔAUC = 0.003-0.004). Altogether the main effect of rs6059655 and in interaction with sex explain 3.4% of the total variance in freckling. Common for both predictive models were DNA variants rs7183877 in HERC2 and rs1448485 in OCA2 that explained 2.1% of the variation with impact mostly on intermediate category (ΔAUC = 0.012 and ΔAUC = 0.008, respectively).
The overall prediction accuracy of the multinomial model reached the value of AUC = 0.814 for heavily freckled, AUC = 0.686 for medium freckled and AUC = 0.771 for non-freckled individuals and after 10-fold cross-validation these values were 0.792, 0.657 and 0.754, respectively. In addition, the total variation explained by this model was 33.0% (Fig. 1). The highest sensitivity of the model was observed for medium freckled individuals (74.3%), being lower for non-freckled (54.5%), and very low for heavily freckled individuals (17.2%), however on the contrary, the highest specificity of freckles prediction was recorded for heavily freckled individuals – 97.8%, 79.4% for non-freckled and 48.8% for medium freckled. Parameters describing prediction accuracy estimated for the multinomial model are shown in Table 3.
Like for binomial model, prediction accuracy parameters improved with the increase of the probability threshold. The threshold of p > 0.6 resulted in substantially higher prediction sensitivity for all three categories, with the highest increase for medium freckled individuals (10.2%) and positively influenced PPV values. Further increasing the threshold had positive effect on accuracy but also significantly increased the proportion of inconclusive results. Detailed results for all tested thresholds are presented in Supplementary Table S5.
In the next step we assessed the accuracy of predictors proposed in the Spanish study in our cohort [
]. The Hernando’ model consists of 5 genetic variables and information about sex. These predictive markers were used for model building using the current dataset and returned AUC = 0.667 after 10-fold cross-validation. Sensitivity of the prediction reached 85.6% and specificity 32.1%. Positive and negative predictive values equalled 69.3% and 55.3%, respectively. Details of this analysis are presented in Table 3.
Ephelides are easily observed externally visible traits and although freckling is a very heritable trait that additive genetic effects have been estimated to account for 91% of the total variance in freckle counts [
], their biology and genetics are still poorly understood. Because freckles are pigmentation spots and are mainly observed in people with pale skin, the association with genes involved in pigmentation seems obvious. The first predictive DNA model for ephelides was developed based on variation in 4 pigmentation genes adjusted for sex and the information about freckles presence in childhood. It allows for prediction of the freckle’s presence in the independent group of a Spanish adults with AUC = 0.756 [
], while the same set of predictors in the Polish cohort of 960 adults reached the cross-validated prediction accuracy of AUC = 0.667. Here we make a step forward in prediction of freckles by conducting prediction modelling including DNA variation from additional 42 pigmentation genes and by classification of freckle phenotype into 3 categories – non-freckled, freckled and heavily freckled. Nineteen DNA variants showed nominal significant association with freckles. The study confirmed the previous information about significantly higher prevalence of freckles presence in females than in males (p-value = 9.5 × 10−5) [
]. Details of association testing are presented in supplementary materials.
The main objective of this study was to assess whether the currently available markers associated with pigmentation can predict freckles with acceptable accuracy. We experimented with several prediction models including Random Forests, CART Tree and Extreme Gradient Boosting but since their outcome did not exceed regression methods (data not shown) and at the same time they are much more complex and difficult to interpret than logistic regression, we decided to use logistic regression to develop the final models. In order to increase statistical power of detection of significant effects, the models have been fitted using the entire set of 960 samples.
Our simplistic binomial model for freckles presence assuming the studied individuals to be classified in two categories – freckled vs. non-freckled includes 12 predictors (Table 2). Among them are 4 predictors overlapping with the Hernando model and these include MC1R variants, rs12203592 in IRF4, rs4911442 in NCOA6, and information about sex [
]. All of them revealed nominal significance (p < 0.05) in the univariate logistic analyses, confirming significant impact on freckles formation in our independent study. MC1R and IRF4 are important pigmentation markers included in the HIrisPlex-S system for eye, hair and skin colour prediction [
]. Surprisingly, MC1R previously reported as a major freckle gene had a medium impact on freckling prediction in our cohort. Detailed analysis showed that 88.9% individuals with R/R genotype had freckles and this high number dropped to 67.3% for R/0 genotype carriers. These results are in line with the red hair colour phenotype distribution in carriers of MC1R variants. In the Polish population R/R genotype was observed in 98.3% red (or red-blond) while R/0 in 15.4% individuals and other phenotypic effects (like red facial hair) of single R or even r variants could have been noticed in some individuals [
]. We speculate that the influence of MC1R on freckles presence in our study may be slightly underestimated as we analysed only 11 variants with 9 of them showing variation being included in the final statistical analyses. It seems that analysis of full MC1R sequence may be important in order to achieve good accuracy of prediction. Hernando sequenced the entire MC1R gene and the study revealed 18 R variants and 13 r variants and thus the information gained from MC1R analysis could be higher compared to our study. The binomial model further includes OCA2, HERC2 and DEF8, which showed main association effects in freckling. These 3 genes have been linked to skin colour determination. Notably, 11 DNA variants of these genes were included in the HIrisPlex-S skin colour prediction model [
]. It is also important to note that the DEF8 gene located in the neighbourhood of MC1R provided independent information in multivariate testing and thus variants from both genes were included in our final model for freckles presence. The last main effect freckles gene included in our binomial model was NPLOC4 that previously has been associated with quantitative eye colour [
]. Finally, 4 SNP-SNP interaction effects were selected for the binomial model confirming importance of genetic interactions in prediction analyses.
Total prediction accuracy obtained for the binomial model was AUC = 0.752. Although sensitivity of prediction achieved relatively high value (84.0%), resulted from 518 individuals correctly predicted as having freckles out of 617 freckled individuals, the specificity was significantly lower (46.4%) – out of 343 non-freckled individuals, 184 were predicted as freckled, showing high rate of false positive predictions rate (Table 3). This result may have biological explanation, as it is known that freckle phenotype in a childhood may disappear in adulthood while the studied group comprised adult individuals (mean age = 30.7, SD = 9.0). This hypothesis also prompts further experimental work but if confirmed, provides additional argument for including age in prediction modelling of appearance traits [
]. We also used the Hernando markers to develop a model based on our dataset and to make predictions. This model revealed lower accuracy in our cohort (AUC = 0.667 after 10-fold cross-validation) compared to our binomial model. Sensitivity of the prediction was comparable with sensitivity obtained for our model (85.6%), but specificity was much lower (32.1%), showing higher false positive rate. It might be considered whether differences in the performance of both models may result from differences in the model building procedures. The Spanish training set included individuals that reported freckles presence in the youth while our training set comprised individuals with freckling recorded in adulthood what may be important because of the above-mentioned phenomenon of disappearance of freckles with age. Indeed, the Spanish model showed slightly decreased accuracy when tested on an independent set of adults compared to the results obtained for the cohort with information on freckles collected in childhood (AUC decrease from 0.809 to 0.756). Although sensitivity remains roughly at the same level (59.8% vs. 58.9%), specificity slightly decreased (89.1% vs. 80.7%) revealing higher false positive rate in adults, which may be a result of disappearance of freckles with age. Additionally, individual traits may have different genetic bases in populations of various ancestries [
], but also closely related populations may show differences in allele frequencies that may affect the significance of certain predictors and, consequently, affect prediction results. Therefore, further studies on a universal freckle prediction model should involve sample sets from various populations in Europe and Asia. Total variation in freckles presence explained by our binomial model was 27.4%, thus lower than reported by Hernando in Spanish population (32.9%), (Fig. 1).
The multinomial model explained more of the total variation in freckle phenotype in this population compared to the binomial model (33.0% vs. 27.4%). We suggest that the multinomial model makes a step forward in prediction of freckles presence. The model predicts individuals with few freckles and those highly freckled based on 14 variables, 10 of which are shared with the binomial model (Table 2). Similar to the binomial model, the most significant predictor of freckles was the epistatic effect between rs1042602 in TYR and rs9894429 within NPLOC4. Importantly, sex was also found to be significant and valuable particularly for prediction of intermediate category. This result evokes situation known from eye colour prediction. Implementation of the information about sex to the prediction model have improved intermediate eye colour prediction in Polish and Spanish populations [
]. It is still unclear how sex can affect pigmentation phenotype, including freckled phenotype, but as the data show this effect can be differentiated in terms of power and population specificity. Although impact of hormones such as estrogen and progesterone has been shown to be a possible factor increasing the risk of pigmented spots in women [
], it is barely possible to have something to do with freckles which arise in childhood, contrary to facial spots appearing as a sign of aging.
Additionally, similar to binomial model, variants located in important pigmentation genes such as IRF4, MC1R, DEF8, HERC2 and OCA2 revealed significant involvement to 3-category freckles prediction, demonstrating again high contribution of pigmentation genes to freckles genetic determination. Using the multinomial model individuals from two extreme categories, non-freckled and heavily freckled were predicted with good accuracies (AUC = 0.754 and 0.792, respectively) while medium freckled with AUC = 0.657. Non-freckled category was predicted with moderately high sensitivity at the level of 54.5% (as a result of 187 out of 343 individuals correctly predicted as non-freckled) but higher specificity – 79.4% (490 out of 617 correctly predicted as having freckles). Intermediate category predictions were more sensitive – 74.3%, however less specific – 48.8%, which is the result of 235 out of 459 non-freckled individuals, incorrectly predicted as medium freckled. However, 84 (35.7%) of these incorrectly predicted individuals were heavily freckled, showing that in these cases general freckles presence was correctly predicted. Lower specificity can be the result of above-mentioned biological phenomenon of freckles disappearing in adulthood. In case of heavily freckled individuals, specificity of prediction provided by the model was very high (97.8%). Nevertheless, only 20 from 116 individuals heavily freckled were predicted as such, thus sensitivity of recognizing highly freckled individuals by this model is low – 17.2%. On the other hand, 84 from 116 heavily freckled individuals (72.4%) mentioned above, were predicted as medium freckled, showing that even when intense freckling is not recognized by the prediction model, medium state of freckles presence will be predicted. Only 10.3% of heavy freckled persons were incorrectly predicted as non-freckled. Although non-freckled and heavily freckled phenotypes may be predicted with fairly good accuracy, the prediction of the intermediate category prompts further research to provide better insight into the genetics of freckles. The phenomenon of more accurate prediction of extreme categories is already known from other pigmentation traits studies. Both eye and skin colour prediction modelling demonstrated that the highest prediction accuracy has been obtained for extreme phenotypes: blue and brown eye colour (AUC = 0.91 and AUC = 0.93, respectively) as well as white and black skin colour (AUC = 0.97 and AUC = 0.96, respectively). At the same time, significantly lower prediction accuracy was obtained for intermediate eye (AUC = 0.73) and skin (AUC = 0.83) shades [
]. Our results are in line with these findings and confirm that intermediate categories are perhaps the most heterogeneous and most difficult in terms of accurate prediction. As expected, because freckles presence similar to hair colour can change with age, prediction accuracies of these traits are comparably lower (AUC = 0.75 and AUC = 0.78 for blond and black hair colour after 10-fold cross-validation, respectively) [
], likewise non-freckled and heavily freckled phenotype (AUC = 0.75 and AUC = 0.79).
We used cross-validation approach for evaluation of performance of the selected predictors while it has been suggested that this method may overestimate accuracy values. In order to verify how this approach could affect the accuracy measures, the whole cohort was randomly divided into 70% of training set and 30% of testing set. The models developed based on the subset of 70% individuals contained similar sets of predictors compared to the whole-data models. The most important predictors including IRF4, DEF8, TYR, MC1R, HERC2, OCA2, OLFM3, and sex overlapped between the two prediction approaches. Differences were observed for the genetic loci located in the 20q11.22 region, rich in variants associated with skin colour in Europeans [
]. rs4911442 in NCOA6 in binomial model and rs1885120 in MYH7B in multinomial model were replaced by the rs2378249 in PIGU, located 137 kb from rs4911442 and 359 kb from rs1885120 in LD (r2 = 0.5) with rs4911442. Moreover, 70:30 approach identified additionally KITLG to be slightly significant in the binomial model (p = 0.037). The new models showed similar accuracy values obtained for the 30% testing set compared to the models developed based on the whole dataset. For the binomial model AUC=0.74 and in case of multinomial model these values were: AUCnon-freckled = 0.73, AUCmedium freckled = 0.59 and AUChevily freckled = 0.77. One of the obvious reasons for the lower prediction accuracy is the number of samples included in the training set which contains almost three hundred samples less compared to the model developed based on the whole dataset. Both approaches (model building based on 70% and 100% of data) have advantages and disadvantages, however it seems that there is no perfect approach. The advantage of 100% of training data provides the greater power to detect significant predictors. On the other hand, it is not possible to objectively assess the performance of such model because of the possible overestimation of internal validation results. Details for the models developed based on 70% of the data and the prediction accuracies are presented in Supplementary Table S6 and Figure S3.
5. Forensic relevance and future prospects
We found that the presence of freckles can be predicted from DNA but with lower accuracy compared to other pigmentation traits. We assessed that the multinomial model correctly classifies an unknown individual to one of the 3 categories in 60.3% cases while considering medium freckled and heavily-freckled together this classification is correct in 70.5% cases. Prediction accuracy is very important aspect of forensic DNA intelligence studies as the incorrect prediction may mislead investigation. The number of correct predictions can be improved by setting a threshold when interpreting prediction probabilities. Based on our data, we propose p > 0.6 as the optimal threshold, which allow to increase prediction accuracy parameters (only slight drop of NPV value was noted for heavily freckled category: -0.93) and decrease the error rate (Supplementary Table S5). This solution apparently makes predictions more reliable, as the prediction error rate is decreasing what is important for the forensic investigations. Taking into account this threshold, the number of correct predictions was elevated up to 67.6% (when 3 categories were predicted) and 78.8% (when freckled categories were considered together). The developed model for freckles presence can supplement the currently available forensic DNA tests for prediction of pigmentation and enable better description of unknown perpetrators in criminal cases and identification of human remains. Notably, among 14 predictors included in the multinomial model, 5 predictors including rs12203592 in IRF4, rs1042602 in TYR, rs6059655 in RALY, MC1R R and MC1R r are present in the HIrisPlex-S models and additional 4 (rs1448485 in OCA2, rs7183877 in HERC2 and rs8049897 in DEF8 and rs12350739 in BNC2) represent different variants of genes considered in HIrisPlex-S. This means that it may be relatively simple to extend the currently available pigmentation markers to enable prediction of freckle presence in forensic genetics. It is clear that the model for the presence of freckles can be further improved and some progress can be made simply by sequencing the entire MC1R gene, which is very variable and contains highly relevant information. Further studies involving whole genome analyses and larger cohorts of individuals shall identify additional genes that explain further portion of variance in freckles presence. Moreover, advanced methods for big genomic data analysis, greater focus on results with lower statistical confidence and consideration of different forms of heterogeneity would improve prediction of freckles and other phenotypes in humans [
The authors wish to thank all volunteers who contributed to this study. This research was supported by the Polish National Centre for Research and Development [grant number DOB-BIO7/17/01/2015 ]. EP and WB are supported by the European Union’s Horizon 2020 research and innovation programme under grant agreement No. 740580 within the framework of the Visible Attributes Through Genomics (VISAGE) Project and Consortium.
Appendix A. Supplementary data
The following is Supplementary data to this article: