Association of human height-related genetic variants with familial short stature in Han Chinese in Taiwan

Human height can be described as a classical and inherited trait model. Genome-wide association studies (GWAS) have revealed susceptible loci and provided insights into the polygenic nature of human height. Familial short stature (FSS) represents a suitable trait for investigating short stature genetics because disease associations with short stature have been ruled out in this case. In addition, FSS is caused only by genetically inherited factors. In this study, we explored the correlations of FSS risk with the genetic loci associated with human height in previous GWAS, alone and cumulatively. We systematically evaluated 34 known human height single nucleotide polymorphisms (SNPs) in relation to FSS in the additive model (p < 0.00005). A cumulative effect was observed: the odds ratios gradually increased with increasing genetic risk score quartiles (p < 0.001; Cochran-Armitage trend test). Six affected genes—ZBTB38, ZNF638, LCORL, CABLES1, CDK10, and TSEN15—are located in the nucleus and have been implicated in embryonic, organismal, and tissue development. In conclusion, our study suggests that 13 human height GWAS-identified SNPs are associated with FSS risk both alone and cumulatively.

signaling, the extracellular matrix, and cancer-associated pathways [14][15][16][17] . The results of genome-wide association studies (GWAS) have thus contributed new information on the susceptibility loci related to human height and provide insights into the extreme polygenic nature of human height.
Short stature is defined as a body height below the 3 rd percentile for the corresponding age and gender of the population. Among people affected by short stature, 37% have familial short stature (FSS), 27% have a constitutional growth delay, 17% have both FSS and a constitutional growth delay, 9% have a systemic disease, 5% have an endocrine condition, and the cause is unknown in the remaining 5% of cases [26][27][28] . Representing the majority of cases of short stature, FSS is characterized by a height below the 3 rd percentile in the population, a normal annual growth rate, normal bone age, family history of short stature, normal onset of puberty, and normal results of biochemical analyses. FSS is also called "genetic short stature" because the growth conditions are normal but the body height is below the 3 rd percentile and there is a family history of short stature. FSS represents a very suitable study object for research on short stature genetics because disease association has been ruled out for the short stature in FSS. FSS is caused only by genetically inherited factors. Thus, we initiated a pilot study to search for genetic loci that may predispose individuals to FSS. We explored the associations of FSS risk in Taiwan with genetic loci that have been related to human height according to previous GWAS, both alone and cumulatively.

