Twelve years of GWAS discoveries for osteoporosis and related traits:advances,challenges and applications

Xiaowei Zhu,Weiyang Bai and Houfeng Zheng ,4

INTRODUCTION

Osteoporosis is a common skeletal disease affecting~200 million people around the world;it is characterized by decreased bone density,bone microstructural damage and a consequent increase in bone fragility.1–2Nearly 22 million women and 5.5 million men were estimated to have osteoporosis in Europe3and 10 million in the United States,and this number continues to rise.4In China,~83.9 million people are estimated to suffer from osteoporosis,and this number,including osteopenia,should increase to~212 million people by 2050.5Bone fragility is a poor outcome of osteoporosis,where long-term therapy and medical management are needed,especially in elderly individuals.1By 2050,it is estimated that~51.1% of worldwide hip fracture cases will be from Asia.6Accordingly,the burden of treatment for osteoporosis and osteoporotic fractures has been rising very rapidly,with an annual cost of$17 billion to treat fractures in the United States.4,7In China,~2.33 million osteoporotic fractures occurred in 2010,costing$9.45 billion,and the annual costs are estimated to double by 2035.8Therefore,as aging-related diseases,osteoporosis and osteoporotic fracture inflict a substantial economic,social,and clinical burden.

Osteoporosis,as a complex disease,is influenced by many factors,including diet(calcium and protein intake),physical activity,endocrine status,coexisting diseases,and genetic factors.1Osteoporosis is mainly characterized by low bone mineral density(BMD),which is highly heritable,with heritability ranging from 50%to 80%.2,9To date,genome-wide association studies(GWASs)(Supplemental Note Box 1)and meta-analyses have discovered many loci that are associated with BMD,osteoporosis,and osteoporotic fractures.10–12Furthermore,next-generation sequencing(NGS)of large-scale samples has also uncovered novel rare/low-frequency variants in susceptible genes/loci for BMD,osteoporosis and fracture.13–14Recently,the approach of Mendelian randomization was widely used to identify the causative risk factors for osteoporosis by using GWAS results.15

In this article,we first reviewed the fruitful discovery achieved by GWASs and meta-analyses for osteoporosis and related traits in the last 12 years(Fig.1).We introduced several newly raised conceptual models,such as omnigenic models and natural selection,which might explain the mystery of missing heritability of complex traits.We then summarized the clinical use of GWAS findings in the bone field,such as the identification of causal clinical risk factors,the development of drug targets,and disease prediction.

Fig.1 Timeline highlighting important milestones during the 12 years of GWAS discoveries for osteoporosis and related traits.Blue boxes indicate the studies from the GEFOS and GENOMOS consortia.The green box indicates the studies focused on the Chinese population.Red boxes indicate GWASs including rare variants.Yellow boxes indicate the UK Biobank-based GWAS.BMD bone mineral density,ESR1 estrogen receptor 1,GWAS genome-wide association study,LRP5 low-density lipoprotein receptor-related protein 5,LRP40 low-density lipoprotein receptor-related protein 4,OPG osteoprotegerin,RANK receptor activator of nuclear factor-kappa β,RANKL RANK ligand,SPTBN1 spectrin beta,nonerythrocytic 1,WES,whole-exome sequencing,WGS whole-genome sequencing,ZBTB40 zinc finger and BTB domain containing 40

GWAS IN THE BONE FIELD

Measurement of bone mass

Most studies have focused on areal BMD(aBMD)obtained from a 2-dimensional projection scan with dual energy X-ray absorptiometry(DXA).14,16–19The T-score is measured in standard deviation(SD),a mathematical term that calculates how much one’s bone mass varies from the average.It defines an individual’s bone mass as normal(above?1 SD),osteopenia(between?1.0 and?2.5 SD)and osteoporosis(below?2.5 SD).20This measurement could be influenced by several different skeletal parameters,such as periosteal expansion,trabecular volumetric BMD(vBMD),cortical BMD,cortical thickness,trabecular number and trabecular thickness.21

Bone mass can also be assessed with other radiological imaging tools,such as quantitative computed tomography(QCT),which has the advantage of revealing unique bone information.Paternoster et al.22performed the first GWAS on cortical vBMD measured by QCT and found that the genetic variant rs1021188 near theRANKLgene was associated with the density of cortical bone,and rs9287237 onFMN2was associated with the trabecular bone fraction,while the other three SNPs(rs271170 onLINC00326nearEYA4,rs7839059 onCOLEC10nearOPGand rs6909279 onCCDC170nearESR1)had previously been reported to be associated with aBMD.However,QCT was not applicable to the WHO definition of osteoporosis that was based on DXA measurement,and QCT was more expensive with a higher dosage of exposure to radiation but might not predict fractures better than DXA measurement.23

An alternative method of estimating bone mass is derived from quantitative ultrasound(QUS).This measurement is quick,safe,and relatively inexpensive and can therefore be assessed in very large sample sizes,such as~500 thousand samples in the UK Biobank.The advantages over DXA make QUS a complementary(not replacement)approach to bone health assessment.QUS consists of the use of two separate ultrasound measurements,speed of sound(SOS)and broadband ultrasound attenuation(BUA),typically at the heel.Measures of estimated BMD derived from ultrasound were moderately correlated with DXA-derived BMD at the hip and spine.24A meta-analysis of GWASs25using heel ultrasound parameters identified a novel locus(rs597319 nearTMEM135)and replicated 6 previously reported loci(ESR1,SPTBN1,RSPO3,WNT16,DKK1andGPATCH1).

Early GWAS design

