Identifying genetic variants associated with disease phenotypes is one of the major approaches to deepen our understanding of the genomic and molecular basis of diseases. Most notably, major effort has been put into finding variants associated with changes in gene expression, termed expression quantitative trait loci (eQTL). Recent advancements in single-cell sequencing technology have also enabled the possibility of mapping eQTLs with individual cellular resolution with single-cell RNA-seq. However, while assay technologies for other biological modalities exist, single-cell QTL studies for molecular phenotypes other than gene expression are rare. In particular, there is no bespoke computational workflow tailored for mapping single-cell chromatin accessibility QTL (sc-caQTL) with single-cell ATAC-seq data. This is partly due to the extreme sparsity and the difference in underlying biology of the data. In this study, we explored the statistical properties of scATAC-seq data and proposed a statistical model for sc-caQTL mapping. We proposed a workflow for sc-caQTL mapping and compared it against existing workflows designed for bulk caQTL mapping. Our results show that using a bespoke single-cell ATAC-seq count model can provide more power in identifying sc-caQTLs than traditional pseudobulk approaches.