ABSTRACT Genome-wide association studies (GWASs) for osteoporotic traits have identified over 1000 associations; however, their impact has been limited by the difficulties of causal gene identification and a strict focus on bone mineral density (BMD). Here, we used Diversity Outbred (DO) mice to directly address these limitations by performing the first systems genetics analysis of 55 complex skeletal phenotypes. We applied a network approach to cortical bone RNA-seq data to discover 72 genes likely to be causal for human BMD GWAS associations, including the novel genes SERTAD4 and GLT8D2 . We also performed GWAS in the DO for a wide-range of bone traits and identified Qsox1 as a novel gene influencing cortical bone accrual and bone strength. Our results provide a new perspective on the genetics of osteoporosis and highlight the ability of the mouse to inform human genetics.

Abstract Genome-wide association studies (GWASs) for bone mineral density (BMD) have identified over 1,100 associations to date. However, identifying causal genes implicated by such studies has been challenging. Recent advances in the development of transcriptome reference datasets and computational approaches such as transcriptome-wide association studies (TWASs) and expression quantitative trait loci (eQTL) colocalization have proven to be informative in identifying putatively causal genes underlying GWAS associations. Here, we used TWAS/eQTL colocalization in conjunction with transcriptomic data from the Genotype-Tissue Expression (GTEx) project to identify potentially causal genes for the largest BMD GWAS performed to date. Using this approach, we identified 512 genes as significant (Bonferroni <= 0.05) using both TWAS and eQTL colocalization. This set of genes was enriched for regulators of BMD and members of bone relevant biological processes. To investigate the significance of our findings, we selected PPP6R3, the gene with the strongest support from our analysis which was not previously implicated in the regulation of BMD, for further investigation. We observed that Ppp6r3 deletion in mice decreased BMD. In this work, we provide an updated resource of putatively causal BMD genes and demonstrate that PPP6R3 is a putatively causal BMD GWAS gene. These data increase our understanding of the genetics of BMD and provide further evidence for the utility of combined TWAS/colocalization approaches in untangling the genetics of complex traits.

ABSTRACT Bone mineral density (BMD) is a strong predictor of osteoporotic fracture. It is also one of the most heritable disease-associated quantitative traits. As a result, there has been considerable effort focused on dissecting its genetic basis. Here, we performed a genome-wide association study (GWAS) in a panel of inbred strains to identify associations influencing BMD. This analysis identified a significant (P=3.1 x 10 −12 ) BMD locus on Chromosome 3@52.5 Mbp that replicated in two seperate inbred strain panels and overlapped a BMD quantitative trait locus (QTL) previously identified in a F2 intercross. The association mapped to a 300 Kbp region containing four genes; Gm2447, Gm20750, Cog6 , and Lhfp . Further analysis found that Lipoma HMGIC Fusion Partner ( Lhfp ) was highly expressed in bone and osteoblasts and its expression was regulated by local expression QTL (eQTL) in multiple tissues. A co-expression network analysis revealed that Lhfp was strongly connected to genes involved in osteoblast differentiation. To directly evaluate its role in bone, Lhfp deficient mice ( Lhfp −/− ) were created using CRISPR/Cas9. Consistent with genetic and network predictions, bone marrow stromal cells (BMSCs) from Lhfp −/− displayed increased osteogenic differentiation. Lfhp −/− mice also had elevated BMD due to increased cortical bone mass. In conclusion, we used GWAS and systems genetics in mice to identify Lhfp as a regulator of osteoblast activity and bone mass.

Abstract Genome-wide association studies (GWASs) have revolutionized our understanding of the genetics of complex diseases, such as osteoporosis; however, the challenge has been converting associations to causal genes. Studies have demonstrated the utility of transcriptomics data in linking disease-associated variants to genes; though for osteoporosis, few population transcriptomics datasets have been generated on bone or bone cells, and an even smaller number have profiled individual cell-types. To begin to evaluate approaches to address this challenge, we profiled the transcriptomes of bone marrow-derived stromal cells (BMSCs) cultured under osteogenic conditions, a popular model of osteoblast differentiation and activity, from five Diversity Outbred (DO) mice using single-cell RNA-seq (scRNA-seq). The goal of the study was to determine if BMSCs could serve as a model for the generation of cell-type specific transcriptomic profiles of mesenchymal lineage cells derived from large populations of mice to inform genetic studies. We demonstrate that dissociation of BMSCs from a heavily mineralized matrix had little effect on viability or their transcriptomic signatures. Furthermore, we show that BMSCs cultured under osteogenic conditions are diverse and consist of cells with characteristics of mesenchymal progenitors, marrow adipogenic lineage precursors (MALPs), osteoblasts, osteocyte-like cells, and immune cells. Importantly, all cells were nearly identical from a transcriptomic perspective to cells isolated directly from bone. We also demonstrated the ability to multiplex single cells and subsequently assign cells to their “mouse-of-origin” using demultiplexing approaches based on genotypes inferred from coding SNPs. We employed scRNA-seq analytical tools to confirm the biological identity of profiled cell-types. SCENIC was used to reconstruct gene regulatory networks (GRNs) and we showed that identified cell-types show GRNs expected of osteogenic and pre-adipogenic lineage cells. Further, CELLECT analysis showed that osteoblasts, osteocyte-like cells, and MALPs captured a significant component of BMD heritability. Together, these data suggest that BMSCs cultured under osteogenic conditions coupled with scRNA-seq can be used as a scalable and biologically informative model to generate cell-type specific transcriptomic profiles of mesenchymal lineage cells in large mouse, and potentially human, populations.

Abstract Osteoporosis, characterized by low bone mineral density (BMD), is the most common complex disease affecting bone and constitutes a major societal health problem. Genome-wide association studies (GWASs) have identified over 1100 associations influencing BMD. It has been shown that perturbations to long non-coding RNAs (lncRNAs) influence BMD and the activities of bone cells; however, the extent to which lncRNAs are involved in the genetic regulation of BMD is unknown. Here, we combined the analysis of allelic imbalance (AI) in human acetabular bone fragments with a transcriptome-wide association study (TWAS) and expression quantitative trait loci (eQTL) colocalization analysis using data from the Genotype-Tissue Expression (GTEx) project to identify lncRNAs potentially responsible for GWAS associations. We identified 27 lncRNAs in bone that are located in proximity to a BMD GWAS association and harbor SNPs demonstrating AI. Using GTEx data we identified an additional 31 lncRNAs whose expression was associated (FDR correction<0.05) with BMD through TWAS and had a colocalizing eQTL (regional colocalization probability (RCP)>0.1). The 58 lncRNAs are located in 43 BMD associations. To further support a causal role for the identified lncRNAs, we show that 23 of the 58 lncRNAs are differentially expressed as a function of osteoblast differentiation. Our approach identifies lncRNAs that are potentially responsible for BMD GWAS associations and suggest that lncRNAs play a role in the genetics of osteoporosis.