It has been established that the variation in BMD is the most important predictor for osteoporosis and fracture.Therefore,GWASs for osteoporosis mainly investigated the effect of genetic influence on BMD.In 2007,Kiel et al.26published the first GWAS,including 1 141 Framingham Heart Study subjects,and they identified 40 SNPs that could potentially be associated with several bone phenotypes(Fig.1).Unfortunately,owing to the small sample size,none of thePvalues exceeded the threshold of genome-wide significance(P<5×10?8).In 2008,two GWASs were published and identified 4 loci associated with BMD(LRP5,

ESR1,OPG,andRANKL);in addition,LRP5,RANK,LRP40,ZBTB40,andSPTBN1were found to be associated with the risk of osteoporotic fracture(Fig.1).16,27Soon afterwards,a deluge of GWASs were conducted on osteoporosis and related traits(Fig.2,Table 1 and Supplemental Table 1).

Fig.2 Genetic loci reported by the GWAS catalog for osteoporosis and related traits.#Fracture occurring at any site,except fingers,toes and skull,after age 18.BMC bone mineral content,BMD bone mineral density,BUA broadband ultrasound attenuation,FN femoral neck,LL lower limbs,OF osteoporotic fractures,SOS speed of sound

GWASs in the East Asian population

It is worth noting that the success of GWASs mainly came from studies performed in Caucasian populations,while only a few GWASs focused on East Asian populations(Fig.1 and Table 1).Yang et al.28performed a case-control GWAS in 700 elderly Chinese Han subjects(350 hip fracture patients and 350 healthy matched controls)and found thatUGT2B17copy number variation was associated with hip fracture.The same GWAS dataset was reanalyzed by Guo et al.29who found that the rs13182402 SNP inALDH7A1at 5q31 was strongly associated with hip fracture(Fig.1).29Kung et al.18conducted a GWAS and meta-analysis of BMD and fragility fractures in Chinese women(Hong Kong population)and found that the intronic SNP rs2273061 in theJAG1gene was strongly associated with the BMD of lumbar vertebrae(Fig.1).The first GWAS of osteoporosis conducted in a Japanese population30found that a common variant(rs7605378 onFONG)at 2q33.1 conferred the risk of osteoporosis in elderly individuals using a total of~6 700 subjects.30Recently,a GWAS with only 254 Japanese patients with inflammatory bowel disease(IBD)found that no SNPs reached genome-wide significance(P<5×10?8)for femoral neck(FN)and lumbar spine(LS)BMD.31Hwang et al.32performed an association study in 1 119 fracture patients and 3 444 controls in Korean and Japanese populations and found a newMECOMlocus associated with osteoporotic fracture[P=3.59×10?8,odds ratio(OR)=1.39]32.Another multistage GWAS meta-analysis identified a novel heel SOS locus(rs2446422 onGLDN)in the Korean population,33and the allele-specific epigenetic modifications of the SNP were confirmed using ENCODE annotations.33

GWAS meta-analysis

The allelic architectures of BMD and osteoporosis are likely to be multifactorial,with each factor imparting a relatively small effect.The identification of these loci with weak effects required studies with comprehensive coverage of the genome and very large sample sizes.The Genetic Factors of Osteoporosis(GEFOS)consortium(www.gefos.org)and the Genetic Markers for Osteoporosis(GENOMOS)consortium(www.genomos.eu)were employed to maximize the samples available for large GWAS meta-analyses,with a consequent increase in statistical power and new locus discovery.In 2009,a meta-analysis(GEFOS-1)of femoral neck and lumbar spine BMD was performed within the GEFOS consortium with 19 195 subjects of northern European descent and identified 20 BMD loci,including 13 novel regions that reached genome-wide significance(P<5×10?8)17(Fig.1 and Table 1).In 2012,the GEFOS consortium released their second-round GWAS meta-analysis results(GEFOS-2);19compared to the first GEFOS meta-analysis,the sample size of the study population increased significantly,which led to the identification of 56 loci associated with BMD(Fig.1 and Table 1).However,not all genome-wide significant results of the first GEFOS meta-analysis could be replicated in the second study because increasing the sample size could also lead to sample heterogeneity.34GEFOS-2 also revealed 6 loci associated with the risk of all types of fractures(FAM210A,SLC25A13,LRP5,MEPE,SPTBN1,and DKK1);19however,the definition of fracture in this study was quite heterogeneous,including hip,spine,wrist and other types of fractures.Therefore,these findings should be interpreted with caution before independent validation in other samples with homogeneous fracture types.34Zhang et al.35conducted a three-stage GWAS meta-analysis,and two novel loci were identified in the pooled sample of males and females(SMOC1)and in the female-specific sample(CLDN14);they also independently confirmed 13 previously reported loci(ZBTB40,

Table 1.continued

Table 1.continued

Table 1.Genome-wide association studies conducted on osteoporosis and related traits

GPR177,FGFRL1,MEPE,MEF2C,ESR1,SHFM1,WNT16,OPG,SOX6,LRP5,AKAP11,andFOXL1).Further gene expression analysis in osteogenic cells implied a potential functional association of theSMOC1andCLDN14genes in bone metabolism35.For fracture risk,the largest GWAS meta-analysis to date,including 25 cohorts from Europe,the United States,East Asia and Australia,identified 15 genetic loci for fracture,all of which also influenced BMD.36

GWASs including rare variants

