Expression quantitative trait loci (eQTL) mapping is a powerful tool for investigating the impact of genetic variants on gene expression. In eQTL mapping, the data is typically high-dimensional whereas the sample size is limited. Sparse learning models such as LASSO have shown their strengths to select associated features in such high dimensional data. However, plain LASSO performs poorly when dealing with extremely high dimensional datasets. In this study, we introduce two novel scalable methods named SLASSO and PLASSO which allow efficient learning for datasets of ultra-high dimension based on "divide and conquer". We provided a multi-round procedure to address the sample size limitation for real applications. We performed extensive simulations on synthetic data to validate our methods and evaluate their performance. Comparing to similar methods, our methods showed similar precision and recall, but outperformed them on scalability, especially for increased data dimensions. We further demonstrated the application of our methods by applying them to a real human genomics data set for eQTL mapping. Our methods are not limited to plain LASSO models, it can be extended to variations of LASSO and many other machine learning models.