Abstract Background The mixed linear model employed for genomic best linear unbiased prediction (GBLUP) includes the breeding value for each animal as a random effect that has a mean of zero and a covariance matrix proportional to the genomic relationship matrix ( $${\mathbf {G}}{gg}$$ G g g ), where the inverse of $${\mathbf {G}}{gg}$$ G g g is required to set up the usual mixed model equations (MME). When only some animals have genomic information, genomic predictions can be obtained by an extension known as single-step GBLUP, where the covariance matrix of breeding values is constructed by combining the pedigree-based additive relationship matrix with $${\mathbf {G}}{gg}$$ G g g . The inverse of the combined relationship matrix can be obtained efficiently, provided $${\mathbf {G}}{gg}$$ G g g can be inverted. In some livestock species, however, the number $$N_{g}$$ N g of animals with genomic information exceeds the number of marker covariates used to compute $${\mathbf {G}}{gg}$$ G g g , and this results in a singular $${\mathbf {G}}{gg}$$ G g g . For such a case, an efficient and exact method to obtain GBLUP and single-step GBLUP is presented here. Results Exact methods are already available to obtain GBLUP when $${\mathbf {G}}{gg}$$ G g g is singular, but these require working with large dense matrices. Another approach is to modify $${\mathbf {G}}{gg}$$ G g g to make it nonsingular by adding a small value to all its diagonals or regressing it towards the pedigree-based relationship matrix. This, however, results in the inverse of $${\mathbf {G}}{gg}$$ G g g being dense and difficult to compute as $$N{g}$$ N g grows. The approach presented here recognizes that the number r of linearly independent genomic breeding values cannot exceed the number of marker covariates, and the mixed linear model used here for genomic prediction only fits these r linearly independent breeding values as random effects. Conclusions The exact method presented here was compared to Apy-GBLUP and to Apy single-step GBLUP, both of which are approximate methods that use a modified $${\mathbf {G}}_{gg}$$ G g g that has a sparse inverse which can be computed efficiently. In a small numerical example, predictions from the exact approach and Apy were almost identical, but the MME from Apy had a condition number about 1000 times larger than that from the exact approach, indicating ill-conditioning of the MME from Apy. The practical application of exact SSGBLUP is not more difficult than implementation of Apy.