Early GWAS design included only common variants[minor allele frequency(MAF)>=5%]and poorly covered low-frequency(1%<=MAF<5%)and rare variants(MAF<1%).Next-generation sequencing technology provides an approach to capture rare and low-frequency variants,which might be identified to be associated with complex traits with large effects.Styrkarsdottir et al.13performed the first whole-genome sequencing(WGS)study for BMD in an Icelandic population and found a rare nonsense mutation(c.376 C>T)withinLGR4that was strongly associated with low BMD and osteoporotic fracture.However,the mutation was not present in the public Exome Variant Server(EVS)database or in the Australian samples(Fig.1)13.In 2015,Zheng et al.14integrated WGS data(n=2 882),whole-exome sequencing data(n=3 549),deep imputation of genotyped data(n=26 534)37,and de novo replication genotyping data(n=20 271)and found that a low-frequency noncoding variant(rs11692564_T,nearEN1,MAF=1.6%)could result in an increased lumbar spine BMD(effect size=+0.20 standard deviation,P=2×10?14)and reduced fracture risk(OR=0.85,P=1×10?11)(Fig.1).Conditional loss ofEn1in a cre/flox mouse model resulted in osteopenia and increased skull bone resorption via an indirect effect sinceEn1was not expressed in osteoclasts.14

Large-scale biobank based GWASs

The UK Biobank(www.ukbiobank.ac.uk)recruited 502 647 individuals aged between 37 and 76 years from all over the country in 2006–2010,and the heel bone quality of the participants was evaluated by quantitative ultrasound SOS and BUA.In 2017,Kemp et al.38conducted a GWAS of 142 487 individuals from the UK Biobank using BMD as estimated by quantitative ultrasound of the heel.They demonstrated that 307 conditionally independent SNPs attained a genome-wide significance level at 203 loci,explaining~12% of phenotypic variance(Fig.1).Next,they investigated the underlying mechanism of these SNPs by four steps(including a.bioinformatic,functional genomic annotation and human osteoblast expression studies;b.gene function prediction;c.skeletal phenotyping of 120 knockout mice with deletions in genes adjacent to lead independent SNPs;d.the analysis of gene expression in mouse osteoblasts,osteocytes and osteoclasts)and suggested that theGPC6gene was a novel determinant of BMD and the pathophysiology of osteoporosis.38

A new study by Morris et al.evaluated genetic determinants of BMD as estimated by heel quantitative ultrasound in 426 824 individuals,identifying 518 loci(301 novel)that reached a genome-wide significance level,explaining 20% of its variance(Fig.1).39They also undertook a meta-analysis of~1.2 million individuals and identified 13 fracture loci(all associated with heel BMD),highlighting the importance of BMD as a determinant of fracture risk.39They found that target genes were enriched in those known to influence bone density and strength from cellspecific features(maximum OR=58,P=1×10?75)and found an increased abnormal skeletal phenotype frequency through the phenotyping of 126 knockout mice with disruptions in predicted target genes.39Finally,DAAM2showed critical effects on bone strength,porosity,composition and mineralization.

While most BMD GWASs analyzed data derived from DXA,these two studies of UK Biobank data used estimated BMDs derived from calcaneus ultrasound.Though Pearson’s correlation coefficients between DXA and QUS parameters showed a moderate association(r=0.42–0.61)24,40and quantitative ultrasound had the ability to predict the occurrence of fractures in older women41and men,42there were some essential differences.The GWASs did not replicate 18 known loci from previous studies utilizing DXA-derived BMD measures,and 6 loci had opposite effects on heel BMD in the study of Kemp et al.38compared to previous DXA BMD studies.Although these differences may be due to various reasons,differences in measurement by QUS and DXA were likely to be causes.

Age-and sex-specific BMD GWASs

Although most studies have focused on adults,GWASs have also performed in younger individuals,including children,43–45teenagers,22and premenopausal women.46–47The first GWAS reported for BMD in children identified theSP7locus,which encodes the transcription factor osterix,as being associated with whole-body BMD,and replication was subsequently achieved in adult lumbar spine BMD.44Recently,Chesi et al.43found that two loci achieved genome-wide significance:rs7797976 withinCPED1in girls and rs7035284 nearMTAPin boys at the distal radius.Actually,signals at theCPED1-WNT16-FAM3Clocus have been previously reported to be associated with BMD at other skeletal sites in adults17and children(skull and total body aBMD)of European ancestry45.Interestingly,this locus was also associated with cortical bone thickness,bone strength,and the risk of forearm fracture in adults;48peak bone mass in premenopausal women;47and BMD and fracture in elderly individuals.19The integration of functional studies inWnt16knockout mice revealed reductions in bone mineral content(BMC),bone area and bone strength.45,48Both natural variation in humans and functional studies inWnt16knockout mice demonstrated thatWNT16was an important determinant of the bone mass at different body sites in children and adults and the risk of fracture,suggesting that this genetic effect acted over the whole lifetime.

To detect genetic variants influencing variability in peak BMD in premenopausal women,Koller et al.conducted a GWAS in 1 524 US Caucasian women(aged 20–45 years)and 669 African American women(aged 20–44 years).A novel gene,CATSPERB,was identified to be significant in femoral neck BMD.46CATSPERBwas not found to be significant in the meta-analyses of GEFOS-2,although the samples from the above study were included in GEFOS-2.19Later,a meta-analysis was carried out restricting samples to premenopausal white women from 4 cohorts(n=4 061,aged 20–45 years),and two loci(WNT16)and(ESR1/C6orf97)were identified to influence the peak bone mass at the lumbar spine and femoral neck.47Only 4 out of the 56 GEFOS-2 loci19were observed to havePvalues below 5×10?5in this meta-analysis.

Although most of the published GWASs on skeletal phenotypes did not have adequate power to test sex-specific genetic effects,49there was suggestive evidence for an interaction between sex and SNP rs1021188(nearRANKL)(P=0.01),with a stronger association in males than females(at age 15,males?6.77 mg·cm3per C allele,P=2×10?6;females?2.79 mg·cm3per C allele,P=0.004).22In the GEFOS-2 study,19two loci(Xp22.31 in men and 8q13.3 in women)were discovered in the sex-stratified meta-analysis;however,only the locus in Xp22.31(nearFAM9B)showed significant heterogeneity(Phet=1.62×10?8),and the imbalance in sample size between women and men and the conservative heterogeneity test limited the ability to identify sex-specific findings.In a study of European American children(n=1 419),50four novel loci(IZUMO3,RBFOX1,SPBT,andTBPL2)were identified to be associated with BMD at the 1/3 distal radius,spine,total hip and femoral neck,two of which were sex-specific loci(SPTBin females andIZUMO3in males).

