Every Breeding Value Is a Regression: How Covariance Connects Falconer to BLUP

No one has ever seen a breeding value. Not in a field, not in a catalogue, not in a selection list. Every number we rank our lines on is a prediction of something that cannot be measured.

Quantitative geneticists mark this with a small hat. The true breeding value is u; what our software returns is û. Because u is a random variable and not a fixed parameter, Henderson insisted that we predict it rather than estimate it. "Estimated breeding value" is a name we keep out of habit.

The hat hides a simple fact that this edition sets out to show: every breeding value is a regression. Falconer wrote it for one plant with one record. Henderson wrote it for a whole programme with millions. The slope is the same in both, and covariance sets it. In Edition 01 we followed covariance from shared alleles to the A and G matrices. Now we watch it turn into a breeding value, and then clear up four words that are often mixed up: BV, EBV, GEBV and BLUP.

1. Act I · Figure 1Heritability is a slopeFalconer's BV, h² as regression, Â = h²(P − μ)
2. Act II · Figure 2Borrowing strength from relativesOff-diagonals become weights; BLUP reduces to Falconer
3. Act III · Figure 3u, û and the four namesA or G, same engine; four names, one model

What Falconer meant, and why heritability is a slope

A breeding value is what a plant passes on, not what it is. A parent does not hand its genotype to its offspring; it hands over a random half of its alleles. What travels between generations is therefore the average effect of each allele: how much, on average, the offspring that receive it differ from those that do not (Fisher, 1918). Falconer defined the breeding value as the sum of these average effects over all loci (Fig. 1a). The same quantity can be measured in the field: twice the deviation of a parent's progeny mean from the population mean, under random mating.

Two things follow from this definition. A breeding value belongs to a population, because average effects depend on allele frequencies. And it can never be observed directly, because we see only the phenotype, P = μ + A + E, where the breeding value A is mixed with everything else.

So we ask a regression question: given a plant's phenotype, what is our best guess of its breeding value? Plot the breeding values of 300 plants against their phenotypes and fit a line (Fig. 1b). The slope of that line is

This gives Falconer's prediction, the first breeding value most of us met:

And the adjusted mean is a BLUE. In a real trial, P is not a single plot. A line is grown in several plots and blocks, and its adjusted mean is what the mixed model returns when lines are fitted as fixed effects: a best linear unbiased estimate, a BLUE. Even the mean μ̂ is a BLUE. Falconer's formula then reads:

Shrinkage is not a correction; it is the regression. A plant one tonne above the mean is predicted at +0.8 when h² is 0.8 and at only +0.2 when h² is 0.2 (Fig. 1c). Part of every phenotypic advantage is environment, and the slope removes exactly that part, on average. The accuracy of the prediction, the correlation between  and A, is h (Fig. 1d): 0.63 for our plant.

Heritability answers one question: how far should I believe a phenotype? It is the slope that turns what we see into what we expect to be passed on.

Borrowing strength from relatives: it's all in the off-diagonals

Plant T is not alone. It has a full sib, a half sib, a cousin, and a neighbour in the trial that is no relation at all (Fig. 2a). Each has a yield record. Which of them should change our view of T?

Regression gives the answer again. We now regress T's breeding value on all five records at once. The slope becomes a vector of weights, and it is built from the same two ingredients as before: how each record co-varies with T's breeding value, and how the records co-vary with one another.

Work it through (Fig. 2b). T's own record gets a weight of 0.370, a little below h², because part of what it says is already said by the relatives. The full sib gets 0.120, the half sib 0.051 and the cousin 0.026. The unrelated plant, the highest yielder in the group at 6.8 t/ha, gets a weight of exactly zero. It is impressive, and it tells us nothing about T. The prediction rises to +0.58 t/ha and its accuracy from 0.63 to 0.67.

Now suppose T was never planted. The same formula, using only the relatives, still returns a prediction: +0.33 t/ha, with accuracy 0.35. This is the logic of the progeny test and, as Edition 01 showed, of genomic prediction. With enough relatives, a parent can be judged more accurately from its offspring than from its own plot (Fig. 2c), above all when heritability is low.

Every breeding value is a regression. Write the same idea for every plant and every record in a programme, add fixed effects for trials, blocks and years, and you have Henderson's best linear unbiased predictor:

