A systematic search of the Pubmed database was conducted to identify relevant English-language studies on NS using the search terms “Netherton syndrome OR SPINK5 OR LEKTI”. All papers published up to 4 October 2023, were reviewed to determine whether the reported patients met the criteria for inclusion in our study.

Previous studies on NS have limitations, such as small sample sizes or inadequate inclusion criteria (Nartisa et al., 2023; Utsumi et al., 2020; Hachem et al., 2006; Sarri et al., 2017), which may introduce bias into the results. To address these issues, we established stringent inclusion and exclusion criteria for participant selection. After applying these criteria, a total of 73 studies were included, yielding 162 patients: 2 patients from our center (Shi et al., 2017; Xu et al., 2017) and 160 patients reported in the literature. Detailed information on all participants is provided in Supplementary Table S1.

The inclusion criteria were as follows: (1) patients presenting with one or more symptoms of NS; (2) patients with confirmed homozygous or compound heterozygous mutations in the SPINK5 gene.

The exclusion criteria were as follows: (1) incorrect assessment of gene mutations, such as cDNA changes not conforming to predicted amino acid alterations or reported cDNA base not matching the corresponding base in the reference sequence; (2) presence of two variants of SPINK5 in a patient where the pathogenic variant could not be determined based on current evidence; (3) identification of only one variant of SPINK5; (4) duplicate records of the same patient across multiple publications; (5) cases involving large deletions or duplications which could not be assigned to certain domains.

