A Novel Pathophysiological Mechanism for Osteoporosis Suggested by an in Vivo Gene Expression Study of Circulating Monocytes*

Bone mineral density (BMD) is a major risk factor for osteoporosis. Circulating monocytes may serve as early progenitors of osteoclasts and produce a wide variety of factors important to bone metabolism. However, little is known about the roles of circulating monocytes in relation to the pathophysiology of osteoporosis. Using the Affymetrix HG-U133A GeneChip® array, we performed a comparative gene expression study of circulating monocytes in subjects with high and low BMD. We identified in total 66 differentially expressed genes including some novel as well as some already known to be relevant to bone metabolism. Three genes potentially contributing to bone metabolism, CCR3 (chemokine receptor 3), HDC (histidine decarboxylase, i.e. the histamine synthesis enzyme), and GCR (glucocorticoid receptor), were confirmed by quantitative real-time reverse transcriptase-PCR as up-regulated in subjects with lower BMD. In addition, significant negative correlation was observed between expression levels of the genes and BMD Z-scores. These three genes and/or their products mediate monocyte chemotaxis, histamine production, and/or sensitivity to glucocorticoids. Our results suggest a novel pathophysiological mechanism for osteoporosis that is characterized by increased recruitment of circulating monocyte into bone, enhanced monocyte differentiation into osteoclasts, as well as osteoclast stimulation via monocyte functional changes. This is the first in vivo microarray study of osteoporosis in humans. The results may contribute to identification of new genes and their functions for osteoporosis and suggest genetic markers to discern individuals at higher risk to osteoporosis with an aim for preventive intervention and treatment.

Osteoporosis is most frequently characterized by low bone mineral density (BMD). 1 The well accepted pathophysiological mechanisms for low bone mass include prolongation of the life span of osteoclasts (1) and early apoptosis of osteoblasts (2) and osteocytes (3,4). Recently, the imbalance between osteoblastogenesis and adipogenesis of bone marrow mesenchymal stem cells was suggested as a new mechanism for osteoporosis (5). However, the above mechanisms were discovered mainly through in vitro studies or in vivo studies using animals. Their comprehensiveness and exhaustiveness is unclear, and their direct relevance to osteoporosis in humans is yet to be established.
Circulating monocytes may serve as progenitors of osteoclasts (6 -9). Monocytes isolated from peripheral blood can differentiate in vitro into the multinuclear cells that express typical osteoclast markers and are able to resorb bone (6,8). All osteoclasts in peripheral skeleton (7,10) and a considerable amount of osteoclasts in the central skeleton (11) come from circulating monocytes. Blood monocytes also produce a wide variety of factors involved in bone metabolism, such as interleukin-1, tumor necrosis factor-␣, interleukin-6, platelet-derived growth factor, transforming growth factor-␤, and 1,25(OH) 2 D 3 (12). Circulating mononuclear cells of osteoporotic patients, when induced in vitro into osteoclasts, demonstrated increased levels of bone resorption activity as compared with normal controls (13).
For both peripheral and central skeleton, the access route of circulating monocytes into the bone microenvironment is through the structure called sinusoid (14) that occupies the center of a basic multicellular unit, the functional unit of bone remodeling. Recruitment of circulating monocytes into bone is modulated mainly by chemotaxis (15), whereby monocyte chemoattractants from bone and their receptors on monocytes, such as chemokine receptors, interact to initiate and maintain a directed migration of monocytes toward bone tissue.
Given the importance of the circulating monocytes in bone metabolism, currently, surprisingly little is known about the in vivo roles of circulating monocytes in the pathophysiology of osteoporosis. Microarrays are a powerful tool to provide a comprehensive picture of cell function as they can assay expression of tens of thousands of genes simultaneously; however, they have not been attempted for in vivo studies in humans to investigate the molecular mechanisms underlying osteoporo-sis. Here, with an in vivo approach directly in humans, we explore other potential mechanisms for osteoporosis, which are mediated by or reflected in circulating monocytes.
We hypothesized that gene expression profiles of circulating monocytes differ between high and low BMD women. These differentially expressed genes (DEGs) are potentially important for the pathogenesis of osteoporosis. To identify these genes, we applied the Affymetrix HG-U133A GeneChip® array to analyze the gene expression profiles of circulating monocytes in subjects with extremely discordant BMD values.

