Computation of restricted maximum likelihood estimates of variance components

The method preferred by animal breeders for the estimation of variance components is restricted maximum likelihood (REML). Various iterative algorithms have been proposed for computing REML estimates. Five different computational strategies for implementing such an algorithm were compared in terms of flops (floating-point operations). These strategies were based respectively on the LDL\u27 decomposition, the W transformation, the SWEEP method, tridiagonalization and diagonalization of the coefficient matrix of the mixed-model equations;The computational requirements of the orthogonal transformations employed in tridiagonalization and diagonalization were found to be rather extensive. However, these transformations are performed prior to the initiation of the iterative estimation process and need not be repeated during the remainder of the process. Subsequent to either diagonalization or tridiagonalization, the flops required per iteration are very minimal. Thus, for most applications of mixed-effects linear models with a single set of random effects, the use of an orthogonal transformation prior to the initiation of the iterative process is recommended. For most animal breeding applications, tridiagonalization will generally be more efficient than diagonalization;In most animal breeding applications, the coefficient matrix of the mixed-model equations is extremely sparse and of very large order. The use of sparse-matrix techniques for the numerical evaluation of the log-likelihood function and its first- and second-order partial derivatives was investigated in the case of the simple sire and animal models. Instead of applying these techniques directly to the coefficient matrix of the mixed-model equations to obtain the Cholesky factor, they were used to obtain the Cholesky factor indirectly by carrying out a QR decomposition of an augmented model matrix;The feasibility of the computational method for the simple sire model was investigated by carrying out the most computationally intensive part of this method (which is the part consisting of the QR decomposition) for an animal breeding data set comprising 180,994 records and 1,264 sires. The total CPU time required for this part (using an NAS AS/9160 computer) was approximately 75,000 seconds