Mutation screening of eight genes and comparison of the clinical data in a Chinese cohort with congenital hypothyroidism

Congenital hypothyroidism (CH) is a common neonatal endocrine disorder, characterized by irreversible intellectual disability and short stature if left untreated. It can be divided into thyroid dysgenesis (TD), including athyreosis, ectopy and hypoplasia, and dyshormonogenesis (DH), also referring to gland in situ (GIS), in which patients have eutopic thyroids with normal size or goiter. This study aims to analyze the clinical and genetic data of 375 Chinese CH patients without DUOX2 and thyroid transcription factor (TTF) variants, and to explore the mutation frequencies of the eight genes and the inheritance pattern of CH. Targeted next generation sequencing (NGS) and statistical analysis were performed for mutation screening on eight CH-related genes and the comparison of clinical data in a cohort of 606 Chinese CH patients from Henan Province. A total of 104 variants were detected in genes required for thyroid formation (TSHR, GLIS3, BOREALIN, NTN1, JAG1 and TUBB1) and thyroid hormone synthesis (TG and TPO) in 83 subjects. Monogenic variants were the most prevalent with a percentage of 75.00% (78/104) followed by oligogenic variants (25.00%, 26/104). No differences were found in various clinical data between patients with and without variants. However, it should be noted that only initial L-T4 dose was statistically different between patients with monogenic variants and oligogenic variants. Our results suggested that apart from Mendelian monogenic inheritance, oligogenic inheritance of CH could not be excluded and also involves other factors, such as penetrance, epigenetic mechanisms and environmental factors.


Introduction
Congenital hypothyroidism (CH) is one of the frequent neonatal endocrine diseases and is characterized by irreversible intellectual disability and physical developmental disorder if left untreated, affecting 1 out of 3000 live births worldwide [1]. The incidence of CH differs from country to country and the estimated incidence is 1:2000-2500 individuals among Chinese populations, higher than that in the world [2]. However, the prevalence of CH has shown an increasing trend in recent years due to the lower thyroid-stimulating hormone (TSH) cut-off value utilized in newborn screening (NBS) strategy in which an increasing number of GIS cases are screened [3]. The female:male ratio of CH is 1.5-2:1 [4].
CH is attributed to either abnormalities in thyroid development (thyroid dysgenesis, TD), including athyreosis, ectopy and hypoplasia, or defective thyroid hormone biosynthesis (dyshormonogenesis, DH), also referring to gland in situ (GIS), in which patients have eutopic thyroids with normal size or goiter [5]. Additionally, CH can be These authors contributed equally: Liangshan Li, Xiaole Li * Fang Wang 18660291711@163.com * Shiguo Liu liushiguo2002@126.com divided into permanent CH (PCH) and transient CH (TCH) depending on the length of clinical course. TD accounts for most cases of PCH, and although the vast majority of TD cases are sporadic, approximately 2% are familial, supporting the involvement of genetic factors in TD [6]. Abnormal thyroid development and thyroid function are more common in the first degree relatives of TD patients in comparison with those in the general population [7]. The formation of thyroid is a complicated process that needs polygenic contributions such as PAX8, GLIS3, NKX2-1/TTF1, FOXE1/TTF2, NKX2-5, TSHR and the latest reported genes TUBB1, CDCA8/BOREALIN, NTN1 and JAG1 [1]. Nevertheless, defects in any of the steps of thyroid hormone synthesis can result in GIS, which is associated with pathogenic variants in TG, TPO, DUOX2, DUOXA2, IYD/DEHAL1, PDS/SLC26A4, NIS/SLC5A5 and SLC26A7 with an autosomal recessive inheritance pattern [8,9]. Although many scholars focused on the study of the relationships between genotype and phenotype, genetic variants in the same gene may be responsible for various CH phenotypes, which indicates that the genotype-phenotype correlations of CH seem to be complex and need to be further established. For example, it is interesting that the known GIS causative gene TPO is also relevant to hypoplasia and similar findings in other CH-related genes have been reported in recent years [10].
Early diagnosis and timely levothyroxine (L-T4) replacement therapy could improve the prognosis of CH children. The affected individuals should be treated with a dose of 10-15 μg/kg/day within the first two weeks of birth [11,12]. However, there is evidence indicating that highdose L-T4 has an impaired effect in cognition when compared with the low-dose L-T4 [13].
In this study, comprehensive mutation screening was performed in a cohort of 606 Chinese CH patients from Henan Province by targeted next generation sequencing (NGS). The gene panel contains eight CH-related causative genes including TG, TPO, TSHR, GLIS3, BOREALIN, NTN1, JAG1 and TUBB1. Genetic screenings for other known pathogenic genes such as DUOX system genes and thyroid transcription factor (TTF) genes have been reported by our team [14,15]. This study aims to analyze the clinical and genetic data of 375 Chinese CH patients without DUOX2 and TTF variants, and to explore the mutation frequencies of the eight genes and the inheritance pattern of CH.