Subjects
The study was approved by the Creighton University Institutional Review Board. All the study subjects signed informed consent documents before entering the project. We established strict criteria to exclude non-genetic factors that might cause BMD variation. Additional criteria were implemented to exclude diseases/conditions that may potentially lead to gene expression changes of circulating monocytes. The detailed criteria are included in the Supplemental Data.
The recruited women were all Caucasians and aged between 47 and 55 with an average age of ϳ51. The narrow age window may minimize the potential age effects on monocyte gene expression (16,17). Among the 19 subjects, 10 had high spine BMD, with an average L1-L4 Z-score of ϩ2.47, and 9 had low spine BMD, with an average L1-L4 Z-score of Ϫ1.43. The Z-score is defined as the number of standard deviations of a BMD measurement above (i.e. a positive Z-score) or below (i.e. a negative Z-score) the age-, gender-, and ethnicity-matched mean BMD. The average Z-score for the high BMD subjects places them in the top 0.68% of the normal population distribution, whereas the Z-score for the low BMD subjects places them in the bottom 7.6% of the normal population distribution. Among the 10 high BMD subjects, five were premenopausal, and five were postmenopausal. Among the nine low BMD subjects, four were premenopausal, and five were postmenopausal. In our study, menopause was defined as cessation of menstruation for at least one year. The recruitment scheme for menopausal status ensures that the DEGs found between the 10 high and 9 low BMD subjects are only related to BMD. Detailed characteristics of the subjects are available in Table I.
We used a method recommended by Affymetrix to evaluate our sample size (e.g. the number of replicates per condition) (18). According to the method, coefficient of variation (CV) is calculated for certain percentiles (e.g. the 25th, 50th, and 75th) of signal intensities for a data set containing 1, 2, 3, . . ., and n replicates. Comparison is made between the calculated CVs of certain percentiles for different numbers of replicates. The optimal sample size is achieved when the value of CV is stabilized, that is, does not change appreciably from n replicates to n ϩ 1 replicates. Under such a situation, it is unlikely that additional replicates will improve significantly the accuracy of the estimation for the standard deviations of the sample, which is ultimately used to determine statistical significance in parametric statistical tests.
Based on the calculations with the raw microarray data in our study, we observed that the CV became stabilized at seven replicates for the 25th and the 50th percentile of signal intensities and at nine replicates for the 75th percentile of signal intensities (Fig. 1). Therefore, our sample size is sufficient and decent to provide adequate statistical power for detecting differential gene expression.

BMD Measurement
BMD (g/cm 2 ) for the lumbar spine (L1-L4) was measured with a Hologic 4500 dual energy x-ray absorptiometer scanner (Hologic Corp., Waltham, MA). The machine was calibrated daily. The coefficient of variation of the dual energy x-ray absorptiometer measurements for BMD was 0.9%. None of our subjects had any known conditions that might artificially increase BMD values, such as osteophytes and facet sclerosis.