SPINK5 gene variants were validated using MutationTaster (https://www.mutationtaster.org/) (Schwarz et al., 2014). All mutation data, including cDNA positions and base changes, were standardized according to transcript number NM_006846.4. The transcript was divided into two regions: the 5′ half (1-1598 bp) and 3′ half (1,599–3195 bp). LEKTI consisted of 15 domains, as follows: Domain-1 (AA 28-66), Domain-2 (AA 91-153), Domain-3 (AA 155-216), Domain-4 (AA 219-285), Domain-5 (AA 291-352), Domain-6 (AA 361-423), Domain-7 (AA 431-489), Domain-8 (AA 490-551), Domain −9 (AA 561-622), Domain-10 (AA 626-688), Domain-11 (AA 701 - 757), Domain-12 (AA 768 - 830), Domain-13 (AA 843-905), Domain-14 (AA 910- 971), Domain-15 (AA 987-1,048) (Consortium, 2023). The domains were grouped into five domain regions based on different bio-reactive fragments: DomainR-1 (Domain-1 to Domain-5), DomainR-2 (Domain-6), DomainR-3 (Domain-7), DomainR-4 (Domain-8 to Domain-9), DomainR-5 (Domain-10 to Domain-15) (Sarri et al., 2017). Variants were categorized based on these domains and domain regions for further analysis.

This study was designed as a case-control study. Data analysis was performed using SPSS version 24.0 (IBM Corporation, Chicago, IL, United States) and GraphPad PRISM version 7.01 (GraphPad Software, San Diego, CA, United States). The distribution of domain regions and mutation types was described using frequency distributions and pie charts. Binary logistic regression (Enter method) was employed to evaluate the correlation between genotypes (independent variable) and phenotypes (dependent variable), with age included as a covariate due to the age-dependent nature of certain phenotypes. Commonly recorded phenotypes, including ILC, erythroderma, hair shaft anomalies, failure to thrive, recurrent/systematic infections, hypernatremia, angioedema, urticaria, and asthma, were selected for analysis to ensure adequate sample sizes and statistical reliability.

For the assignment of values to independent variables, variants in the 5′ half and 3′ half of the SPINK5 gene were coded as 1 and 2, respectively. To assess the association between the location of individual domains and phenotypes, variants in different domains were coded according to the domain number. Additionally, variants in DomainR-5 and other domains were assigned values of 1 and 0, respectively. Similarly, variants in domain regions and non-domain regions were assigned values of 1 and 0, respectively. The presence and absence of a phenotype were coded as 1 and 0, respectively. To evaluate the overall effect of two allele variants, we multiplied the values of the two allele variants as independent variables, including Halfx, DomainR-5x, and Domainx.

To investigate genotype-phenotype correlation in patients with homozygous variants, the value of a single allele variant was used, as both alleles were identical. Odds ratios (ORs) with 95% confidence intervals (CIs) were calculated to assess the strength of associations. Contingency tables and Chi-square test were used to investigate the relationship between homozygous condition of fatal mutations and death during infancy. A P-value <0.05 was considered statistically significant.

A total of 162 participants were included in this study. Among them, 100, 142, 134, 73, 95, 56, 35, 44, and 50 patients had clear clinical records of ILC, erythroderma, hair shaft anomalies, failure to thrive, recurrent/systemic infections, hypernatremia, angioedema, urticaria, and asthma, respectively. Detailed data are presented in Table 1. In total, 324 allele variants and 75 different mutations were identified in the 162 patients. Among them, 73 patients carried heterozygous variants, and 89 carried homozygous variants. The distribution of the 324 allele variants in the SPINK5 gene is shown in Figure 1A. The majority of variants (87.04%) were located within domain regions, with the 5′-terminus (DomainR-1) and 3′-terminus (DomainR-5) representing mutation hotspots (Figure 1B). The most common mutation types observed in patients were nonsense, splicing, and deletion mutation (Figure 1C). Detailed genotype information is provided in Table 2.

A previous study suggested that variants located upstream of the SPINK5 gene are associated with more severe phenotypes (Sarri et al., 2017). To further investigate this, we defined 5′ half as upstream and 3′ half as downstream. Our analysis confirmed this trend, showing that patients with variants located in the 5′ half of the gene tended to experience failure to thrive [P < 0.05, OR (95% CI): 0.605 (0.384–0.954); Figure 2A]. Furthermore, a detailed analysis of the correlation between phenotypes and variants in different domains revealed that patients with variants in upstream domains were more likely to experience failure to thrive [P < 0.05, OR (95%): 0.987 (0.977–0.988)]. When analyzing domain regions, variants outside of DomainR-5 were more strongly associated with failure to thrive [P < 0.05, OR (95%): 0.250 (0.066–0.942); Figure 2B], suggesting that mutations in specific domain regions have a significant impact on clinical phenotypes. Additionally, patients with variants in domains were more likely to have ILC [P < 0.01, OR (95%): 6.385 (1.722–23.679); Figure 2C]. No significant correlations were found between individual domains and other phenotypes; detailed information is provided in Supplementary Table S2.

To further validate these findings, we performed the same analysis in patients with homozygous variants. The results revealed that patients with homozygous variants located in the 5′ half of the gene [P < 0.05, OR (95% CI): 0.121 (0.020–0.720); Figure 3A] or upstream domains [P < 0.05, OR (95%): 0.769 (0.618–0.956)] or outside of DomainR-5 [P < 0.05, OR (95%): 0.121 (0.020–0.720); Figure 3B] were more likely to experience failure to thrive. However, although a tendency was observed, no significant correlation was found between ILC and whether the variants were located in domains [P > 0.05, OR (95% CI): 8.339 (0.573–131.243); Figure 3C]. The lack of correlation may be attributed to the small number of patients with variants outside domains (N = 3). No significant correlations were found between individual domains and other phenotypes; detailed information is provided in Supplementary Table S3.

The most severe phenotype observed in this study was death during infancy. Fourteen patients died during infancy, all of whom carried homozygous fatal variants, including c.153delT, c.1431-12G>A, c.1111C>T, c. 1887 + 1G>A, and c. 995delT. Notably, most of these fatal variants were located at the 5′ terminus of DomainR-5.

Furthermore, 18 patients with homozygous fatal variants and 12 patients with heterozygous fatal variants were identified. None of the patients with heterozygous fatal variants died during infancy, and this difference was statistically significant (P < 0.001). These findings suggest that homozygous fatal variants are strongly associated with premature death. Detailed data are presented in Table 3.