Objectives: In this paper, we introduce the use of generalized linear mixed models (GLMM) as a better alternative to traditional statistical methods for studying factors associated to the prevalence of degenerative joint disease (DJD) in bioarchaeological contexts.Materials and Methods: DJD prevalence was assessed for the appendicular joints and the spine of a Spanish population dated from the 15th to the 18th century. Data were analyzed using contingency tables, logistic regression models, and logistic GLMM.Results: In general, results from GLMMs find agreement in other methods. However, by being able to analyze the data at the level of individual bones instead of aggregated joints or limbs, GLMMs are capable of revealing associations that are not evident in other frameworks.Discussion: Currently widely available in statistical analysis software, GLMMs can accommodate a wide array of data distributions, account for hierarchical correlations, and return estimates of DJD prevalence within individuals and skeletal locations that are unbiased by the effect of covariates. This gives clear advantages for the analysis of bioarchaeological datasets which can lead to more robust and comparable analyses across populations.