GWAS-determined genetic variants related to human height and FSS risk in the Han Chinese population of Taiwan.
To identify the susceptible single nucleotide polymorphisms (SNPs) associated with FSS risk, we chose genetic SNPs that have been related to human height according to published GWAS, which were considered to be the closest genes potentially associated with FSS. There were 1,033 genetic SNPs reported in 14 human height GWAS conducted to date (Table S1). Of these, 34 SNPs were ultimately identified as candidate susceptibility SNPs for FSS according to the selection criteria [call rates of genotyping for both FSS and controls >95%, genotypes for genetic variants in controls conforming to Hardy-Weinberg equilibrium (p > 0.05), and p-value for the additive model <5.00E-05 (0.05/1,033)]. The 34 SNPs in 15 closest genes identified in human height GWAS are summarized in Table 1 and Supplementary Table S2. As shown in Table 1, the minor allele frequencies (MAFs) of these 34 genetic SNPs were similar to the MAFs in Han Chinese from Beijing (CHB) NCBI GRCh38.p2 assembly) from the online database of the dbSNP website (HAPMAP-CHB; http://www.ncbi. nlm.nih.gov/projects/SNP/index.html). The pair-wise linkage disequilibrium (LD) between the 34 SNPs was also estimated (Supplementary Tables S2-S5 and Figure S1).
The excluded SNPs were those with strong LD (D′ > 0.8). According to the LD results, we then re-selected these SNPs for calculation of the cumulative effect. The re-selected SNPs included 13 SNPs in the 13 closest genes: rs1046934 in TSEN15, rs3791679 in EFEMP1, rs3771381 in ZNF638, rs10935120 in CEP63, rs7632381 in ZBTB38, rs13131350 in LCORL, rs6845999 in HHIP, rs4240326 in ANAPC10, rs6470764 in GSDMC, rs12338076 in QSOX2, rs4842838 in ADAMTSL3, rs258324 in CDK10, and rs4308051 in CABLES1. The association was statistically evaluated by the odds ratio (OR) and its corresponding 95% confidence interval (CI) in the additive genetic model. As shown in Table 2, these 13 genetic SNPs showed a statistically significant association with FSS risk (p-value < 0.00005).
Cumulative effects of the 13 human height genetic SNPs on FSS risk. In accordance with the potential additive inheritance model for the individual SNP association analysis (Table 2), we subdivided the three genotypes of each SNP into three groups: risk allele homozygote (the risk genotype is coded as "2"), risk allele heterozygote (the risk genotype is coded as "1"), and Non-risk allele homozygote (the non-risk genotype is coded as "0") (Table S6, Figures S2 and S3). We calculated the multi-locus genetic risk score (GRS) for each individual by summing the number of risk alleles (0/1/2) for each of the 13 SNPs weighted by their estimated effect sizes (log OR) (Table S6 and Figure S3). Estimates of the association between the weighted GRS divided into quartiles (Q1-Q4) and FSS were calculated using a logistic regression model; individuals in genetic risk score Q1 served as the reference. As shown in Table 3, there was a continuous increase in FSS risk with the cumulative weighted GRS quartile numbers (p < 0.001; Cochran-Armitage trend test). That is, compared with individuals in Q1, those identified in Q2, Q3, and Q4 had a significant association with increased FSS risk, with an increasing risk for each quartile (Table 3). These results suggested a cumulative effect of these 13 SNPs on FSS risk. Functional analysis. We applied these 13 SNPs re-selected from the LD results to possible pathway mapping and further functional validation to explore possible functional relations among the genes affected by these 13 SNPs. The network analysis identified a single cluster of 35 genes that included 11 associated genes discovered in this study ( Figure S4) related to functions of embryonic development, organismal development, and tissue development. Six of these genes are all located in the nucleus: ZBTB38, ZNF638, LCORL, CABLES1, CDK10, and TSEN15.

Discussion
Human height GWAS have led to the discovery of thousands of novel genetic variants [14][15][16][17][18][19][20][21][22][23][24] . In this work, we initiated our FSS study on a single ethnic group (Han Chinese) and investigated the associations of FSS risk with genetic variants previously related to human height in GWAS, both alone and cumulatively. We identified a continuous increase in FSS risk with increasing genetic risk score quartiles, suggesting a cumulative effect of these 13 SNPs on FSS risk. These data suggest that FSS represents a suitable research object for investigating short stature genetics. To our knowledge, this is the first study reporting the genetic profile of FSS.
The LD analysis revealed several SNPs in strong LD that have been linked to adult height in previous genetic studies, mapping to genes such as COLGALT2 and TSEN15, EFEMP1 (implicated in bone and cartilage development pathways, including the constitutions of the extracellular matrix), ZBTB38, LCORL, HHIP, ANAPC10, GSDMC, QSOX2, ADAMTSL3, and CABLES1 14,[16][17][18][19][20]25 . Other SNPs were found to be independent, located near the ZNF638 and CEP63 genes. CEP63 recruits the cell cycle regulatory protein CDK1 to the centrosome, and has been associated with regulation of mitotic entry, centrosome amplification, and genome maintenance, thereby affecting cell proliferation. Only one SNP was identified for chromosome 16, related to CDK10.
Network cluster analysis revealed functions of these 13 genes that are mainly related to embryonic development, organismal development, and tissue development. Two of the genes found to be associated with FSS belong to the zinc finger protein family: ZBTB38 and ZNF638. There were seven relevant SNPs in ZBTB38: rs6440003, rs6763931, rs724016, rs7632381, rs1344672, rs9825379, and rs10513137 that have previously been reported in human height GWAS 12,16,18,19,25 . Moreover, genetic variants in ZBTB38, including rs724016, rs1582874, rs11919556, rs6440006, rs7612543, rs62282002, and rs18651435, have also been implicated in idiopathic short stature (ISS) 29 . ZBTB38 encodes a zinc finger-and BTB domain-containing protein 38. It belongs to the family of zinc finger proteins and thus serves as a zinc finger transcriptional activator that binds methylated DNA; such potential suppression of the transcription of methylated regions might affect adult stature by regulating the production of insulin-like growth factor-II 30 . Therefore, ZBTB38 may affect adult height, ISS, and FSS by regulating the protein expression of IGF-II. In addition, ZBTB38 also regulates the transcription of tyrosine hydroxylase, the rate-limiting enzyme in catecholamine synthesis 31 .
Moreover, the SNP rs3771381 in ZNF638, encoding zinc finger protein 68, was also associated with FSS risk. ZNF638 genetic variants have been reported in human height GWAS 25 that binds to cytidine-rich sequences in double-stranded DNA 33 . This zinc finger protein is a member of the large class of transcription factors that participates in skeletal development 34 and regulates adipocyte differentiation 35 . ZNF638 has been reported to be a transcriptional coregulator of the early regulator of adipogenesis via induction of PPARγ in cooperation with CCAAT/enhancer binding proteins (C/EBPs) 35,36 . LCORL encodes a ligand-dependent nuclear receptor corepressor-like protein and is a transcription factor. Polymorphisms in this gene are associated with measures of skeletal frame size and adult height in human and in livestock animals 24,25,[37][38][39][40] . SNPs affecting the transcription or translation of LCORL may result in up-or down-regulation of the expression of genes involved in growth 38 , and may affect bone formation, skeletal frame size, as well as body height; thus, further functional investigations of these variants are warranted. One SNP (rs13131350) within this gene was associated with FSS in our study. This SNP was also reported in a meta-analysis study of GWAS of adult height in East Asians 25 . Regulation of LCORL protein expression may affect bone formation, skeletal frame size, as well as body height; thus, further functional investigations are warranted.
There were two cell cycle regulation-related genes-CABLES1 and CDK10-associated with FSS in our study. According to the present study, four SNPs in CABLES1 are associated with FSS in Han Chinese in Taiwan: rs4800452, rs4369779, rs4308051, and rs8094261, which have also been reported in human height GWAS 14,18,23,25 . CABLES1 encodes a Cdk5 and Abl enzyme substrate protein 1, and serves as a candidate tumor suppressor that negatively regulates cell growth by inhibiting cyclin-dependent kinases 41,42 . Loss of CABLES1 protein expression has been observed at high frequency in human colorectal, lung, ovarian, and endometrial hyperplasia and in endometrial cancers [42][43][44][45] . CABLES1 protein overexpression induces apoptosis and inhibits cell growth by stabilizing p21 and decreasing cyclin-dependent kinase 2 kinase activity, thereby affecting cell proliferation 46 .
The other cell cycle-regulatory gene is CDK10. We found one SNP (rs258324) located in CDK10 associated with FSS, which was also reported in a human height GWAS 25 . CDK10 encodes cyclin-dependent kinase 10 and belongs to the CDK subfamily of the family of Ser/Thr protein kinases, which are necessary for cell cycle rs ID Gene Chr.
progression. CDK10 plays important roles in cellular proliferation and in regulation of the G2-M phase of the cell cycle 47,48 . Therefore, CABLES1 and CDK10 may affect adult height and FSS by regulating cell survival, cell cycle, and cellular proliferation. We identified one SNP (rs1046934) in TSEN15 associated with FSS. This SNP has been reported in a human height GWAS study 14 . TSEN15 encodes a highly conserved subunit of a tRNA-splicing endonuclease that catalyzes the removal of introns from tRNA precursors 49,50 . tRNA splicing is a fundamental process for cell growth and division and is highly conserved among vertebrates. Mutations in TSEN15 cause neurogenetic disorders, including pontocerebellar hypoplasia and progressive microcephaly [50][51][52] . These findings are suggestive of TSEN15′s importance in brain development. However, the roles of TSEN15 in adult height and FSS remain unclear and further functional characterization is needed.
In conclusion, our study suggests that 13 SNPs related to human height (according to previous GWAS) are associated with FSS risk alone and cumulatively. This is the first report on the genetic profile of FSS. Such information may be helpful for further research on short stature in the Han Chinese population, and for future functional studies dealing with growth and development. Confirmatory studies with a larger sample size and with functional characterization of these candidate genes would be warranted in the near future.

Methods
Study participants. The case population consisted of 978 individuals with FSS sequentially enrolled from the Children's Hospital, China Medical University, Taichung, Taiwan, from August 1999 to September 2014. All participants were unrelated Han Chinese children who were residents of central Taiwan and fulfilled the diagnostic criteria of FSS 26,28 . The criteria for FSS in our association study were (1) height less than the 3 rd percentile of the corresponding population, (2) normal annual growth rate, (3) bone age appropriate for chronologic age, (4) a family history of short stature, (5) normal onset of puberty, and (6) normal results of physical examination. The controls consisted of 1,129 subjects who were randomly selected from the Taiwan Biobank (http://www.twbiobank.org.tw/new_web/index.php). The criteria for controls for inclusion in our association study were (1) no history of FSS diagnosis, (2) body height above that of the top 75th percentile of the general population in Taiwan; and (3) age <61 years. All the participating FSS cases and controls were of Han Chinese origin, which constitutes 98% of the Taiwanese population. The ethnic background was assigned in accordance with the results of self-reported questionnaires. All research was performed in accordance with the relevant guidelines and regulations. The study protocol was approved by the institutional review board and the ethics committee of Human Studies Committee of China Medical University Hospital. Written informed consent was obtained from the participants, their parent, or legal guardian in accordance with the institutional requirements and the Declaration of Helsinki principles.
SNP selection, genotyping, and quality control. We initially included genetic SNPs listed as supplementary files from previous human height GWAS [14][15][16][17][18][19][20][21][22][23][24] . The genetic resources for the human height GWAS are shown in Table S7. The repeated SNPs were removed and the resulting 1,033 SNPs from 14 GWAS of human height are listed in Table S1. Genotyping of these 1,033 SNPs was performed in both cases and controls. Genomic DNA was extracted from peripheral blood leukocytes according to standard protocols (Qiagen Genomic DNA Isolation Kit; Qiagen, Taichung, Taiwan). The DNA concentration and quality were confirmed on a NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific, Taichung, Taiwan). SNPs of each sample were genotyped using a custom-designed VeraCode GoldenGate Genotyping Assay System (Illumina; http://www.illumina. com/) [53][54][55] . Genotypic data were quality controlled, and SNPs were excluded from further analysis if (1) the MAF was less than 5% in the controls (the Han Chinese population); (2) the total call rate was less than 95% for both cases and controls; or (3) the SNP significantly departed from Hardy-Weinberg equilibrium proportions for the controls (p < 0.05). Then SNPs that passed the correction test were considered to be significant (P = 0.05/1,033 = 0.00005).

Statistical analysis.
Hardy-Weinberg equilibrium for genotype distributions of controls in this study was evaluated by the goodness-of-fit χ 2 test (Table 1).
Alleles were identified by genotyping and direct counting. The difference of allelic frequency in the additive model between the cases and controls were measured by ORs with 95% CIs using logistical regression models ( Table 2). Deviation from the additive model was also tested for 34 selected SNPs and a p-value for the DOMDEV test <1.47 × 10 −3 was considered significant (p = 0.05/34) (Table S2). All statistical analyses were performed in the SPSS software, v12.0 for Windows (IBM, Armonk, NY, USA).
For the analysis of haplotype blocks, we used the Lewontin D′ and R 2 measure to evaluate the intermarker coefficient of linkage disequilibrium in both controls and FSS cases 56 . The confidence interval for LD was estimated using a resampling procedure and was used to construct the haplotype blocks 57,58 .