Materials and methods
Patients A total of 375 Chinese CH patients (169 males, 206 females) were enrolled from the Third Affiliated Hospital of Zhengzhou University between January 2012 and May 2018. The diagnostic criteria for CH were an increased TSH level (reference range, 0.47-7.63 mIU/L) and a decreased free thyroxine (FT4; reference range, 11.5-22.74 pmol/L) level in venous blood samples in accordance with TSH value ≥9 μU/mL in the neonatal heel prick samples. The thyroid phenotypes of the subjects were identified on the basis of thyroid ultrasound and radionuclide scanning. Thyroid function was reevaluated after discontinuation of L-T4 for 1 month at 2-3 years of age. If the TSH and FT4 levels were abnormal, the patients were diagnosed with PCH and continued L-T4 replacement therapy. If the TSH and FT4 levels were normal, the patients continued to discontinue L-T4 and were reviewed regularly. Finally, patients with normal thyroid function were diagnosed with TCH. This study was approved by the Ethics Committee of the Third Affiliated Hospital of Zhengzhou University and the Affiliated Hospital of Qingdao University. Whole blood samples were collected from the 606 CH subjects and 100 healthy controls of Chinese Han population (mean age, 6 years; male to female ratio, 1:1.08) for mutation analysis after written informed consent was obtained from the parents or guardians.

Targeted NGS
Genomic DNA was extracted from blood clots using a Genomic DNA extraction kit (Tiangen Biotech Co., Ltd., Beijing, China). All coding exons and exon-intron boundaries of eight CH causative genes TG (NM_003235.

Variant filtering
Raw reads were filtered to obtain clean reads, which were then aligned to the human reference genome (hg19) with the Burrows-Wheeler Aligner (BWA) software. Picard software was utilized to mark duplicate reads. The single nucleotide polymorphisms (SNPs) and insertions and deletions (InDels) were detected and filtered using Genome Analysis Toolkit

Statistical analyses
All statistical analyses were completed by using IBM SPSS Statistics Version 22.0. The measurement data are shown by median and interquartile range (IQR) and the comparison was performed by Mann-Whitney U test, while counting data are expressed in percentage (%) and analyzed by Pearson chi-square (χ2) test. A P value < 0.05 indicated that the difference was statistically significant.

Summary of the 104 variants identified by targeted NGS
Due to the facts that DUOX2 is the most common causal gene for Chinese CH population [2,16] and our study focused on the eight CH-related genes, we excluded patients with variants in DUOX2 and TTF genes, and then analyzed the variants of the eight genes in the remaining 375 patients. A total of 104 variants were found in 83 out of the 375 subjects (22.13%, 83/375), of which the most recurrent variant was c.1349 G > A/p.R450H in TSHR. All the 83 patients were heterozygous for these variants, except one patient who was a homozygote. The majority of these variants occur within the causative genes for TD, including TSHR (28.85%, 30/104), TUBB1 (4.81%, 5/104), GLIS3  Table 1). The pathogenicity of these variants was evaluated based on the ACMG criteria [17]. As shown in Fig. 3, TG was the most frequently mutated gene (32.43%, 24/74), followed by TSHR (20.27%, 15/74), GLIS3 (14.86%, 11/74) and TPO (10.81%, 8/74). Out of the 74 variants, 69 were missense variants, 4 were nonsense variants and 1 was a splicing variant. Missense variants were distributed in seven genes except BOREALIN, whereas nonsense variants were prevalent in TG, followed by TPO. It is noteworthy that the only splicing variant occurred in the BOREALIN gene.  in GIS-related and TD-related genes were associated with individuals accompanied by GIS or TD.

