Calculation of heterosis and recombination loss

In this post, I demonstrate the calculation of heterosis and recombination loss coefficients using an example pedigree.

Heterosis (hybrid vigor) and recombination loss are critical concepts in animal and plant breeding, particularly in crossbreeding systems.

Heterosis is the favourable gene effect in crossbred offspring, calculated as

$$ h_{AB}^k = p_A^s p_B^d + p_B^s p_A^d, $$

where $h_{AB}^k$ is the heterosis coefficient of progeny k for the A and B breed combination, $p_A^s$ is the breed A proportion in sire s (the sire of k), and $p_B^d$ is the breed B proportion in dam d (the dam of k). $p_B^s$ and $p_A^d$ are self-explanatory.

Recombination loss is the unfavourable gene effect due to the breakdown of the parental epistatic gene complex, calculated as

$$ r_{AB}^k = \left( p_A^s + p_A^d \right) \left( p_B^s + p_B^d \right) - h_{AB}^k, $$

where $r_{AB}^k$ is the recombination loss coefficient of progeny k for the A and B breed combination.

Here, I show how to calculate heterosis and recombination loss given a pedigree and breed compositions. Consider the following pedigree.

ID	SIRE	DAM	A	B	O
1	0	0	1	0	0
2	0	0	0	1	0
3	1	0	0.5	0	0.5
4	1	2	0.5	0.5	0
5	1	2	0.5	0.5	0
6	0	0	0.5	0.5	0
7	3	4	0.5	0.25	0.25
8	5	6	0.5	0.5	0
9	7	0	0.25	0.125	0.625
10	0	8	0.25	0.25	0.5

Two assumptions were made to generate breed proportions (A, B, and O for other/unknown) for animals with unknown parents:

Animals with both parents unknown are considered purebred (including O) or F1 crossbred between two known (excluding O) breeds.
Animals with one unknown parent, receive 0.5 O contribution from the unknown parent.

ped.init <- data.frame(
    ID = 1:10,
    SIRE = c(0,0,1,1,1,0,3,5,7,0),
    DAM = c(0,0,0,2,2,0,4,6,0,8),
    A = c(1,0,0.5,0.5,0.5,0.5,0.5,0.5,0.25,0.25),
    B = c(0,1,0,0.5,0.5,0.5,0.25,0.5,0.125,0.25),
    O = c(0,0,0.5,0,0,0,0.25,0,0.625,0.5)
)
(ped <- ped.init)

   ID SIRE DAM    A     B     O
1   1    0   0 1.00 0.000 0.000
2   2    0   0 0.00 1.000 0.000
3   3    1   0 0.50 0.000 0.500
4   4    1   2 0.50 0.500 0.000
5   5    1   2 0.50 0.500 0.000
6   6    0   0 0.50 0.500 0.000
7   7    3   4 0.50 0.250 0.250
8   8    5   6 0.50 0.500 0.000
9   9    7   0 0.25 0.125 0.625
10 10    0   8 0.25 0.250 0.500

The first step is data preparation: bringing the sire and dam breed proportions ($s_A$, $s_B$, $d_A$, $d_B$) into each animal row. Here I use R, but the same logic can be implemented in other programming languages.

brdcols <- c("A", "B", "O")
stopifnot(all(rowSums(ped[,brdcols]) == 1)) # Check that breed proportions sum to 1
brdcols <- brdcols[brdcols != "O"] # Breed O is not needed
stopifnot(length(brdcols) > 1)
stopifnot(brdcols %in% colnames(ped))

ped <-merge(ped, ped[,c("ID", brdcols)], by.x="SIRE", by.y="ID", all.x=TRUE)
# Rename columns after merging
colends.x <- grepl("\\.x", colnames(ped))
colends.y <- grepl("\\.y", colnames(ped))
colnames(ped)[colends.x] <- substr(colnames(ped)[colends.x], 1, nchar(colnames(ped)[colends.x])-2)
colnames(ped)[colends.y] <- paste0("s", substr(colnames(ped)[colends.y], 1, nchar(colnames(ped)[colends.y])-2))