Bone size/geometry GWASs

Bone size(BS)is also an important factor that influences bone geometry and bone strength.To date,a limited number of GWASs for BS have been conducted compared to studies on BMD.In Table 1,we summarized the current GWASs and meta-analyses in bone size/geometry.In early studies,GWASs of the bone area of the hip or lumbar spine using DXA did not find significant loci,26,51–53possibly due to the small sample size.Recently,a study of a hip shape model(HSM)derived from statistical shape modeling of DXA scans found 8 loci associated with hip shape,54and another GWAS meta-analysis identified 22 significant loci(P<5.0×10?8)for hip bone size.55Styrkarsdottir et al.reported a large GWAS of bone size using a simple parameter from DXA scans,the bone area,56they found that 8 loci for the lumbar spine area,5 loci for the total hip area,4 loci for the intertrochanteric area,3 loci for the trochanter area,and 1 locus for the femoral neck area satisfied the criteria of genome-wide significance(Table 1).

GWAS FINDINGS CANNOT PERFECTLY EXPLAIN THE VARIANCE IN BONE MASS

GWASs on osteoporosis and related traits have made great achievements in the past 12 years and have highlighted many genes/loci and related biological pathways that contribute to the pathophysiology of osteoporosis and/or fracture,such as the RANK-RANKL-OPG and WNT signaling pathways.These pathways are functionally relevant to bone metabolism and endochondral ossification,and their contribution to osteoporosis has been well established.57However,at the same time,similar to other complex traits,the variance in bone mass could not be fully explained by GWAS findings.

Missing heritability and beyond

Unlike Mendelian diseases that are caused by mutations in coding regions,most of the associated SNPs for osteoporosis and related traits are found in noncoding intergenic and intronic regulatory regions.62Therefore,the greatest challenge was to understand the functional consequences of these SNPs and to accurately elucidate the biological mechanism by which these genes and SNPs act.To date,only a small fraction of SNPs/genes and their functional mechanisms have been successfully characterized,63and these variants or regions could be transcription factor binding sites that regulate or affect gene expression.62

Polygenicity and negative selection

GWASs of osteoporosis and related traits often identified a number of SNPs that had significant p values but showed very low disease odds ratios(ORs).For example,in the GEFOS-2 study,1913 of 14 SNPs associated with any low-trauma osteoporotic fracture had ORs<1.10.A recent GWAS involving 53 184 fractures and 373 611 controls39identified 14 association SNPs,all of which had ORs<1.10.In a GWAS of Chinese fractures,the highest OR of rs13182402 in theALDH7A1gene was 2.25.29Despite their statistical significance,the ORs were small and explained little about the genetic contribution to fracture.

Over the last few years,a commonly accepted explanation for the small OR was that osteoporosis was caused by a large number of interacting genes,each with a small effect size and additive increment to disease risk,called“polygenic inheritance”.It is known that common diseases have a polygenic genetic architecture.64Thus,perhaps in many cases,the so-called problem of missing heritability might be synonymous with high polygenicity(defined as the total number of genetic loci or alleles with nonzero effects contributing to a phenotype).61,65The classic polygenic model consists of contributions to disease risk from both common and rare variants.61In 2018,using UK Biobank data,Zeng et al.confirmed that negative selection played a predominant role in shaping the relationship between effect size and MAF for complex traits.66They found that 23 out of the 28 studied complex traits(including heel BMD)showed significant signatures of natural selection,and the genetic variants associated with heel BMD were under negative selection,with a moderate estimate(S?=?0.381),where S?reflected the strength of selection on the trait-associated SNPs.66

More recently,O’Connor et al.redefined polygenicity as the effective number of independently associated SNPs(Me).For the 33 complex traits they studied,the“Me”estimates for common SNPs ranged from 500 to 20 000,with a‘Me’estimate for heel BMD of~800.67This implied that most common SNPs were associated with complex traits and that heritability was spread evenly across the genome.67They found that functionally important regions in the genome had higher polygenicity and higher heritability,but low-frequency SNPs had lower polygenicity than common SNPs on average.The conclusion was that negative selection not only constrained the effect sizes of common variants on average but also flattened their distribution across the genome.67

Polygenicity and omnigenicity

Recently,Boyle et al.68proposed the“omnigenic model”in which gene regulatory networks were fully interconnected;that is,all genes expressed in disease-related cells were considered to affect disease phenotype,but most heritability could be explained by effects on genes outside core pathways.This model tried to answer 2 questions:(1).Why do the lead hits from GWASs for any given trait contribute so little to heritability?(2).Why does so much of the genome contribute to heritability?The key feature of this model was the classification of genes as“core”(direct roles in disease)or“peripheral”(essentially all other expressed genes can transregulate core genes).In fact,the“omnigenicity model”is one scenario of the“polygenicity model”,in which the“polygenicity”is partitioned into different parts.

In the latest point of view,they defined the“core gene”as the only gene from which the gene product(protein or RNA for a noncoding gene)had a direct effect—not mediated through the regulation of another gene—on cellular and organismal processes,leading to a change in the expected value of a particular phenotype;“peripheral genes”were defined as those expressed in relevant cell types that could affect the phenotype only indirectly through regulatory effects on core genes.69This model assumed that the relationship between each core gene and the expected phenotype value was a linear function of the gene expression level;moreover,each core gene was likely affected by large numbers of weak trans(peripheral)variants,and most trait heritability was mediated through trans effects.69

