Introduction

Osteoporosis is a widespread systemic disorder characterized by decreased bone mass and increased risk of fracture. In Europe, the incidence of osteoporosis is increasing exponentially.1 Genetic factors have a role in the regulation of bone mass.2 Changes in bone turnover and the occurrence of fractures are associated with a number of single-nucleotide polymorphisms (SNPs) in genes coding for bone-related proteins, but disappointingly the effect sizes have been lower than anticipated.3 Osteoporosis is a multi-factorial disease that results from a combination of genetic, environmental and other risk factors, such as tobacco use, sedentary lifestyle, low calcium intake, low body weight, certain medications (glucocorticoids), and so on, that affect peak bone mass and rate of bone loss.

The P2X7 purinergic receptor belongs to a family of ligand-gated cation channels.4 In contrast to other P2X receptors, prolonged activation with ligand induces formation of a pore permeable to cations up to 900 Da.5 P2X7 is expressed primarily by cells of the hematopoietic and immune system,6, 7, 8, 9, 10 osteoclasts,11 osteoblasts12, 13 and osteocytes.14 Activation of P2X7 has been shown to induce membrane blebbing (in osteoclasts),15 metalloproteinase activation (in lymphocytes),16 and following prolonged activation, apoptosis and cell death (osteoclasts).17 In addition, P2X7 activation is a major stimulus for formation and activation of the inflammasome 18 and release of interleukin (IL)-1a and IL-1b.10, 19, 20 This is particularly relevant in relation to post-menopausal bone loss, where estrogen depletion results in increased production and release of IL-1, and subsequently osteoclast formation and bone resorption.21

It has been demonstrated that P2X7 is involved in mechanically induced signalling between osteoblasts and osteoclasts,14 as well as with osteoclast survival.12, 17, 22 Furthermore, it has been shown to couple to NF-κB in osteoclasts.23 In a P2X7 knock-out model, P2X7-null mice have decreased bone mass24, 25 and decreased responses to mechanical stimulation of the bones,26 supporting a role for the P2X7 receptor in skeletal mechanotransduction. Thus, defects in normal P2X7 receptor function in humans could potentially have detrimental effects on bone mass, quality and skeletal resistance to mechanical stress.

The P2X7 receptor gene, P2RX7, is located on chromosome 12q24. It is highly polymorphic with 40 coding variants reported in the latest build 132 of dbSNP database with more than a dozen variants characterized for their effect on function.27 At least eight non-synonymous variants in P2RX7 confer loss-of-function in the P2X7 signalling pathway. Three are uncommon (1–3% minor allele frequency (MAF)) but give rise to severe functional defects.28 In particular, one of these three variants, Arg307Gln, which is located in the ATP binding pocket, dramatically alters receptor function by reducing the affinity of agonist binding.29 Only two variants (His155Tyr and Ala348Thr) have been shown to confer gain-of-function. However, while the inheritance of the Gln460Arg minor allele in isolation causes insignificant loss-of-function, the 460-Arg is nearly always co-inherited with 155Tyr and 348Thr and thus is a marker of the major gain-of-function haplotype of P2RX7 in Caucasians.30

In this study, we performed an association analysis of non-synonymous variants in P2RX7 and changes in bone mineral density (BMD) and fracture rates over 10 years in post-menopausal women. A risk model based on the functional effects of the individual variants on receptor function showed a clear association to the regulation of bone mass. Thus, the P2X7 receptor is important in regulating BMD in the ageing human skeleton.

Participants and methods

Study participants

The study population consisted of healthy, post-menopausal Danish women included in the Danish Osteoporosis Prevention Study (DOPS). They were enrolled 3–24 months after menopause and aged 45–58 years at inclusion. The recruitment procedure and study design have been described previously.31 In brief, DOPS is a prospective multi-center trial investigating the preventive effect of hormone replacement therapy (HRT) on the development of osteoporotic fractures. A total of 2016 women were included and they were allocated to a randomized and a non-randomized arm. In both arms, study participants were allocated to either HRT or non-HRT with a randomized open design in the randomized arm while participants included in the non-randomized arm were given the possibility of choosing between HRT and non-HRT. From the total number of 2016 women, 1795 gave informed consent to donate DNA, which was stored at −80 °C. However, as genotyping was performed in duplicate at two sites, and individuals where discrepant results were found, were excluded from the analyses a total of between 1651 and 1694 individuals were available for analysis for the individual variants.