ped <-merge(ped, ped[,c("ID", brdcols)], by.x="DAM", by.y="ID", all.x=TRUE)
# Rename columns after merging
colends.x <- grepl("\\.x", colnames(ped))
colends.y <- grepl("\\.y", colnames(ped))
colnames(ped)[colends.x] <- substr(colnames(ped)[colends.x], 1, nchar(colnames(ped)[colends.x])-2)
colnames(ped)[colends.y] <- paste0("d", substr(colnames(ped)[colends.y], 1, nchar(colnames(ped)[colends.y])-2))
# Set missing to 0
ped[is.na(ped)] <- 0
# Reorder pedigree rows
ped <- ped[order(ped$ID),]
# Set parents' breed for crossbred founders
ped[ped$SIRE==0 & ped$DAM==0 & ped$A==0.5 & ped$B==0.5, c("sA","dB")] <- 1
ped

   DAM SIRE ID    A     B     O  sA   sB  dA  dB
1    0    0  1 1.00 0.000 0.000 0.0 0.00 0.0 0.0
2    0    0  2 0.00 1.000 0.000 0.0 0.00 0.0 0.0
4    0    1  3 0.50 0.000 0.500 1.0 0.00 0.0 0.0
7    2    1  4 0.50 0.500 0.000 1.0 0.00 0.0 1.0
6    2    1  5 0.50 0.500 0.000 1.0 0.00 0.0 1.0
3    0    0  6 0.50 0.500 0.000 1.0 0.00 0.0 1.0
8    4    3  7 0.50 0.250 0.250 0.5 0.00 0.5 0.5
9    6    5  8 0.50 0.500 0.000 0.5 0.50 0.5 0.5
5    0    7  9 0.25 0.125 0.625 0.5 0.25 0.0 0.0
10   8    0 10 0.25 0.250 0.500 0.0 0.00 0.5 0.5

You may choose to reorder columns. sA and dA for animal 1, and sB and dB for animal 2, could be set to 1. No edit is required, as these animals have no heterosis and no recombination to lose. Animal 6 is an F1 crossbred with both parents unknown. sA and dB for animal 6 are set to 1.

The second step is calculating the heterosis coefficients.

ped$hetAB <- ped$sA * ped$dB + ped$sB * ped$dA
ped

   DAM SIRE ID    A     B     O  sA   sB  dA  dB hetAB
1    0    0  1 1.00 0.000 0.000 0.0 0.00 0.0 0.0  0.00
2    0    0  2 0.00 1.000 0.000 0.0 0.00 0.0 0.0  0.00
4    0    1  3 0.50 0.000 0.500 1.0 0.00 0.0 0.0  0.00
7    2    1  4 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00
6    2    1  5 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00
3    0    0  6 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00
8    4    3  7 0.50 0.250 0.250 0.5 0.00 0.5 0.5  0.25
9    6    5  8 0.50 0.500 0.000 0.5 0.50 0.5 0.5  0.50
5    0    7  9 0.25 0.125 0.625 0.5 0.25 0.0 0.0  0.00
10   8    0 10 0.25 0.250 0.500 0.0 0.00 0.5 0.5  0.00

The third and final step is calculating the recombination loss coefficient.

ped$recombAB <- (ped$sA + ped$dA)*(ped$sB + ped$dB) - ped$hetAB
ped

   DAM SIRE ID    A     B     O  sA   sB  dA  dB hetAB recombAB
1    0    0  1 1.00 0.000 0.000 0.0 0.00 0.0 0.0  0.00    0.000
2    0    0  2 0.00 1.000 0.000 0.0 0.00 0.0 0.0  0.00    0.000
4    0    1  3 0.50 0.000 0.500 1.0 0.00 0.0 0.0  0.00    0.000
7    2    1  4 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
6    2    1  5 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
3    0    0  6 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
8    4    3  7 0.50 0.250 0.250 0.5 0.00 0.5 0.5  0.25    0.250
9    6    5  8 0.50 0.500 0.000 0.5 0.50 0.5 0.5  0.50    0.500
5    0    7  9 0.25 0.125 0.625 0.5 0.25 0.0 0.0  0.00    0.125
10   8    0 10 0.25 0.250 0.500 0.0 0.00 0.5 0.5  0.00    0.250