Experimental Procedures
Monocyte Isolation-Monocyte isolation from 30 ml of whole blood was performed with a monocyte-negative isolation kit (Dynal Biotech Inc., Lake Success, NY) following manufacturer's recommendation. The kit contains a highly optimized antibody mix, blocking reagent, and Depletion Dynabeads® to deplete T cells, B cells, and natural killer cells from mononuclear cells, leaving monocytes untouched and free of the surface-bound antibody and beads. Our flow cytometry results showed an average monocyte purity of 86% Ϯ 3% in the samples.
Total RNA Extraction-Total RNA from monocytes was extracted using Qiagen RNeasy Mini kit (Qiagen, Inc., Valencia, CA).
Preparation of cRNA and Gene Chip Hybridization-Preparation of cRNA, hybridization, and scanning of the HG-U133A GeneChip® was performed according to the manufacturer's protocol (Affymetrix, Santa Clara, CA).
Real-time RT-PCR-Real-time RT-PCR was performed in two steps, i.e. the first step of reverse transcription for synthesis of cDNA from total RNA and the secondary step of real-time quantitative PCR.
Reverse transcription reactions were performed in a 30-l reaction volume, containing 3 l of 10ϫ PCR Buffer II (100 mM Tris-HCl, pH 8.3, 500 mM KCl, 0.01% w/v gelatin), 6 l of 2.5 mM MgCl 2 , 8 l of dNTPs, 1 l of murine leukemia virus reverse transcriptase, 1 l of RNase inhibitor, 1 l of oligo(dT), 0.5 g of total RNA, and water to 10 l. All of the above reagents were supplied by Applied Biosystems (Foster City, CA). Reaction conditions were as follows, 10 min at 20°C, 60 min at 42°C, 5 min at 95°C.
Real-time quantitative PCR was performed in a 50-l reaction volume using standard protocols on an Applied Biosystems 7000 sequence detection system. Briefly, 5 l of cDNA was mixed with 25 l of TaqMan universal PCR master mix (2ϫ), 2.5 l of Assays-on-Demand TM gene expression assay mix (contains forward and reverse primers and labeled probe), 2.5 l of human GAPDH (20ϫ), and 15 l of water. The thermocycling conditions were as follows: 2 min at 50°C, 10 min at 95°C, 50 cycles of 15 s at 95°C plus 1 min at 60°C. The thermal denaturation protocol was run at the end of the PCR to determine the copy number of products that were present in the reaction. All reactions were run in triplicates and included "no template" controls for each gene.
The TaqMan gene expression assays all have amplification efficiencies very close to 1 (19). Thus, it was not necessary for us to perform validation experiments to test the equality of the amplification efficiencies between the target genes and the reference gene (GAPDH).

Data Analyses
We used a method by Golub et al. (20) to preselect a subset of genes from the total of around 14,500 genes in the array for further analyses. The preselection process minimizes the multiple-testing problem by excluding from further analysis the genes with relatively constant expression across the samples. It was based on the expression values generated with the MAS 5.0 software (Affymetrix). We selected a subset of 8,623 probe sets from a total of 22,283 for around 14,500 genes. Subsequent data analyses were performed on these selected probe sets using Bioconductor packages (21). The Bioconductor Affy package (22) was used to process the probe-level data (raw data) and transform the data into expression data using the robust multiarray average algorithm (23).
Based on the expression data generated with the robust multiarray average algorithm, the Bioconductor Multtest package was used to identify DEGs between high versus low BMD subjects. To account for multiple testing, the Benjamini and Hochberg stepwise procedure (24), which controls for the false discovery rate, was implemented to generate adjusted p values.

