Implemented here are procedures for fitting hierarchical generalized linear models (HGLM). It can be used for linear mixed models and generalized linear mixed models with random effects for a variety of links and a variety of distributions for both the outcomes and the random effects. Fixed effects can also be fitted in the dispersion part of the mean model. As statistical models, HGLMs were initially developed by Lee and Nelder (1996) . We provide an implementation (Ronnegard, Alam and Shen 2010) following Lee, Nelder and Pawitan (2006) with algorithms extended for spatial modeling (Alam, Ronnegard and Shen 2015) .