Based on this model,most variants that contributed to heritability tended to be spread across the whole genome,and genes with specific functions for osteoporosis or related traits could only explain little heritability.This might explain why some loci/genes identified by GWASs for BMD were considered to have no contribution to the pathophysiology of osteoporosis and/or fracture,while some genes that had known functional relevance to bone metabolism and endochondral ossification tended to be core genes.For example,LRP5,which encodes low-density lipoprotein receptor-related protein 5,could function as a coreceptor together with the seven-transmembrane-spanning Frizzled for Wnt proteins to regulate intracellular signal transduction by β-catenin,70–71and the activation of the Wnt pathway results in cytoplasmic β-catenin accumulation.Consequently,β-catenin translocates to the nucleus and in turn regulates osteoblast proliferation and differentiation,thus determining bone mass.72Osteoblasts produce RANKL following the binding of RANKL to RANK on the surface of osteoclastic precursors,and subsequently,NF-κB is activated and translocates into the nucleus and interacts with NFATc1 to trigger osteoclastogenic gene transcription.73OPG,a member of the tumor necrosis factor(TNF)receptor superfamily(TNFRS),also known as TNFRS member 11B(TNFRS11B),can bind toRANKLto prevent its coupling withRANKand inhibit the maturation of osteoclasts as a result of reducing bone resorption.Notably,these genes related to the bone metabolism pathway,such asLRP5,16RANKL,22ESR1,25,38–39BMP4(bone morphogenetic protein 4),38,74andWNT16,were identified in GWAS signals.38,45,48

首先，進一步發(fā)揮好哈薩克小說在思考人生、珍視大自然等方面獨特的人文教育功能，并在小說中突出體現人性境界提升、理想人格塑造以及個人與社會價值實現的人文教育理念。如艾克拜爾·米吉提的小說《車禍》，通過對區(qū)域人物與事件的敘述描寫，從細微的心理刻畫上，在美心美行的浸潤中探究人與人之間的美好人性。

CLINICAL RELEVANCE OF GWAS FINDINGS

The ultimate goal of genetic study is to translate the discoveries into clinical practice.GWAS discoveries for osteoporosis and related traits in the past 12 years are undoubtedly more fruitful than previous linkage analyses and candidate gene association analyses,and hundreds of loci(thousands of SNPs)have been identified that are significantly and robustly associated with osteoporosis and related traits(Fig.2,Table 1 and Supplemental Table 1).However,it is still too early to understand the function of novel proteins identified by GWASs.This review is not meant to describe novel discovered loci and their interactions.We assumed that there are three ways in which GWAS findings could provide important clinical insight for osteoporosis.First,GWAS results could be employed to investigate the causal risk factors for osteoporosis by using the Mendelian randomization approach.Second,new drug targets and anti-osteoporotic therapeutics should be investigated.Despite the small effect of common variants identified by GWASs,it should be noted that the effect size of the genetic variant on molecular phenotypes could be large,and the drug effect on targets could also be magnified(e.g.,statins).75Third,genetic information could be applied to“personalized”medicine,for example,disease prediction and risk stratification,leading to the overall improvement in disease prevention or intervention.

Mendelian randomization approach to link clinical risk factors to osteoporosis and fracture

The identification of causative risk factors is essential for the prevention and treatment of osteoporosis,and a better understanding of causality could be conducive to further prevention strategies and clinical trials and to providing targets for effective lifestyle and drug intervention.15,76Observational studies have identified associations of potential risk factors(for example,smoking,low body mass index(BMI),low vitamin D level,earlier age at menopause and physical inactivity)with fracture risk.However,because of confounding factors and reverse causality,bias might be introduced into observational studies,thereby reducing their reliability.The gold standard for evidence for causal effects could come from well-conducted randomized control trials(RCTs),but RCTs are resource-intensive and examine mainly shortterm exposures.In addition,not all risk factors can be investigated by RCTs.15

Recently,Mendelian randomization(MR)analyses have been widely used to illustrate the causal effect between exposures and outcomes using large-scale GWAS summary statistic data.77–79MR is a type of analytical approach that takes genetic variants associated with a risk factor(e.g.,calcium)as instrumental variables(IVs)to examine the causality between exposure and outcome(e.g.,BMD).80–81Since the genetic alleles are randomly assorted during conception,MR analyses are less susceptible to confounding factors;additionally,MR analyses are robust to reverse causation bias because genotypes are unlikely to be affected by disease.Further information can be found in Supplemental Note Box 2.Three main assumptions must be applied when conducting Mendelian randomization analyses.82First,the genetic variants should be strongly associated with exposure(the relevance assumption);second,the genetic variants should be independent of factors that confounded the exposureoutcome relationship(the independence assumption);and third,the genetic variants affect the outcome only through the exposure(the exclusion restriction assumption)(Fig.3,Panel A).This approach has advantages over traditional observational studies by minimizing confounding bias.To date,the MR approach in the bone field has been applied predominantly to assess causal relationships between different factors and BMD,osteoporosis and fracture(Fig.3,Panel B and Table 2).Among these factors,vitamin D level,inflammatory disease,obesity and diabetes were frequently investigated.

Fig.3 Mendelian randomization in bone field.Panel A:Principal of Mendelian randomization.Panel B:The causality between the clinical risk factors and osteoporosis from the current literature.Red boxes indicate the causal relationship.Black boxes indicate the noncausal relationship.Blue boxes indicate controversial results.BMD bone mineral density,CAD coronary artery disease,IBD inflammatory bowel disease,DBP vitamin D binding protein,PsA psoriatic arthritis,T1D type 1 diabetes,T2D type 2 diabetes,TSH thyroid stimulating hormone