Now take one plant with one record and no relatives (Fig. 2d). The covariance between its breeding value and its record is σ²_A, the variance of the record is σ²_P, and the fixed effect is μ. The matrix formula collapses to σ²_A/σ²_P × (y − μ) = h²(y − μ). Falconer's prediction is not an old rule of thumb that BLUP replaced. It is BLUP, with one record.

Falconer and Henderson wrote the same regression. Falconer had one record and a slope called h². Henderson had all the records and a slope built from covariance.

u, û and the covariance between them

If û is a regression prediction, it should behave like one. Plot the true breeding values of our 300 plants against their predictions (Fig. 3a). Three properties of BLUP appear.

The slope of u on û is 1. When the prediction says +0.4, the true value is +0.4 on average. This is what unbiased means for a predictor, and it is the reason breeders can compare EBVs at face value. The predictions are shrunk. Var(û) is 0.16 while Var(u) is 0.37, and Cov(u, û) equals Var(û). Reliability is the ratio, r² = Var(û)/σ²_A, here 0.43. A prediction near zero may be a mediocre plant, or simply a plant we know little about.

A or G, same engine

The difference between an EBV and a GEBV is not the method. It is the covariance matrix. Use the pedigree (A) and you get PBLUP and EBVs; use markers (G) and you get GBLUP and GEBVs; combine the two (H) and you get single-step GBLUP. The equations are untouched.

The difference shows best inside a family. A cross produced 200 sibs; 150 were grown in the trial and 50 were not (Fig. 3b). Every untested sib has the same pedigree, so every one receives the same EBV, +0.68 t/ha. A cannot see which sib inherited the better half of each parent. G can. The GEBVs of the 50 untested sibs correlate 0.63 with their true breeding values, and the ten sibs with the highest GEBV are, on average, 0.30 t/ha better than the family mean. With A alone, picking the best sib would be a coin toss.

Four names, one model

Now the four words can be put in their places (Fig. 3c). Read the figure as a metro line. Records go in, a mixed model is fitted, and two kinds of number come out: BLUEs for the fixed effects and BLUPs for the random ones. When the random effect is the additive genetic value, with covariance A or G, its BLUP has a name of its own: the EBV or the GEBV. The true breeding value, u, sits at the end of the line. It is never reached, only predicted.

So when you fit GBLUP, what comes out for the breeding values is a BLUP, and calling it the BLUP of the breeding value is exactly right: GBLUP output = BLUP of BV = GEBV. Two cautions keep the words honest. First, the random effect decides what the BLUP is a BLUP of. Fit lines with I instead of A or G, as in Route 1 of Edition 01, and the BLUP is of a genotypic value, not a breeding value; with a dominance or Gaussian kernel it is a non-additive or total genetic value. Second, the arrow runs one way. Every BLUP of a breeding value is an EBV or GEBV, but not every GEBV is a BLUP. A classical selection index with known means is the best linear predictor (BLP). Bayesian methods with non-normal priors, such as BayesA, BayesB, BayesCπ and the Bayesian LASSO, return posterior means that are not linear in the data. RKHS and machine-learning models take other routes. All of them give GEBVs; none of them gives a BLUP.

From average effects to BLUP, in one animated figure

Press play and watch a breeding value being built from the average effects of five genes, then hidden behind a phenotype. Watch the regression line find h² in a cloud of 300 plants, and use the h² switch to see the shrinkage change. Then see relatives turn into weights for plant T, and the same numbers come out of the mixed model as BLUPs. Switch A to I in the last step and Falconer's h²(y − μ) comes back. The figure ends on the metro map of Figure 3c, drawn station by station.

Open the animation full screen · Five steps, about 80 seconds. Numbers match Figures 1 and 2.

From Falconer to BLUP by hand, in base R

No packages. We predict plant T from its own record, then add the four relatives of Figure 2, take them away again, and finally solve Henderson's equations with μ estimated at the same time.

# The Crop Science Brief, Edition 02: Every Breeding Value Is a Regression
# Base R only. Yield in t/ha. mu = 5, h2 = 0.4, phenotypic variance = 1.

mu   <- 5
h2   <- 0.4
s2A  <- 0.4                 # additive genetic variance
s2e  <- 0.6                 # residual variance  (s2A + s2e = 1)