Analyses of oligogenic variants
Eleven subjects (2.93%, 11/375) exhibited oligogenic inheritance pattern and 26 oligogenic variants (25.00%, 26/ 104) were identified in these cases. The number of genes combined with the eight genes were as follows: TG (n = 2), TPO (n = 1), TSHR (n = 4), GLIS3 (n = 3), BOREALIN (n = 0), NTN1 (n = 2), JAG1 (n = 2) and TUBB1 (n = 0). The most prevalent oligogenic variants combination was in TSHR with other genes such as GLIS3, NTN1, JAG1 and TG, which was identified in 7 cases. Patients with GIS or TD harbored the combination of variants in GIS-causative genes and TDcausative genes. However, it is worth noting that the combined variants in genes involving thyroid development were found only in patients with GIS. In the 11 oligogenic cases, 8 had double missense variants in two alleles of two distinct genes and 3 GIS cases with normal-sized glands carried three or four variants in two different genes ( Table 1).

Analyses of phenotype-genotype correlations
The clinical and biochemical data of the 606 subjects (median age: 6 years, ranging from 4 to 18 years) have been summarized in our previous study (Supplementary  Table 2) [15], while the clinical characteristics of the 375 patients without DUOX2 and TTF variants were shown in Table 2. In the 375 CH-affected patients, the proportion was similar between males and females, and was approximately 1:1.22. The number of patients with GIS (n = 252) was significantly larger than those with TD (n = 123). The results of Mann-Whitney U test revealed that no significant differences were found regarding current age, birth weight, age at diagnosis, TSH value at screening, diagnostic serum TSH value, diagnostic FT4 value, initial L-T4 dose, maximum adjusted dose, current dose of L-T4 and treatment duration (All P values > 0.05) between patients with and without variants ( For patients with variants in the eight genes, the differences were not statistically significant in terms of current age, birth weight, age at diagnosis, TSH value at screening, diagnostic serum TSH value, diagnostic FT4 value, maximum adjusted dose, current dose of L-T4 and treatment duration (All P values > 0.05) between patients with monogenic variants and oligogenic variants (Table  3). However, the initial L-T4 dose in the monogenic variants group was higher than that of the oligogenic variants group (P = 0.009). The number of mutated genes seemed to have no statistically significant correlation with the severity of CH (Pearson χ2 = 0.000, P = 1.000). Among the 83 subjects with variants in the eight genes, monogenic variants were observed in 91.67% (22/24) of the TD group and 84.75% (50/59) of the GIS group, but the difference was not statistically significant (Pearson χ2 = 0.236, P = 0.627). Figure 1 shows the distribution of the 104 variants in patients with different thyroid morphology. Twenty-six variants were identified in cases with TD, and 78 variants were found in cases with GIS. Notably, variants in GIS-associated genes (TG and TPO) and TD-associated genes (TSHR, GLIS3, BOREALIN, NTN1, JAG1 and TUBB1) have been detected in both subjects with GIS and TD. The distribution of our CH  cohort based on different number of mutated genes was presented in Fig. 4. The mutation spectrum of the GIS cases showed that the most common mutated genes were TSHR and TG, while the most frequently defective genes in the TD cases were TG, followed by GLIS3.

Discussion
In this study, we aimed to investigate the clinical and molecular findings in a cohort of 375 Chinese CH patients without DUOX2 and TTF variants. A total of 104   [20]. However, the specific mechanism of oligogenic inheritance of CH remains unclear and needs further elucidation.
Our study revealed that the coincidence of variants in GIS and TD genes could be found in both GIS and TD patients. Moreover, oligogenic variants in TD-related genes could present in GIS patients with normal glands or goiter. A previous study showed that one CH proband with a ectopic thyroid had digenic TG and JAG1 variants, another CH proband with GIS was found to be the carrier of digenic variants in PAX8 and GLIS3 [21]. Likewise, the coexistence of a PAX8 variant p.P25L and a TG variant p.R321X was identified in a CH-affected individual with hypoplasia [22]. In addition, one female subject with hypoplasia harboring both a monoallelic TSHR variant p.R450H and a monoallelic DUOX2 variant p.V779M was described by Abe et al. in 2018 [23], whereas Sasivari et al. subsequently reported a male GIS patient with a normally sized gland who was heterozygous for p.Q202TfsTer99 in DUOX2 and p.C41S in TSHR [24].
Eleven patients carried oligogenic variants, but we found no significant correlations between the number of mutated genes and disease severity or morphological phenotype of the thyroid, which diverged from the conclusion drawn by Yamaguchi et al. that patients with oligogenic variants tended to be more severely affected than those with monogenic variants [25]. Similarly, one previous study has also provided relevant evidence to support this conclusion, where 2 PCH patients respectively had p.E879K in DUOX2 combined with p.P1012L in TG, and p.K530X in DUOX2 combined with p.R109Q in TSHR [16]. Nevertheless, Oliver-Petit et al. reported 8 PCH cases with GIS and goiter, who carried digenic variants in either DUOX2 and TPO, TG and DUOX2 or TG and TPO, and the CH severity differed among them, ranging from mild to severe [26].
The prevalence of TG variants (8.00%, 30/375) in our cohort was higher than that in the other seven genes, which was consistent with one previous finding that TG was the most commonly mutated gene in 49 GIS patients [27]. Surprisingly, we found that 25 out of the 30 patients carried single heterozygous TG variants, which was contrary to the generally believed autosomal recessive inheritance pattern and seemed unlikely to induce CH alone. However, there are also studies consistent with our findings that only single heterozygous TG variants were identified in certain CH patients [22,26]. Whether single heterozygous TG variants could cause CH is unclear and remains poorly explored. Double TG variants and variants in combination with other gene variants were detected in 5 individuals. The prevalence of TSHR variants (5.87%, 22/375) was second only to TG, which was in line with the mutation rate of TSHR (5.91%, 13/220) in a Chinese cohort of 220 patients [28]. The identified hotspot variant p.R450H in previous studies was detected in 8 cases, accounting for 36.36% (8/22) of the total patients with TSHR variants. The mutation frequency of GLIS3 was 3.73% (14/375), lower than previously reported by Zhang et al. (19/353, 5.4%) [29]. In our study, the most common variant in GLIS3 was p.H824R, which was found in 3 subjects. Among the 14 patients, 6 had digenic variants in GLIS3 and other genes (TG, TSHR, JAG1) and 8 carried monoallelic GLIS3 genetic alterations. However, heterozygous variants may not be responsible for CH. TUBB1, BOREALIN, NTN1 and JAG1 are newly identified causal genes for TD [30][31][32]. In this study, we performed large-scale CH population-based mutation screening for these genes by targeted NGS for the first time.
The correlations between phenotype and genotype have not yet been fully established. In the present study, we found that variants in GIS genes (TG and TPO) could result in TD, Fig. 4 Distribution of our CH cohort based on different number of mutated genes in turn, patients with variants in genes associated with thyroid development (TSHR, GLIS3, BOREALIN, NTN1, JAG1 and TUBB1) could develop CH with GIS, which suggested considerable phenotypic variability of CH. Most of the patients carrying monoallelic TSHR variants presented with GIS, but only one with variant p.F525S had athyreosis. Despite the fact that monoallelic and biallelic variants in BOREALIN could give rise to TD, our results showed that defect in BOREALIN may cause GIS [33]. Two GIS patients carried heterozygous p.E107D and p.L452P in TSHR, respectively [34], and two siblings with a homozygous TG variant p.W637X had a small thyroid volume [35]. In 2018, Stoupa et al. reported the first hypoplasia case with a homozygous TPO variant p.S292F in Tunisia [10]. In addition, patients with the same variant may show diverse phenotypes. For example, for two patients with monoallelic variant p.R846W in TPO, they manifested athyreosis and GIS with a enlarged thyroid gland, respectively. In 2016, Fu and colleagues identified a heterozygous TSHR variant p.V689G in two unrelated cases with decreased or normal-sized glands [36].
Although this is the largest cohort for mutation screening on 8 CH-related genes to the best of our knowledge, the following limitations should be considered: (a) Functional studies of novel variants are needed to fully elucidate their contributions in the occurrence of CH. (b) Pedigree validation is required to clarify the origin of the variants.
In summary, our findings implied that the underlying oligogenic mechanism may contribute to CH with the involvement of other factors, such as penetrance, epigenetic mechanisms and environmental factors. Further studies are needed to elucidate the complex genotype-phenotype correlations.

Data availability
The data analyzed during this study are available from the corresponding authors upon reasonable request.