Vitamin D and BMD/fracture.Vitamin D,by improving intestinal calcium absorption,has pivotal roles in bone heath.Vitamin D insufficiency was reported as a risk factor for several common diseases and conditions,including osteoporosis and osteoporotic fracture.83However,the influence of vitamin D on the etiology of low bone mass and osteoporosis is unclear due to inconsistent results from clinical studies.79Leong et al.84investigated the causal relationship between vitamin D-binding protein(DBP)levels and BMD in the Canadian Multicentre Osteoporosis Study(CaMos)using individual-level data,and the results demonstrated a strong causal relationship between serum DBP and 25OHD levels;however,serum DBP had no causal effect on femoral neck BMD or osteoporosis(Table 2).Furthermore,Li et al.85found no evidence for a causal effect of vitamin D levels on BMD(total hip,FN and LS)in Chinese postmenopausal women using four SNPs,GC-rs2282679,NADSYN1-rs12785878,CYP2R1-rs10741657 andCYP24A1-rs6013897,as candidate instrumental variables in the MR analyses.Recently,using data from the GEFOS consortium and UK Biobank,Larsson et al.86found that vitamin D levels had no effect on BMD(FN,LS,heel)(N=32 965).Recently,a study also showed a lack of a causal relationship between vitamin D levels and fracture risk by using 37 857 cases and 227 116 controls from the GEFOS Consortium,UK Biobank,EPIC-Norfolk study and 23andMe36(Fig.3,Panel B and Table 2).Similarly,as a provider of protein,micronutrients and dairy calcium,milk was recommended by some dietary guidelines,particularly for bone health.However,MR studies using a SNP(rs4988235)located upstream of the lactase gene as an instrumental variable found that milk consumption had no causal effect on BMD87or fracture36(Fig.3,Panel B and Table 2).

Diseases and BMD/fracture.To date,diseases such as type 2 diabetes(T2D)and inflammatory diseases have been studied for their effects on osteoporosis or fracture(Table 2).Trajanoska et al.36found that IBD was not a causal factor for fracture risk in 185 057 cases and 377 191 controls.More recently,with 432 513 samples from the UK Biobank dataset,Xia et al.found that psoriatic arthritis might be a risk factor for low BMD,but the link was not genetically determined.Psoriasis without arthritis is not a risk factor for osteoporosis.88

Table 2.continued

Table 2.continued

Table 2.continued

Table 2.Mendelian randomization studies in the bone field

By using SNPs as IVs[32 SNPs strongly associated with type 2 diabetes(T2D),30 SNPs associated with fasting glucose and 4 SNPs associated with 2-h glucose(2hGlu)],Ahmad et al.89found that a genetically increased risk of T2D and a genetically increased risk of fasting glucose both had weak effects on increasing femoral neck BMD,but no significant trends were observed for the effect of T2D and fasting glucose on lumbar spine BMD.89Furthermore,Trajanoska et al.36found that T2D and fasting glucose were not causal for fracture in 185 057 cases and 377 191 controls.The study also reported no causal effect of type 1 diabetes(T1D)and coronary artery disease(CAD)on fracture.36

Other factors and BMD/fracture.Fat mass might be a causal decisive factor of bone mass,but the evidence was contradictory.90–93By using variants of two loci[FTO(fat mass and obesity-associated gene)and MC4R(melanocortin 4 receptor)]strongly associated with fat mass and obesity,Timpson et al.evaluated the relation between fat mass and bone outcomes in~5 000 children at a mean age of 9.9 years from the Avon Longitudinal Study of Parents and Children(ALSPAC)cohort and suggested that fat mass was the causal pathway for bone mass in children.94In 2016,a study investigated whether adiposity was causal for BMD at the skull,upper limbs and lower limbs,pelvis and lumbar spine in 5 221 children from ALSPAC using 32 SNPs(strongly associated with BMI),and the results suggested that adiposity was causally related to increased BMD at all sites except the skull.95The relationship between obesity and BMD was also investigated in adults,and it was found that obesity might be causally related to BMD at the femur but not at the spine.96In addition,the MR approach has been used to show a positive causal association between serum estradiol concentrations and femoral neck BMD,lumbar spine BMD and heel BMD.86Other studies demonstrated that earlier menopause and late puberty were causal factors for increasing fracture risk.36,97However,urate,78–98thyroid stimulating hormone(TSH),36,99homocysteine,36alcohol consumption100and smoking status36,100were not identified as causal factors for BMD or fracture by the MR approach(Fig.3,Panel B).Notably,it was demonstrated that genetically decreased BMD was the only clinical risk factor with evidence for an effect on fracture risk among 15 clinically identified fracture factors.36More recently,Cerani et al.undertook an MR study and found that a standard deviation increase in genetically derived serum calcium(0.13 mmol·L?1or 0.51 mg·dL?1)was not associated with increased estimated BMD(426 824 subjects,P=0.92)or a reduced risk of fractures(76 549 cases and 470 164 controls;P=0.85).101

Therapeutic targets for osteoporosis