To determine the association between the individual genotypes and BMD and fracture rate in HRT naïve women and the response to HRT, the 5 and 10 years longitudinal analyses were performed separately for the two groups. At 5/10 years, 876/768 in the HRT-naïve group and 418/121 in the HRT-group were evaluable. The main reason for the decrease in evaluable individuals is that study participants were allowed to change treatment during the course of the study. All procedures involving patients and patient materials were approved by the Danish Ethics Committee (approval #1990/1821), and all participants were recruited after obtaining informed consent.

Genotyping

Non-synonymous P2RX7 variants, which were selected for a known effect on P2X7 receptor function, were analyzed using TaqMan allelic discrimination assays (Department of Endocrinology at Aarhus University Hospital, Aarhus, Denmark) and a homogeneous mass extension assay (HME) at the Australian Genome Research Facility (St Lucia, Queensland). Samples that failed HME were re-analyzed using restriction enzyme digestion of appropriate PCR products or by Taqman assay as previously described.32, 33 Homogenous mass extension method failed in 3–4% of the samples and for these samples genotyping was repeated using real-time PCR. To determine accuracy, genotyping results were compared between the two genotyping methods. Discrepancies were found for 79 individuals at one or more variants. There was a discrepancy rate between the two genotyping methods of <1% and all discrepant results were removed from the final analysis. In 40 individuals, genotyping failed for eight or more variants, and these subjects were excluded from the analyses. Haplotype analysis was performed using Phase 2.0 software.34

Biochemistry

At baseline, blood samples were collected in the morning after an overnight fast. Serum was isolated and stored at −80 °C for later analysis. Bone-specific alkaline phosphatase (U/l) was lectin-precipitated and analyzed using spectrophotometry with an intraassay coefficient of variation (CV) of 8% and an interassay CV of 25%. Serum osteocalcin was analyzed by radio-immuno assay as described previously.31 Intraassay CV was 5% and interassay CV was 10%. Hydroxyproline (OHP), was measured using fasting second void urine and analyzed by spectrophotometry according to the manufacturer's instructions (Organon Teknika, Boxtel, The Netherlands). Values were expressed as a ratio relative to creatinine (Cr) excretion (OHP/Cr, μmol/mmol).

Bone mass measurements

At inclusion, and after 5 and 10 years, BMD (g/cm2) was assessed using dual-energy X-ray absorptiometry for all study participants. BMD was determined at the lumbar spine, the hip and the femoral neck using cross-calibrated Hologic QDR-1000/W and Hologic QDR-2000 densitometers (Hologic Inc., Waltham, MA, USA).

Fracture determination

At inclusion, and after 5 and 10 years, X-rays of the vertebral spine were obtained as lateral projections covering T4 to L5. The presence of vertebral fractures was evaluated by a trained radiologist and a fracture was defined as >20% reduction in the height at the anterior, middle or posterior of the vertebrae. Of the 1680 having a vertebral X-ray performed at the 10 years visit, 143 were diagnosed as having one or more vertebral fractures. After 10 years, genotypes were available for 799 individuals in the non-HRT group included in the fracture study, of which 71 had sustained one or more vertebral fractures. Of these, 51 subjects had one vertebral fracture and the remaining 20 had >1. In the HRT group, genotypes were available for 124 individuals. Six had one or more vertebral fractures diagnosed on X-ray, and of these three had one fracture and three had more than one fracture.

Statistical analysis

Distribution of genotype frequencies was tested for adherence to Hardy–Weinberg equilibrium (HWE) using the χ2-test and Fisher's exact test. Differences in BMD or rate of bone loss were tested by one-way analysis of variance (ANOVA) after testing for normal distribution of the data. Results are expressed as mean±SD. Differences in covariates between genotypes were tested by ANOVA and Kruskal–Wallis test. The influence of covariates on baseline BMD values were tested using the general linear model procedure for covariance analysis (ANCOVA) and BMD values were corrected for age, body mass index (BMI), and time following menopause. For fracture analysis, differences in fracture incidences were tested using the χ2-test. For all analyses, P-values <0.05 were considered statistically significant. All statistic calculations were performed using the SPSS statistical software, version 11.5 (SPSS Inc., Chicago, IL, USA).

Results

P2X7 genotypes

The study population was genotyped for 12 different non-synonymous SNPs in the gene coding for the P2X7 receptor. The variants were selected based on published effects on receptor function (Table 1). Furthermore, the MAFs are shown in Table 1. All variants were in HWE, except for the Gly150Arg (P<0.001).

Table 1 List of P2X7 receptor single-nucleotide polymorphisms for which the participants in the Danish Osteoporosis Prevention Study have been genotyped

Baseline characteristics at menopause of the P2X7 genotypes

First, age, height, BMI and bone markers were compared for the individual genotypes. No statistically significant differences were found between genotypes. Next, baseline values of BMD were compared. Interestingly, significant differences between genotypes were found for both the Val76Ala and the Gln460Arg polymorphism (Supplementary Material Table). However, after correcting for multiple analyses, only the Gln460Arg polymorphism remained significantly associated to BMD at the femoral neck (P=0.01) with the homozygotes for the variant alleles having the higher BMD. Although not significant, this trend was also present for the two other sites of measurement.

P2X7 genotype and bone loss over 5 and 10 years after menopause

An accelerated bone loss takes place after menopause and the rate of bone loss was expressed as change of BMD per year over 0–5 and 0–10 years for each genotype of the 12 functional variants in P2RX7. For women not treated with HRT at any time, BMD loss at the hip was greatest in those heterozygous for Arg307Gln (Figure 1a and Table 2). A similar trend was observed in both lumbar spine and femoral neck at 0–5 and 0–10 years. In subjects receiving HRT treatment, rates of bone loss were much reduced in all P2X7 genotypes. However, even in this HRT group, the variant allele Arg307Gln, was still associated with more bone loss over 0–5 years than those who were wild type for this allele (Figure 1b). Owing to low numbers in the individual genotype groups, it was not possible to make similar analyses for the 0–10 year bone loss. These results suggest a major effect of the P2X7 loss-of-function 307-Arg allele in accelerating bone loss in post-menopausal women.

Figure 1
figure 1

Relationship of P2X7 SNPs with changes in BMD. (a) Post-menopausal change in BMD in non-HRT individuals for the total loss-of-function SNP Arg307Gln. Individuals carrying the mutant allele have an increased rate of bone loss. This is significant for the total hip (TH) 10 years after menopause, and the same trend is present at the LS and FN (top panel). Furthermore, a similar trend is seen after the first 5 years after menopause (bottom panel). The number of individuals with the two genotypes was 860 GG and 16 GA after 5 years and 758 GG and 15 GA after 10 years. No study participants were homozygous for the variant A allele. (b) The presence of Arg307Gln variant, increases the rate of bone loss at the femoral neck (0–5 years) despite HRT treatment. (c) Changes in BMD associated with the Ile568Asn variant over 5 years post menopause. Individuals homozygous for the variant allele (AA) have significantly accelerated bone loss in the femoral neck, with the same trend seen at TH and LS. *P-value <0.05; LS, lumbar spine; FN, femoral neck; NS, not significant; TH, total hip.

Table 2 Change in bone mineral density (BMD) for individuals not on hormone replacement therapy

Another functional variant in P2RX7 had significant effects on the rate of bone loss. The Ile568Asn allele was significantly associated to the rate of bone loss from 0 to 5 years after menopause in the non-HRT group (Figure 1c) with two subjects homozygous for the variant allele (AA) clearly having the highest rate of bone loss of all the P2X7 genotypes. However, analysis of the rate of bone loss from 0 to 10 years after menopause showed that this variant was not significantly associated to the bone loss, although numerically the trends were the same as from 0 to 5 years (data not shown). To demonstrate that the findings were not due to differences in basic parameters at baseline for the subgroups not on HRT at 5 and 10 years, we compared the baseline values for age, height, BMI and BMD, and found no differences between genotypes for those two subgroups.

As the frequency of the rare allele for the four polymorphisms Val76Ala, Gly150Arg, Thr357Ser and Ile568Asn are very low, we performed the above analyses again after grouping individuals carrying the rare allele (heterozygotes and homozygotes for the rare allele). However, no additional associations could be demonstrated (not shown).

P2X7 genotype and fracture risk after 5 and 10 years after menopause

We examined the distribution of vertebral fractures for the 12 P2X7 variants. In a previous analysis of our data, the common Glu496Ala variant was associated with vertebral fracture incidence.17 In addition to this, we found a significant association of the gain-of-function polymorphism Ala348Thr to fracture incidence, where individuals homozygous for the variant allele (AA) and thereby with gain of receptor function were partially protected against vertebral fractures as the incidence was reduced by almost 70% compared with subjects homozygous for the wild-type allele (GG) (Figure 2).

Figure 2
figure 2

Association of fracture incidence with Ala348Thr variant in non-HRT subjects. A χ2-test was performed to test for differences in fracture incidences between the three genotype groups 0–10 years post menopause. Significantly fewer fractures are associated with homozygous mutant (AA) group than the wild type (GG) or heterozygous (GA) group (P=0.04). The number of individuals with the genotype is shown below the bars. The number of individuals having sustained a fracture is shown above each bar.

To determine whether the detected differences could be attributed to differences in confounding factors such as age, time since menopause, daily calcium intake, intake of vitamin D, BMI, comparisons by ANOVA were performed among the three genotypes for the Ala348Thr polymorphism. No differences between groups for any of the confounding factors could be demonstrated. A similar analysis between fracture rate and Ala348Thr in the HRT group failed to show significant association (data not shown).

Again, we grouped individuals carrying the rare allele (heterozygotes and homozygotes for polymorphisms Val76Ala, Gly150Arg, Thr357Ser and Ile568Asn). Now, as previously reported,17 we found an association between reduced risk of vertebral fracture and the rare allele of the Ile568Asn polymorphism. For the remaining three rare polymorphisms, no association to fracture risk could be demonstrated.

P2X7 risk model analyses – bone loss 5 and 10 years after menopause

The results above suggest that several of the uncommon variants confer a predisposition to accelerated loss of BMD. To further test this hypothesis, we examined a risk model based on combining subjects into one of three risk groups according to the functional effects (Table 1) of the P2RX7 variants. The groupings are as follows: (A) A ‘high-risk’ group consisting of individuals carrying one or more of the variant alleles for any of the five total loss-of-function variants (151+1g->t, Gly150Arg, Arg276His, Arg307Gln and Ile568Asn). As several of these alleles exert a dominant-negative effect if expressed together with a gain-of-function allele, we allowed inclusion of a co-inherited gain-of-function allele. (B) A ‘low-risk group’ consisting of individuals homozygous for the Gln460Arg polymorphism, as we have recently within the research consortium shown that it is an indicator of the most pronounced increase in function of the receptor in primary cells.30 These individuals are homozygous for the wild-type allele of all five total loss-of-function variants detailed in group A and thus do not carry any of the mutant alleles. (C) An ‘intermediate-risk’ group consisting of all other individuals.

As for the individual variants, we compared the baseline values of age, height, BMI and markers of bone metabolism for the three risk-groups. No statistically significant differences were found between genotypes (Table 3). Next, baseline values of BMD at the lumbar spine, femoral neck and total hip were compared using analysis of covariance, with correction for logBMI, age and time since menopause, and again no significant differences between groups could be detected (Table 3).

Table 3 Baseline characteristics of the three risk-groups based on P2X7 genotype

We then tested whether there were any differences in rate of bone loss between the three risk-groups for the untreated individuals. We found a highly significant difference between the risk-groups for BMD loss in the lumbar spine from baseline to 10 years after menopause. Individuals belonging to the high-risk group lost more than twice as much bone per year as the individuals in the low-risk group with the intermediate-risk group showing bone loss at a rate between these two extremes (Figure 3). Interestingly, the trend was the same for 0 to 5 years, although not significant, while a significant difference in bone loss between risk groups from 5 to 10 years was found (P=0.02) clearly confirming that intact P2X7 receptor function is important for maintaining a strong and healthy skeleton after menopause. Haplotype analysis of P2X7 genotypes identified 17 haplotypes with population frequencies >1.0 % (Table 4). Grouping the haplotypes based on the haplotype blocks containing the 1068, 1096, 1405 and 1513 polymorphisms as described by Stokes et al30 showed similar haplotype frequencies between the Australian cohort and the DOPS cohort. Aggregating patients into a group of haplotypes carrying ‘high-risk’ loss-of-function variants showed the same accelerated BMD loss as for the Group A defined above on the basis of SNPs (data not shown).

Figure 3
figure 3

The rate of post-menopausal bone loss at the LS is associated with P2X7 functional groupings. In the absence of HRT treatment women who have one or more P2X7 loss-of-function SNP (high-risk group) have significantly greater loss in BMD than women who have P2X7 gain-of-function variants (low-risk group) and an intermediate group consisting of all other individuals. **P<0.01.

Table 4 The table shows the 17 most frequent haplotypes (H1–H17) with a frequency of more than 1% in the cohort together with their frequency and the predicted functional effects according to reference35

Finally, the fracture risk between the risk-groups was determined using χ2-test. Although no significant differences could be detected (data not shown), it should be noted that the small number of individuals in each group reduced the power of the analysis.

Discussion

This 10-year study shows that a SNP in the P2RX7 gene conferring severe loss-of-function on the P2X7 receptor is associated with accelerated bone loss in women following menopause. This variant (Arg307Gln) is located in a coding region of the gene substituting an uncharged glutamine in place of arginine-307 in the ATP-binding site of the receptor thereby abolishing activation by agonist.29 Subjects who were heterozygous for this Arg307Gln variant had >40% greater bone loss at the hip over the 0–10 years interval than subjects who were wild type at this position. Estrogens have a strong anabolic effect on bone mass and in this study estrogens reduced or even reversed bone loss in all P2RX7 subgroups. However, the effect of the Arg307Gln genotype to accelerate bone loss at the femoral neck was evident even in subjects taking HRT continuously over the 0- to 5-year period of study. Many of these subjects subsequently stopped HRT so that numbers at 10 years became too small for analysis of this effect. This observation suggests that the functional effects of polymorphisms in the P2RX7 and estrogen receptor (ESR1) genes have independent effects on bone mass. The association of low bone mass with the Arg307Gln variant has been replicated in the Aberdeen Prospective Osteoporosis Screening Study who studied 506 post-menopausal women over 6- to 7-year follow-up.35 Subjects with the Arg307Gln variant had lower BMD of the lumbar spine both at baseline and at follow-up, suggesting this variant of the P2X7 receptor may be exerting an effect on bone mass even before the menopause.35 The Arg307Gln has also been associated with other bone-related events. A recent study showed an association of the Arg307Gln polymorphism with increased failure rate for total hip replacement,36 although this is not necessarily the result of decreased bone mass. Another polymorphism (Ile568Asn) in P2X7 showed a significant association with changes in bone mass over 0–5 years. Subjects who were heterozygous for Ile568Asn had more bone loss than those with wild-type Ile-568, but as previously reported the difference did not reach significance.17 Two subjects, however, were homozygous for 568-Asn, which totally prevents trafficking of P2X7 to the cell surface,37 and these two suffered major loss of BMD. Extensive genotyping in 3430 Caucasian subjects has defined two haplotypes (P2X7-2 and -4), which exhibit gain-of-function over wild type (P2X7-1). Both P2X7-2 and P2X7-4 contain the Ala348Thr polymorphism and this variant is critical for the gain-of-function effect.30 Analysis showed a lesser fracture rate for subjects with one or two of this gain-of-function allele.

Interestingly, we were able to identify three groups of individuals in the cohort with highly different risk of developing osteoporosis based on their combined genotypes and the knowledge of functional effects of these genotypes on P2X7 receptor function. The groups with the lowest risk of developing osteoporosis were the ones homozygous for the variant allele of the Gln460Arg polymorphism. However, in vitro studies have shown that the variant allele in isolation reduces receptor function slightly.45 In contrast cells transfected with the full haplotype block has significantly increased function.30 Thus, the variant allele of Gln460Arg is therefore likely to be co-inherited with other gain-of-function polymorphisms in the P2X7 receptor, which subsequently are associated with low risk of bone loss.

Recently, it was shown that the P2X4 receptor can be co-expressed together with the P2X7 receptor and that P2X4 and P2X7 subunits can form heteromeric channels.38 If this is the case in vivo, P2X4 receptor subunits might theoretically replace ‘defect’ P2X7 receptor subunits in the cell membrane thereby reducing the true effect of a genetic defect in the gene encoding the P2X7 receptor. Furthermore, the P2X4 gene is situated close to and downstream for the P2X7 receptor on the same chromosome. Thus, polymorphisms in the P2X7 gene might be in linkage with polymorphisms in the P2X4 receptor gene. The associations we find might therefore also be a consequence of gain- or loss-of-function polymorphisms in the P2X4 gene. Thus, it is highly warranted to investigate the association of P2X4 polymorphisms with bone status and also the interaction between P2X4 and P2X7 polymorphisms.

In this study, loss-of-function polymorphisms in the human P2RX7 gene have produced analogous changes to deletion of the P2RX7 gene in mice, which gave rise to animals with reduced BMD in weight-bearing limbs because of reduced periosteal bone formation and increased trabecular bone resorption.24 Similar skeletal changes have been described in humans with disuse of weight-bearing limbs39 raising the possibility that the purinergic P2X7 receptor forms part of a signaling pathway, which transduces mechanical strain into an osteoblastic response and improved bone strength. Studies of osteoblasts and osteocytes in vitro have shown that fluid shear stresses increase the production of prostaglandin E2 (PGE2) in wild type but not in P2X7 –/– cells26 and PGE2 has been shown to have major anabolic effects on bone formation. However, the first study to demonstrate a direct link between P2X7 receptor activation and osteogenesis was the study by Dixon's group showing that activation of P2X7 receptors in osteoblasts gives rise to the lipid mediators, lysophosphatidic acid and PGE2,40 both of which have been reported to regulate osteogenesis. Mechanical stimulation leads to a release of nucleotides, which activate P2X7-mediated apoptosis of osteoclasts41 and thus lead to increased bone strength; an effect that would augment the P2X7-induced osteogenesis.

We find that one of the polymorphisms (Gly150Arg) is not in HWE. As we have genotyped the samples with two different methods independently, genotyping error is not very likely. Other reasons could be high mutation rate, selection toward one of the genotypes, assortative mating, and so on. However, none of these seem very likely, as only this polymorphism is not in HWE. Therefore, this is likely to be a chance finding.

In conclusion, one uncommon P2X7 receptor variant, Arg307Gln, exerted a dominant-negative effect on P2X7 function and heterozygotes were associated with greater bone loss over 5 or 10 years than wild types. A second uncommon variant, Ile568Asn, was present in homozygous dosage in two subjects who exhibited the highest rate of bone loss in our study. Our results for the Arg307Gln variant are replicated in an independent cohort studied over 6 years in Aberdeen, UK.35 and show that uncommon genetic variants of the P2X7 receptor are an important risk factor in osteoporosis of post-menopausal women.