There may be more than two breeds involved. There are n(n - 1)/2 heterosis combinations for n breeds. Let’s expand the code to accommodate more than two breeds via two “for” loops.

ped <- ped.init
brdcols <- c("A", "B", "O")
stopifnot(all(rowSums(ped[,brdcols]) == 1)) # Check that breed proportions sum to 1
brdcols <- brdcols[brdcols != "O"] # Breed O is not needed
stopifnot(length(brdcols) > 1)
stopifnot(brdcols %in% colnames(ped))

ped <-merge(ped, ped[,c("ID", brdcols)], by.x="SIRE", by.y="ID", all.x=TRUE)
# Rename columns after merging
colends.x <- grepl("\\.x", colnames(ped))
colends.y <- grepl("\\.y", colnames(ped))
colnames(ped)[colends.x] <- substr(colnames(ped)[colends.x], 1, nchar(colnames(ped)[colends.x])-2)
colnames(ped)[colends.y] <- paste0("s", substr(colnames(ped)[colends.y], 1, nchar(colnames(ped)[colends.y])-2))

ped <-merge(ped, ped[,c("ID", brdcols)], by.x="DAM", by.y="ID", all.x=TRUE)
# Rename columns after merging
colends.x <- grepl("\\.x", colnames(ped))
colends.y <- grepl("\\.y", colnames(ped))
colnames(ped)[colends.x] <- substr(colnames(ped)[colends.x], 1, nchar(colnames(ped)[colends.x])-2)
colnames(ped)[colends.y] <- paste0("d", substr(colnames(ped)[colends.y], 1, nchar(colnames(ped)[colends.y])-2))
# Set missing to 0
ped[is.na(ped)] <- 0
# Reorder pedigree rows
ped <- ped[order(ped$ID),]

brdcols2 <- brdcols
for(brd1 in brdcols)
{
    brdcols2 <- brdcols2[-1]
    for(brd2 in brdcols2)
    {
        # Set parents' breed for crossbred founders
        ped[ped$SIRE==0 & ped$DAM==0 & ped[,brd1]==0.5 & ped[,brd2]==0.5, paste0(c("s", "d"), c(brd1, brd2))] <- 1
        # Calculate heterosis
        sA <- paste0("s", brd1)
        dA <- paste0("d", brd1)
        sB <- paste0("s", brd2)
        dB <- paste0("d", brd2)
        ped[,paste0("het", brd1, brd2)] <- ped[,sA] * ped[,dB] + ped[,sB] * ped[,dA]
        # Calculate recombination loss
        ped[,paste0("recomb", brd1, brd2)] <- (ped[,sA] + ped[,dA]) * (ped[,sB] + ped[,dB]) - ped[,paste0("het", brd1, brd2)]
    }
}
ped

   DAM SIRE ID    A     B     O  sA   sB  dA  dB hetAB recombAB
1    0    0  1 1.00 0.000 0.000 0.0 0.00 0.0 0.0  0.00    0.000
2    0    0  2 0.00 1.000 0.000 0.0 0.00 0.0 0.0  0.00    0.000
4    0    1  3 0.50 0.000 0.500 1.0 0.00 0.0 0.0  0.00    0.000
7    2    1  4 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
6    2    1  5 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
3    0    0  6 0.50 0.500 0.000 1.0 0.00 0.0 1.0  1.00    0.000
8    4    3  7 0.50 0.250 0.250 0.5 0.00 0.5 0.5  0.25    0.250
9    6    5  8 0.50 0.500 0.000 0.5 0.50 0.5 0.5  0.50    0.500
5    0    7  9 0.25 0.125 0.625 0.5 0.25 0.0 0.0  0.00    0.125
10   8    0 10 0.25 0.250 0.500 0.0 0.00 0.5 0.5  0.00    0.250

I hope you found this demo useful.