RESULTS
DEGs Suggested by Microarray Analyses-All the raw microarray data are available. 2 We also submitted the data to the NCBI Gene Expression Omnibus data repository under accession number GSE2208. Based on the analysis with the MAS 5.0 software (Affymetrix), on average, 39.69% of a total of 22,283 probe sets in the array were called "present" for our samples.
Of the preselected 8,623 probe sets, we identified 67 (corresponding to 66 genes) differentially expressed between the high and low BMD subjects (Benjamini and Hochberg stepwise procedure-adjusted p Ͻ 0.05). These genes are listed in the order of the adjusted p values in Table 1 of the Supplemental Data. Forty-six genes are up-regulated and 20 are down-regulated in the low BMD subjects as compared with the high BMD ones. Since the majority of the genes listed in the table have not been confirmed with real-time RT-PCR, the significance of their differential expression (i.e. the p values) should be interpreted with caution.
DEGs Validated through Real-time RT-PCR-The number of monocytes we obtained from the 30 ml of blood of each subject was limited, and the total RNA was scarce. Therefore, our criteria for selection of genes for further real-time RT-PCR validation were stringent and comprehensive. They were based not only on the statistical significance of the differential expression of a gene as detected by microarray (i.e. the Benjamini and Hochberg stepwise procedure-adjusted p Ͻ 0.05) but also on its potential functional relevance to bone metabolism. Accordingly, we selected six genes for real-time RT-PCR validation, which are CCR3 (chemokine (C-C motif) receptor 3), HDC (histidine decarboxylase), GCR (glucocorticoid receptor), RAMP3 (receptor (calcitonin) activity-modifying protein 3), SCYE1 (small inducible cytokine subfamily E, member 1), and FCER1A (Fc fragment of IgE, high affinity I, receptor for ␣-polypeptide). These genes and/or their products, through modulation of monocytes, are known to participate in bone monocyte recruitment (CCR3, HDC) (15,25), osteoclastogenesis (GCR, HDC, CCR3) (26 -31), osteoclast regulation (HDC, RAMP3) (32,33), and/or production of proinflammatory cytokines important to osteoporosis (SCYE1, FCER1A) (34,35).
Based on the relative gene expression (i.e. 2 Ϫ⌬C T , where ⌬C T ϭ (C T Target Gene Ϫ C T GAPDH )) for each sample obtained through real-time RT-PCR, the Student's t test was employed to compare expression for each gene between the high and low BMD subjects. As our subjects were also stratified by menopausal status, comparison of gene expression between high and low BMD subjects was also performed among premenopausal as well as postmenopausal women. Since the major aim of this study was focused on identifying genes related to BMD variation, factorial analysis for identifying differential gene expression due to menopause and/or interaction between BMD and menopause was not performed here and was pursued elsewhere (36). The t tests confirmed the differential expression for three genes. CCR3 and HDC were up-regulated, by 3-(p ϭ 0.0002) and 4-fold (p ϭ 0.022), respectively, in the low versus the high BMD subjects, and GCR was 2-fold up-regulated in the premenopausal low versus high BMD subjects (p ϭ 0.027).
Although there was significant difference in expression between the high and the low BMD groups for the three genes, there is no clear threshold level for their expression, which might be associated with the low BMD subjects. Therefore, the expression levels may represent continuous variables, which negatively correlate with BMD.
The reason for only three out of the six candidate genes being confirmed by real-time RT-PCR might be related to the absolute expression levels of the genes. Genes with moderate expression levels tend to have higher correlation between expression quantification results obtained by microarray and RT-PCR (37). Indeed, in our study, among the six selected genes, the three confirmed ones had expression levels within the middle range (i.e. from 800 to 1,300 according to microarray quantification using the MAS 5.0 software), whereas the other three had more extreme expression levels (i.e. less than 300 or more than 3,000). Statistically, this may also be due to problems associated with multiple testing and statistical power for replication (see detailed elaboration under "Discussion").
We chose GAPDH as the endogenous control gene for the real-time RT-PCR because its expression levels in microarrays were very consistent and stable across all the samples. In the HG-U133A array, there are three independent probe sets, among a total of 22,283, for the GAPDH gene. The three probe sets for the GAPDH gene were ranked by their CVs of expres-sion at the upper 0.2th, 0.3th, and 1.7th percentile, respectively, higher than all other well known housekeeping genes, including the actin and the 18S gene. Therefore, the GAPDH gene has the most stable expression as compared with other well known endogenous control genes and is also ranked among the most consistently expressed genes on the array. DISCUSSION It has already been suggested that functional changes of osteoclasts (1), osteoblasts (2), osteocytes (3,4), and mesenchymal stem cells (5) underlie the pathogenesis of osteoporosis. Although circulating monocytes are an important cell type for bone metabolism (6 -13), their role in the pathophysiology of osteoporosis has not been systematically explored in vivo. In this study, we scanned expression of ϳ14,500 genes of circulating monocytes from human subjects using the microarray approach and identified, among the others, three important genes related to BMD variation. The potential functions of the three genes allowed us to propose a model (Fig. 2) suggesting a potential and novel monocyte-mediated pathophysiological mechanism for osteoporosis.
Chemotaxis is a major mechanism whereby circulating monocytes migrate into the bone microenvironment (15). Given that CCR3 ligand, RANTES (Chemokine, CC motif, ligand 5), was found in both osteoblasts and osteoclasts (38 -40), CCR3 may play an important role in the recruitment of monocytes into bone. In addition, in dendritic cells that are very similar to monocytes in terms of phenotype and origin, the CCR3 gene mediates in vitro chemotactic response (41), a process essentially the same as that by which monocytes migrate into the bone microenvironment. The CCR3 gene was also found to be up-regulated during the receptor activator of NF-B ligandinduced osteoclastogenesis (28), which suggests its potential role in promoting monocyte differentiation into osteoclasts. Thus, the up-regulation of CCR3 in circulating monocytes of low BMD women as discovered in the present study may facilitate both the access of monocytes to the bone microenvironment and their subsequent differentiation into osteoclasts. These two concurrent processes may lead to increased bone resorption.
HDC is the only enzyme for synthesis of histamine in the human body (42). With the HDC gene up-regulated, monocyte histamine production may increase in subjects with lower BMD. A high level of circulating histamine was associated with osteoporosis and abnormal bone histomorphometric indices of resorption (43)(44)(45). Histamine receptor antagonists were found to have treatment effects on osteoporosis (46,47). HDC knockout mice (suffering from histamine deficiency) demonstrated increased bone formation and decreased bone loss induced by ovariectomy (42). Histamine may also have effects on the kidney by decreasing synthesis of calcitriol and leading to consequent secondary hyperparathyroidism (42). Locally, histamine can stimulate osteoclast activity (33) and increase osteoclast number (33), production of granulocyte-macrophage colonystimulating factor and interleukin-6 by mononuclear cells (29,30), and expression of tumor necrosis factor-␣ by monocytes (31), which are important cytokines promoting osteoclastogenesis (48,49). Histamine may also stimulate monocytes to produce MCP-1 (monocyte chemoattractant protein-1) and CCR2 (chemokine receptor-2) (25), two important factors mediating chemotaxis of monocytes (50 -52), thereby increasing monocyte recruitment into bone.
Higher expression of the GCR gene in subjects of lower BMD may increase sensitivity of their monocytes to glucocorticoids. Glucocorticoid treatment is known to promote secondary osteoporosis (53). Glucocorticoids may serve to prime monocytes/ macrophages toward differentiation into osteoclasts (26,27). A polymorphism of the GCR gene has been associated with a higher sensitivity to glucocorticoids and lower BMD (54).
As illustrated in Fig. 2, the three genes and/or their products may have similar functions in bone remodeling. For example, both CCR3 and HDC may increase bone monocyte recruitment, whereas both HDC and GCR may promote osteoclastogenesis. In such a way, the genes may act synergistically, thereby resulting in a higher rate of monocyte recruitment into bone, an increased number of monocytes committed to differentiation into osteoclasts, and, subsequently, lower BMD.
It has already been suggested that there is no increase in the number of circulating osteoclast precursors in osteoporotic versus normal subjects (13). Thus, an increased accessibility of these precursors (i.e. monocytes) to the bone tissue, as suggested in our study, may represent an important and novel mechanism for osteoporosis. This hypothesis is supported by studies on chemokines and chemokine receptors, the factors responsible for monocyte recruitment, which have shown their important roles in BMD regulation and bone remodeling (52,55).
Our study for the first time suggested monocyte production of histamine as a pivotal factor underlying osteoporosis. As shown in Fig. 2, histamine affects bone metabolism at multiple levels, locally in osteoclasts (33), mononuclear cells (29,30), and monocytes (25,31) and systemically in the kidney (42). Although the primary cell type for histamine production is mast cells, they are normally not found in the circulation and thus are probably clinically unimportant for osteoporosis, except in hyperplastic conditions, such as mastocytosis (43)(44)(45). In monocytes, the expression of the HDC gene and production of histamine has already been observed in several studies (56 -58). Histamine produced by monocytes was also implicated in the pathogenesis of chronic inflammatory diseases, such as atherosclerosis (59). Recently, in a study showing the effects of histamine on trabecular bone loss in ovariectomized rats, the phagocyte/monocyte lineage was suggested as a cellular source of histamine causing the bone loss in the animals (46). Our findings, together with the above, suggested histamine production from circulating monocytes as an important factor modulating bone loss and proposed histamine receptor antagonists as a potential drug for treatment of osteoporosis.
According to our real-time RT-PCR results, there was no significant difference in the GCR gene expression between the low and the high BMD subjects in the postmenopausal females, although such difference was detected in the premenopausal females. This may suggest an interaction between menopause and BMD for the GCR gene expression. Indeed, by fitting the real-time RT-PCR data of the gene into a general linear model, a significant interaction (p ϭ 0.016) between BMD and menopause was detected (data not shown). Consistent with this, the GCR gene was down-regulated due to menopause (p ϭ 0.03) only in low BMD females.
In our study, only three out of six selected genes suggested by microarrays were confirmed with real-time RT-PCR. Although this may reflect a high false-positive discovery rate associated with microarray technology, statistically it may also be attributed to the multiple-testing procedure in the process of microarray gene identification. Generally, due to multiple-testing procedures in a microarray experiment, the power required to suggest a gene differentially expressed by microarrays is considerably higher than that required to confirm the differential expression of that gene in a subsequent experiment. If, for example, for a list of 10 genes, the power to detect each one of them individually is 10%, then the power required to suggest, in a microarray experiment, any one of them (i.e. at least one of them) is increased to about 65% (i.e. 1-(1-10%) 10 ). However, for such a gene suggested by microarrays, the power decreases to essentially 10% again when it comes to confirm it in a subsequent focused replication experiment, either in an independent sample, or as in our study, in the same sample but on a different platform (e.g. RT-PCR, Northern blotting, etc.). This drastic drop of power makes the confirmation of the differential expression of a gene suggested by microarrays difficult and may contribute to a moderate confirmation rate in our study.
Although our data suggested a significant correlation between the expression levels of the three genes and BMD Zscores, the difference in expression levels of monocyte HDC, GCR, and CCR3 gene among women, whose BMD Z-scores span from normal to severely osteoporotic, might be moderate. In normal conditions, circulating monocytes are mature cells of finite differentiation stage, functionally stable and homogeneous in phenotype until entering peripheral tissues (60). Therefore, we expected that, unless pathologic conditions occur, such as immune diseases and inflammation, their expression change should be generally small in magnitude. Although in our study the low BMD subjects were not severely osteoporotic, most of the high BMD subjects had very high Z-scores, which made the difference in Z-score between the two groups (high versus low BMD) as large as 3.90. Even with the two groups so extremely segregated in terms of BMD, the gene expression fold change for the three genes between the two groups was only from 2 to 4.
Minor functional differences of monocytes, as detected by microarrays, might magnify exponentially in the bone microenvironment, where many specific bone-related factors may interact or synergize. For example, in the bone tissue, an increase in monocyte histamine production may locally induce MCP-1 and CCR2 production (25), leading to monocyte accumulation. As more monocytes accumulate, local histamine levels may increase, which in turn, may induce production of more MCP-1 and CCR2. With such cycles, a tiny increase in monocyte histamine production may be enhanced locally in the bone to physiologically significant levels necessary to initiate a higher rate of monocyte recruitment and osteoclastogenesis and, consequently, bone resorption as illustrated in Fig. 2.
The present study reports data about differential expression of genes at the mRNA level only. In future functional genomics studies of circulating monocytes, the results may be validated also at the protein level. Another interesting direction is to genetically manipulate the expression of the three genes in freshly isolated monocytes using such approaches as RNA interference and viral overexpression and look at the ability of the cells to attach to bone, undergo osteoclastogenesis, or resorb bone ex vivo. This may further test the functions of the genes as suggested by the present study.
The present study for the first time provides evidence for the functional difference of circulating monocytes between low and high BMD women. The model proposed here (Fig. 2) suggests how the three genes, CCR3, HDC, and GCR, may potentially contribute to increased bone resorption, making them a new group of functionally related skeletal genes in monocytic cell lineages. To the best of our knowledge, this is the first in vivo microarray study in humans to suggest another potential pathophysiological mechanism underlying osteoporosis.