Despite the small effect size of common variants identified by GWASs,most of the osteoporosis agents in use(or undergoing trials)target pathways related to the GWAS-discovered BMD genes,and genetic information might significantly improve the search for drug targets and increase the success rate of preclinical and clinical trials.102Moreover,it is well recognized that the effect size of association is not well correlated with clinical relevance,as many FDA-approved medications target proteins linked to common variants identified by GWASs.102–104An example of success in the field was the use of GWAS data for drug repositioning studies.Sanseau et al.105found that among the publicly relevant disease-related GWAS loci,155 out of the 991 loci(15.6%)were related to drug development.Among them,the drug indications of 63 targeted proteins matched the corresponding GWAS traits,indicating that the pathogenic genes excavated by GWASs had a higher probability of being directly used as drug targets.105For example,theIL12B(interleukin 12B)gene found in psoriasis GWASs encodes the target of ustekinumab,a newly proven drug for psoriasis.In addition,the gene was considered to be related to Crohn’s disease,and thedevelopment of related drugs was in a phase II clinical trial.105Another example was denosumab,which is a drug marketed for the treatment of osteoporosis in postmenopausal women,targeting the geneTNFSF11(tumor necrosis factor superfamily,member 11),also known asRANKL.Denosumab is a RANKL inhibitor that functions by preventing the development of osteoclasts.Recently,it was speculated that the drug might have a therapeutic effect on Crohn’s disease,asTNFSF11was found to be significantly associated with Crohn’s disease in GWASs.106The current drugs that are available for the treatment of osteoporosis and their most likely targets are listed in Table 3.Five antiosteoporosis therapeutics currently approved or in advanced clinical trials were supported by GWAS data.It was reasonable to believe that the findings of GWASs could be potentially powerful in the identification of anti-osteoporosis drug targets and drug repositioning.

Table 3.Present and potential near-term osteoporosis drug targets that have been linked to changes in BMD by GWAS.Table adapted from178

TheSOST(sclerostin)gene was found to be strongly associated with BMD by GWASs;107SOSTproduces sclerostin,which is a key Wnt pathway regulator that is preferentially expressed by osteocytes.Sclerostin acts by binding to the Wnt coreceptor LRP5/6 by competing with Wnt protein;as a consequence,sclerostin blocks the accumulation of β-catenin in the cytoplasm,inhibits the differentiation and proliferation of osteoblasts,enhances osteoclastogenesis and causes bone loss.108–109Given the inhibitory effect of sclerostin on osteoblast function and bone formation,blocking the activity of sclerostin to activate this pathway seems to be a potential strategy in the treatment of osteoporosis.Romosozumab(AMG785/CDP-7851),a monoclonal humanized antibody to sclerostin,was evaluated for its efficacy.Compared with the traditional bone resorption inhibitor alendronate and the bone formation promoter teriparatide,the greatest feature of romosozumab was its ability to reverse postmenopausal osteoporosis in women with hormone deficiency.110–112Saag et al.compared the effect between romosozumab(210 mg monthly administered subcutaneously)and alendronate(70 mg weekly)for 12 months,followed by open-label alendronate 70 mg weekly for another 12 months in postmenopausal women with osteoporosis and a fragility fracture113.After 24 months,a lower risk of fractures,including clinical fractures(27% lower),hip fractures(38%),new vertebral fractures(48% lower)and nonvertebral fractures(19%),was observed in the romosozumab-toalendronate group than in the alendronate-to-alendronate group.113A phase III clinical trial was conducted to estimate the effect of romosozumab(n=206)versus teriparatide(n=209)on osteoporosis in postmenopausal women who took oral bisphosphonate for at least 3 years,and it was found that romosozumab(210 mg once monthly)had a greater effect on hip BMD than subcutaneous teriparatide(20 μg once daily)114.Another trial recruited 7 180 postmenopausal women who had osteoporosis,and the subjects were randomly assigned to receive subcutaneous injections of romosozumab(at a dose of 210 mg)or placebo monthly for 12 months;thereafter,both groups received denosumab 60 mg every 6 months twice.115At the end of the initial 12 months,romosozumab had decreased the incidence of new vertebral fractures and nonvertebral fractures by~73% and 24%,respectively.115At 24 months,a 75% lower risk of vertebral fractures was seen in the romosozumab group after the transition to denosumab.115

On 9 April 2019,the US Food and Drug Administration(FDA)approved romosozumab for the treatment of osteoporosis in postmenopausal women at high risk of fracture,with a boxed warning highlighting the risk of cardiovascular adverse events and a postmarketing requirement to assess the cardiovascular safety of romosozumab.116On 28 June 2019,the European Medicines Agency(EMA)recommended the refusal of the marketing authorization for Evenity(romosozumab)because the results suggested that patients given Evenity had an increased risk of serious effects on the heart and circulatory system,such as heart attacks or strokes.117In addition,there were more deaths in patients aged over 75 years who were given the medicine.As it was unclear why the medicine appeared to increase the risk of heart and circulatory problems,measures to reduce the risk could not readily be put in place.117Dramatically,after re-examining initial opinions,the EMA noted that the medicine showed convincing evidence of benefit in women with severe osteoporosis,with better effect than alendronate,and it was suggested that only women who had no history of heart attack and stroke could take the medicine.117On 17 October 2019,the EMA recommended that marketing authorization be granted but for a restricted indication in postmenopausal women with severe osteoporosis at high risk of fracture.117

Dickkopfs(DKKs)are secreted proteins composed of two cysteine-rich domains with four homologous forms(DKK-1~4)in vivo.DKK-1 inhibited the Wnt/β-catenin signaling pathway by directly binding to LRP5/6 and formed a complex with Kringen,a transmembrane protein containing a Kringle domain,which increased endocytosis and decreased LRP5/6 content,thus leading to the inactivation of the Wnt pathway.118–119DDK1 is closely related to bone mass,120–121and similar to sclerostin monoclonal antibodies,monoclonal antibodies to DKK-1 increase trabecular mass and density in mice122and restore bone density in osteoporotic mice and rhesus monkeys.123Monoclonal antibodies to DKK-1 included BHQ880124and PF04840082,125but both were in the preclinical phase.

Prediction of osteoporosis and fracture

