Imaginings of a Livestock Geneticist

Recursive Method to Generate A
(Cholesky decomposition)

Similar to the previous section, Recursive Method to Create A, this section describes a method to generate A based on constructing a lower triangular matrix, referred to as L, that is outlined in Henderson (1976). When utilizing a relationship matrix in the mixed model equations the A-1 is needed, instead of A. When this paper was written computing the inverse for large matrices was not possible due the computing power being much lower than it is at the current time. As a result, methods to recursively generate A-1 were needed instead of inverting the matrix based on Gauss elimination methods. Methods to generate A-1 will be discussed in the next two sections, but both of these methods and advances in algorithms to generate A-1 all originate with the following method.
    A triangular matrix is a matrix with either the upper or lower elements of a matrix having a value of 0 and have appealing computational properties, when inverting the matrix. Dr. Henderson doesn't explicitly say it in the paper, but L is the result of doing a Cholesky decomposition of A (try "chol(A)" in R). For example, in the code outlined below only the diagonal and lower elements of matrix are computed when generating L. Once L is generated, A is calculated based on the following:

$$ A = LL' $$

The diagonal elements of L, referred to as D, are based on the inbreeding value of the parents and are needed when generating A-1 in an inbred population. As outlined in the paper D values can be generated based on the inbreeding value of the parents if they are known. These numbers can be easily verified by calculating D based on the inbreeding values of the parents derived from the Recursive Method to Create A section. The formulas are outlined below (p = parent; s = sire; d = dam):

Similar to the previous method, the pedigree has to be sorted so that parents come before progeny. Due to the pedigree being sorted, the recursive method has a very simple form. Look at the 'Recursive Method to Create G Inverse" sections R Code to see the impact of using a pedigree that isn't sorted. Lastly, if animals are numbered from 1 to the very last animal then the sire and dam values can be used to index where the respective elements are located within the algorithm. Algorithms to sort and renumber a pedigree is outlined in this section.

The following pedigree file from Henderson (1976) can be utilized with the R code above. The columns are animal, sire and dam and the pedigree is already ordered so that parents come before progeny. Lastly the animals go from 1 to the total number of animals. Outlined below is what L looks like at the end of each iteration of the for loop.

L: Loop Iteration 1
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 0.0 - - - - -
3 0.0 0.0 0.0 - - - -
4 0.0 0.0 0.0 0.0 - - -
5 0.0 0.0 0.0 0.0 0.0 - -
6 0.0 0.0 0.0 0.0 0.0 0.0 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 2
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.0 0.0 0.0 - - - -
4 0.0 0.0 0.0 0.0 - - -
5 0.0 0.0 0.0 0.0 0.0 - -
6 0.0 0.0 0.0 0.0 0.0 0.0 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 3
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.5 0.0 0.866 - - - -
4 0.0 0.0 0.0 0.0 - - -
5 0.0 0.0 0.0 0.0 0.0 - -
6 0.0 0.0 0.0 0.0 0.0 0.0 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 4
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.5 0.0 0.8660 - - - -
4 0.5 0.5 0.0 0.7071 - - -
5 0.0 0.0 0.0 0.0 0.0 - -
6 0.0 0.0 0.0 0.0 0.0 0.0 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 5
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.5 0.0 0.8660 - - - -
4 0.5 0.5 0.0 0.7071 - - -
5 0.5 0.25 0.4330 0.3536 0.7071 - -
6 0.0 0.0 0.0 0.0 0.0 0.0 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 6
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.5 0.0 0.8660 - - - -
4 0.5 0.5 0.0 0.7071 - - -
5 0.5 0.25 0.4330 0.3536 0.7071 - -
6 0.75 0.25 0.0 0.3536 0.0 0.7071 -
7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
L: Loop Iteration 7
1 2 3 4 5 6 7
1 1.0 - - - - - -
2 0.0 1.0 - - - - -
3 0.5 0.0 0.8660 - - - -
4 0.5 0.5 0.0 0.7071 - - -
5 0.5 0.25 0.4330 0.3536 0.7071 - -
6 0.75 0.25 0.0 0.3536 0.0 0.7071 -
7 0.625 0.25 0.2165 0.3536 0.3536 0.3536 0.6374

References