Comparison of PQL and Laplace 6 estimates of hierarchical linear models when comparing groups of small incident rates in cluster randomised trials

The variances of the random components in hierarchical generalised linear models (HGLMs) with binary outcomes have been reported to have a considerable downward bias when estimated with the commonly used penalised quasilikelihood (PQL) technique. The more recently proposed Laplace 6 approximation promises to reduce this bias. This study compares the performance of these two techniques when estimating the parameters of a particular HGLM. This comparison is performed via Monte Carlo simulations in which the difference between two groups of proportions, modelled after those appearing in many epidemiological cluster randomised interventions, are tested using this model. The Laplace 6 approximation does reduce the bias mentioned above, but at the price of a higher mean square error. The results of this study suggest that the optimal solution involves using a combination of these two techniques. This combination is illustrated by analysing a data set from a real cluster randomised intervention.