This article presents a new investigation to the linear mixed model $\by = \bX \bbe + \bZ\bga + \bve$ with fixed effect $\bX\bbe$ and random effect $\bZ\bga$ under a general assumption via some novel algebraic tools in matrix theory, and reveals a variety of deep and profound properties hidden behind the linear mixed model. We first derive exact formulas for calculating the best linear unbiased predictor (BLUP) of a general vector $\bph = \bF\bbe + \bG\bga + \bH\bve$ of all unknown parameters in the model by solving a constrained quadratic matrix-valued function optimization problem in the L\”owner partial ordering. We then consider some special cases of the BLUP for different choices of $\bF$, $\bG$, and $\bH$ in $\bph$, and establish some fundamental decomposition equalities for the observed random vector $\by$ and its covariance matrix.