# ---- STEP 1. Falconer: predict a breeding value from one plant's own record ----
P    <- 6.0                 # plant T yielded 6.0 t/ha
Ahat <- h2 * (P - mu)       # regression of A on P: slope = Cov(A,P)/Var(P) = h2
Ahat                        # 0.4
sqrt(h2)                    # accuracy of this prediction = h = 0.63

# ---- STEP 2. Add four relatives: the covariance does the work ----
# Plants: T (target), FS (full sib), HS (half sib), C (cousin), U (unrelated)
A <- matrix(c(1,     0.5,  0.25, 0.125, 0,
              0.5,   1,    0.25, 0.125, 0,
              0.25,  0.25, 1,    0,     0,
              0.125, 0.125,0,    1,     0,
              0,     0,    0,    0,     1), nrow = 5,
            dimnames = list(c("T","FS","HS","C","U"), c("T","FS","HS","C","U")))
y <- c(T = 6.0, FS = 6.4, HS = 5.6, C = 5.3, U = 6.8)

V    <- A * s2A + diag(5) * s2e     # Var(y): covariance of the five records
covu <- A["T", ] * s2A              # Cov(u_T, y): each record with T's breeding value
w    <- covu %*% solve(V)           # regression weights = Cov(u,y') Var(y)^-1
round(w, 3)                         # T 0.370, FS 0.120, HS 0.051, C 0.026, U 0.000
uhat <- w %*% (y - mu)              # BLUP formula with mu known
round(uhat, 2)                      # 0.58
round(sqrt(covu %*% solve(V) %*% covu / s2A), 2)   # accuracy rises to 0.67

# ---- STEP 3. Take the relatives away: A becomes I, and we are back to Falconer ----
I5   <- diag(5)
V0   <- I5 * s2A + I5 * s2e
w0   <- (I5[1, ] * s2A) %*% solve(V0)
round(w0, 3)                        # 0.4 on T, 0 on everyone else: the weight IS h2
w0 %*% (y - mu)                     # 0.4 = h2 * (P - mu)

# ---- STEP 4. T never planted: predict it from relatives only ----
r    <- 2:5                         # FS, HS, C, U
wr   <- covu[r] %*% solve(V[r, r])
round(wr, 3)                        # FS 0.190, HS 0.081, C 0.041, U 0.000
round(wr %*% (y[r] - mu), 2)        # 0.33, with accuracy:
round(sqrt(covu[r] %*% solve(V[r, r]) %*% covu[r] / s2A), 2)   # 0.35

# ---- STEP 5. Henderson's mixed-model equations: estimate mu at the same time ----
X   <- matrix(1, 5, 1)              # fixed effect: overall mean
Z   <- diag(5)                      # each record belongs to one plant
lam <- s2e / s2A
LHS <- rbind(cbind(t(X) %*% X, t(X) %*% Z),
             cbind(t(Z) %*% X, t(Z) %*% Z + solve(A) * lam))
RHS <- c(t(X) %*% y, t(Z) %*% y)
sol <- solve(LHS, RHS)
round(sol, 3)   # first value = BLUE of mu, the rest = BLUPs of breeding values
# With mu known (steps 2-4) the predictor is the BLP, the classical selection index.
# Estimating mu jointly makes it the BLUP. Here mu comes from only five plants (6.02),
# so every prediction moves; with thousands of records the two agree closely.

Three works behind this edition

1. Falconer DS, Mackay TFC (1996). Introduction to Quantitative Genetics, 4th edn. Longman.Average effects, breeding value and heritability as a regression (chapters 7 and 10): Figure 1.
2. Henderson CR (1975). Best linear unbiased estimation and prediction under a selection model. Biometrics 31: 423–447.BLUE for fixed effects, BLUP for random ones, from one set of equations: Figure 2d.
3. Robinson GK (1991). That BLUP is a good thing: the estimation of random effects. Statistical Science 6: 15–32.The clearest account of what BLUP is, and is not, across statistics and breeding: Figure 3c.

The Crop Science Brief now comes every two weeks. Next time a list of EBVs or GEBVs crosses your desk, ask three questions: what covariance was used, how reliable is each value, and which model produced it?

Your turn: in your programme, which word is used most loosely: EBV, GEBV or BLUP? Tell us in the comments.

If this was useful, share it with one colleague and subscribe to The Crop Science Brief.

#CropScience #PlantBreeding #QuantitativeGenetics #BLUP #GenomicSelection #PlantNura