Exercise: Matrix calculations in Julia
6. Exercise: Matrix calculations in Julia#
\[
\newcommand{\E}{\text{E}}
\newcommand{\mbf}{\mathbf}
\newcommand{\bs}{\boldsymbol}
\newcommand{\Var}{\text{Var}}
\newcommand{\Cov}{\text{Cov}}
\newcommand{\e}{\frac{1}{\sigma^2_e}}
\newcommand{\f}{\frac{1}{\sigma^2_{\alpha}}}
\]
This exercise will provide an opportunity to practice your matrix calculation skills using Julia in the context of examining how \(\mbf{A}_i\) can be expressed as:
\[\begin{split}
\begin{align}
\begin{bmatrix}
\mbf{A}_{i-1} & \mbf{A}_{i-1}\mbf{q}_i \\
\mbf{q}'_i\mbf{A}_{i-1} & \frac{1}{2}(2 + a_{m_if_i})
\end{bmatrix}
&=
\begin{bmatrix}
\mbf{I} & \mbf{0} \\
\mbf{q}'_i& 1
\end{bmatrix}
\begin{bmatrix}
\mbf{A}_{i-1} & \mbf{0} \\
\mbf{0} & a_{ii} - \mbf{q}'_i\mbf{A}_{i-1}\mbf{q}_i
\end{bmatrix}
\begin{bmatrix}
\mbf{I} & \mbf{q}_i \\
\mbf{0}& 1
\end{bmatrix} \\
&=
\mbf{LDL}'.
\end{align}
\end{split}\]
This identity is key to deriving the efficient algorithm to invert the additive relationship matrix. It can also be used to derive a recurrsive algorithm to invert any positive-definite matrix.
Compute the matrix:
\[\begin{split}
\mathbf{B} =
\begin{bmatrix}
\mbf{I} \\
\mbf{q}'
\end{bmatrix}
\mathbf{A}_3
\begin{bmatrix}
\mbf{I} & \mbf{q}
\end{bmatrix}
\end{split}\]
for the pedigree given in class.
Verify that
\[\begin{split}
\mathbf{B} =
\begin{bmatrix}
\mbf{A}_3 & \mbf{A}_3\mbf{q}\\
\mbf{q}'\mbf{A}_3 & \mbf{q}'\mbf{A}_3\mbf{q}
\end{bmatrix}
\end{split}\]
Compute the difference between \(\mbf{A}_4\) and \(\mbf{B}\)
Verify that
\[\begin{split}
\begin{align}
\mbf{A}_4 &=
\mbf{B} +
\begin{bmatrix}
\mbf{0} & \mbf{0}\\
\mbf{0} & a_{44} - \mbf{q}'\mbf{A}_3\mbf{q}
\end{bmatrix}\\
&= \mbf{B} +
\begin{bmatrix}
\mbf{0} \\
1
\end{bmatrix}
(a_{44} - \mbf{q}'\mbf{A}_3\mbf{q})
\begin{bmatrix}
\mbf{0} & 1
\end{bmatrix}\\
&= \mbf{B} + \mbf{C}
\end{align}
\end{split}\]
where
\[\begin{split}
\mbf{C} =
\begin{bmatrix}
\mbf{0} \\
1
\end{bmatrix}
(a_{44} - \mbf{q}'\mbf{A}_3\mbf{q})
\begin{bmatrix}
\mbf{0} & 1
\end{bmatrix}
\end{split}\]
Verify that
\[\begin{split}
\begin{align}
\mbf{B} + \mbf{C}
&=
\begin{bmatrix}
\mbf{I} & \mbf{0} \\
\mbf{q}'& 1
\end{bmatrix}
\begin{bmatrix}
\mbf{A}_3 & \mbf{0} \\
\mbf{0} & a_{44} - \mbf{q}'\mbf{A}_3\mbf{q}
\end{bmatrix}
\begin{bmatrix}
\mbf{I} & \mbf{q} \\
\mbf{0}& 1
\end{bmatrix} \\
\end{align}
\end{split}\]
Answer to 1.
using LinearAlgebra
### A_2
A2 = [1 0
0 1]
q = [0.5, 0.5]
amf = A2[1,2]
### A_3
A3 = [A2 A2*q
q'A2 1 + 0.5amf ]
q = [0.5, 0.5, 0]
amf = A3[1,2]
### A_4
A4 = [A3 A3*q
q'A3 1 + 0.5amf ]
### Computing B
P = [I
q']
B1 = P*A3*P'
4×4 Matrix{Float64}:
1.0 0.0 0.5 0.5
0.0 1.0 0.5 0.5
0.5 0.5 1.0 0.5
0.5 0.5 0.5 0.5
Answer to 2.
B2 = [A3 A3*q
q'A3 q'A3*q]
4×4 Matrix{Float64}:
1.0 0.0 0.5 0.5
0.0 1.0 0.5 0.5
0.5 0.5 1.0 0.5
0.5 0.5 0.5 0.5
B1 - B2
4×4 Matrix{Float64}:
0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0