maximum Likelihood (ML), Restricted Maximum Likelihood (REML) and Bayesian methods are often preferred over other methods for estimating variance components in animal breeding. Iterative computing stategies are required for obtaining estimates with unbalanced data and models with at least two variance components. If iteration converges and the converged value is within the parameter space, it is commonly accepted as the “estimate”. However, if the likelihood or posterior density function is not unimodal, the converged value may correspond to a local but not global maximum, and not be the desired estimate. Bimodal likelihood functions have been found in ML. Analytical proofs and numerical results are presented for two–variance–component models showing that the likelihood or posterior density function is always unimodal for REML and a Bayesian method but sometimes bimodal over the permissible parameter space for ML and another Bayesian method.