Background: The study aimed to utilize machine learning (ML) approaches and genomic data to develop the prediction model for bone mineral density (BMD), and to identify the best modeling approach for BMD prediction. Method: The genomic and phenotypic data of Osteoporotic Fractures in Men Study (n=5,130), was analyzed. Genetic risk score (GRS) was calculated from 1,103 associated SNPs for each participant after a comprehensive genotype imputation. Data were normalized and divided into a training set (80%) and a validation set (20%) for analysis. Random forest, gradient boosting, neural network, and linear regression were used to develop prediction models for BMD separately. The 10-fold cross-validation was used for hyperparameter optimization. Mean square error and mean absolute error were used to assess model performance. Results: When using GRS and phenotypic covariates as the predictors, the performance of all ML models and linear regression in BMD prediction is similar. However, when replacing GRS with the 1,103 individual SNPs in the model, ML models performed significantly better than linear regression, and the gradient boosting model performed the best. Conclusion: Our study suggested that ML models, especially gradient boosting, can improve BMD prediction in genomic data.