One of the goals of genetic study is to improve the value of clinical application,for example,to predict osteoporosis or fracture risk from GWAS findings.Studies have shown that at least 150 loci with an OR value of 1.5 or 250 loci with an OR value of 1.25 were required for the prediction of disease risk.126This suggested that any single locus could not be useful in clinical prediction,regardless of the size of the effect.However,theoretical and empirical studies have suggested that profiling multiple variants that are associated with bone phenotypes could improve the accuracy of fracture prediction and classification beyond that obtained by conventional clinical risk factors,such as the Fracture Risk Assessment Tool(FRAX).127

Polygenic risk scoring.Polygenic risk scoring was one primary approach for disease risk prediction.In a semisimulation study for fracture,it was shown that a profiling of up to 25 genes/variants(each with a relative risk of 1.10–1.35 and frequency ranging from 0.25 to 0.60)in the presence of clinical risk factors—with or without BMD—could achieve an AUC of 0.80.128Ho-Le et al.took 62 BMDassociated SNPs to define the predictive value of genetic profiling for fracture prediction in 557 men and 902 women and found that individuals with a greater polygenic risk score(PRS)had a lower femoral neck BMD(P<0.01);each unit increase in PRS was associated with a hazard ratio of 1.20 for fracture,and this association was independent of age,prior fractures and falls.129.However,polygenic risk scoring remained limited due to the linkage disequilibrium(LD)pruning of SNPs(prioritizing the most significant associations up to an empirically determinedPvalue threshold,and pruning the SNPs based on LD).130To remediate this issue,recent developments in machine learning may be a novel strategy.131–132Machine learning approaches adapted a set of sophisticated statistical and computational algorithms to make predictions by mathematically mapping the complex associations between a set of risk SNPs to complex disease phenotypes.133The optimal predictive ability for the target disease was obtained by mapping the pattern of selected features in the training genotype data,and at the end of the training stage,the model with the maximum predictive ability of the training dataset was selected for validation.131,134Machine learning has been applied to the prediction of diseases or traits,such as inflammatory bowel disease,135Alzheimer’s disease,136cancers,137–138heart failure139and height.140

Machine learning methods.Through the analysis of 341 449 individuals from the UK Biobank,Forgetta et al.tested whether machine learning methods could provide a clinically relevant genomic prediction of quantitative ultrasound speed of sound(SOS)—a risk factor for osteoporotic fracture.141In the Model Selection Set,age,sex and BMI explained 4.0% of the variance in SOS;the addition of the remaining FRAX clinical risk factors increased the variance explained to 4.8%,whereas when polygenic risk scores across differentPvalue thresholds were added,the variance explained increased to at most 18.5%.141Surprisingly,the machine learning algorithm improved the explained variance in SOS to a maximum of 25.0%.Then,they selected the top model(the machine learning algorithm selected 21 717 activated SNPs with aPvalue≤10?4)from the Model Selection Set to test for its correlation with the SOS in the validation set and found that the model could explain 23.2% of the variance in the measured SOS.Subsequently,they evaluated the associations among SOS,genomically predicted SOS(gSOS),BMD and fracture and found that decreased SOS and fracture were both strongly associated with increased odds of incident fracture(gSOS had the highest risk per SD)in the univariate model.However,in multivariate models,gSOS was more strongly associated with major osteoporotic fracture than SOS or BMD.141For fracture prediction,gSOS outperformed FRAX clinical risk factors alone.The machine learning algorithms provided better predictions than traditionally used polygenic risk scores.These findings suggested that genetic profiling of BMD-associated genetic variants could improve the accuracy of fracture prediction over and above that of clinical risk factors alone.

Perspective

Despite fruitful GWAS discoveries in the bone field,most of these GWAS participants were of European descent.In fact,if we extended to other complex traits,~79%of GWASs were conducted in European populations according to the GWAS catalog.Martin et al.142systematically evaluated the polygenic risk prediction accuracy in Japanese,British and African-descent individuals on the basis of using independent GWASs of equal sample sizes from BioBank Japan(BBJ)and UK Biobank,including 17 quantitative anthropometric and blood panel traits and five disease endpoints;they demonstrated that prediction accuracy was consistently higher with GWAS summary statistics from ancestry-matched summary statistics.The condition of genetic resources and analyses overwhelmingly centered on individuals of European ancestry would lead to imbalances in the subsequent translatability of findings.To realize the full and equitable potential of the polygenic risk score,it was encouraged that more GWASs and sequencing studies on osteoporosis,BMD and fracture should be carried out in additional ethnic populations,such as the Chinese population,which made up~20% of the global population.Fortunately,the cost of wholegenome sequencing and genotyping has dramatically decreased,making the utility of genetic variants more affordable and practical.In addition,the prioritization of the recruitment and analysis of diverse cohorts would become smooth with an increasingly globalized and connected research community.143

In summary,the achievement of GWASs is unprecedented in the understanding of how genetic variants influence osteoporosis and fracture.In the future,by mining large databases with detailed characterization of relevant phenotypes,more causal genes/mutations will be identified.In addition,large-scale genetic data could provide a new way to identify new drug targets and could be translated into precision treatment options to prevent and treat osteoporosis and fracture.

ACKNOWLEDGEMENTS

This study was supported by the National Natural Science Foundation of China(81871831 and 32061143019).The funding agencies had no role in the study design,data collection and analysis or the decision to publish or prepare the manuscript.We thank the peer reviewers for their thorough and helpful review of this manuscript.

We also thank the High-Performance Computing Center at Westlake University for the facility support and technical assistance.

ADDITIONAL INFORMATION

Supplementary informationThe online version contains supplementary material available at https://doi.org/10.1038/s41413-021-00143-3.

Competing interests:The authors declare no competing interests.