class: center, middle, inverse, title-slide .title[ # 2. Linear Algebra ] .subtitle[ ## AKA: The important foundation for everything that follows! ] .author[ ### ] --- ### Linear Algebra: Single Goal (maybe not realized in this lecture, but in this week) Understand why these equations are foundational, and used in geometric morphometrics (GM) and multivariate statistics: `$$\mathbf{AZ}=\mathbf{UDV}^T$$` or `$$\mathbf{ZA}=\mathbf{VDU}^T$$` --- ### Overview: - Scalars and Vectors - Vector addition and subtraction - Vector multiplication - Matrices (collections of vectors) - Special Matrices - Matrix inversion - Special Matrix Properties - General Linear Model Preview - Decompositions - Summary --- .center[ # Part I. Basics ## (vector and matrix operations) ] --- ### Scalars and Vectors Scalars `$$a = 5$$` `$$b = 2$$` `$$c = -3$$` Vectors `$$\mathbf{a} = \begin{pmatrix} 5 \\ 3 \\ -2 \end{pmatrix}$$` `$$\mathbf{a}^T = \begin{pmatrix} 5 & 3 & -2 \end{pmatrix}$$` The superscript, `\(^T\)`, means vector transpose. Note, `\(\left(\mathbf{a}^T \right)^T = \mathbf{a}\)` --- ### Scalars and Vectors (cont.) Vector addition/subtraction (consistent orientation is important) .pull-left[ `$$\mathbf{a} = \begin{pmatrix} 5 \\ 3 \\ -2 \end{pmatrix}$$` `$$\mathbf{b}^T = \begin{pmatrix} -5 & 0 & 1 \end{pmatrix}$$` ] .pull.right[ `$$\mathbf{a}^T + \mathbf{b}^T= \begin{pmatrix} 0 & 3 & -1 \end{pmatrix}$$` `$$\mathbf{a} - \mathbf{b} = \begin{pmatrix} 10 \\ 3 \\ -3 \end{pmatrix}$$` ] --- ### Vector Multiplication #### Not this! `$$\mathbf{a}^T = \begin{pmatrix} 5 & 3 & -2 \end{pmatrix}$$` `$$\mathbf{b}^T = \begin{pmatrix} -5 & 0 & 1 \end{pmatrix}$$` `$$\mathbf{a}^T \times \mathbf{b}^T= \begin{pmatrix} -25 & 0 & -2 \end{pmatrix}$$` But it is kind of that, in part... --- ### Vector Multiplication (cont.) Vector multiplication is both multiplication and summation of scalars. Before performing vector multiplication, there has to be consistency with **inner dimensions**. - Dimensions `\(\mathbf{a}: 3 \times 1\)` - Dimensions `\(\mathbf{a}^T: 1 \times 3\)` - Dimensions `\(\mathbf{b}: 3 \times 1\)` - Dimensions `\(\mathbf{b}^T: 1 \times 3\)` For any vector product, arrange the dimensions for the attempted product and see if inner dimensions match; e.g., `$$\mathbf{a} \times \mathbf{b}: 3 \times \color{blue} {1 \times 3} \times 1$$` Does not match, so this product is not possible --- ### Vector Multiplication (cont.) Vector multiplication is both multiplication and summation of scalars. Before performing vector multiplication, there has to be consistency with **inner dimensions**. - Dimensions `\(\mathbf{a}: 3 \times 1\)` - Dimensions `\(\mathbf{a}^T: 1 \times 3\)` - Dimensions `\(\mathbf{b}: 3 \times 1\)` - Dimensions `\(\mathbf{b}^T: 1 \times 3\)` For any vector product, arrange the dimensions for the attempted product and see if inner dimensions match; e.g., `$$\mathbf{a}^T \times \mathbf{b}: 1 \times \color{blue} {3 \times 3} \times 1$$` **This matches!** Multiplication is possible. --- ### Vector Multiplication (cont.): Vector Inner-products If inner dimensions match, the multiplication is possible. The product has dimensions equal to the outer dimensions of the match. There are two types of products: - Inner-product: the vector product results is a single scalar, i.e., the outer dimensions are `\(1 \times 1\)` - Outer-product: the vector product results in a series of scalar products with number and arrangement defined by outer dimensions --- ### Vector Multiplication (cont.): Vector Inner-products (cont.) #### Vector Inner-Product For `\(n \times 1\)` vectors `$$\mathbf{a}^T \mathbf{b} = \sum_{i=1}^n a_ib_i$$` For example, `$$\mathbf{a}^T = \begin{pmatrix} 5 & 3 & -2 \end{pmatrix}$$` `$$\mathbf{b} = \begin{pmatrix} -5 \\ 0 \\ 1 \end{pmatrix}$$` `$$\mathbf{a}^T \mathbf{b} = \left(5\times -5 \right) + \left(3 \times 0 \right) + \left(-2 \times 1 \right) = -27$$` ###### *Helpful mnemonic: "Run Amok Computational Demon!" for row-across and column-down pattern (will be helpful for matrices)* --- ### Vector Multiplication (cont.): Vector Outer Products If inner dimensions match, the multiplication is possible. The product has dimensions equal to the outer dimensions of the match. There are two types of products: - Inner-product: the vector product results is a single scalar, i.e., the outer dimensions are `\(1 \times 1\)` - Outer-product: the vector product results in a series of scalar products with number and arrangement defined by outer dimensions --- #### Vector Outer-Product For `\(n \times 1\)` vectors `$$\mathbf{a}\mathbf{b}^T = \begin{pmatrix} a_1b_1 & a_1b_2 & a_1b_3 & \cdots & a_1b_n\\ a_2b_1 & a_2b_2 & a_2b_3 & \cdots & a_2b_n\\ a_3b_1 & a_3b_2 & a_3b_3 & \cdots & a_3b_n\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ a_nb_1 & a_nb_2 & a_nb_3 & \cdots & a_nb_n\\ \end{pmatrix} = \mathbf{M}_{n \times n}$$` --- ### Vector Multiplication (cont.): Vector Outer Products (cont.) ### Important Notes `$$\mathbf{a}\mathbf{b}^T = \begin{pmatrix} a_1b_1 & a_1b_2 & a_1b_3 & \cdots & a_1b_n\\ a_2b_1 & a_2b_2 & a_2b_3 & \cdots & a_2b_n\\ a_3b_1 & a_3b_2 & a_3b_3 & \cdots & a_3b_n\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ a_nb_1 & a_nb_2 & a_nb_3 & \cdots & a_nb_n\\ \end{pmatrix} = \mathbf{M}_{n \times n}$$` - `\(\mathbf{M}\)` is called a matrix, with dimensions defined - `\(\mathbf{M}\)` comprises row vectors and column vectors - because `\(\mathbf{M}\)` has the same number of rows and columns, it also has a diagonal vector - the sum of the diagonal vector, called the **trace**, is the inner-product of `\(\mathbf{a}\)` and `\(\mathbf{b}\)` --- ### Mike's simple definition of a matrix ### *A matrix is a template for algorithms used on vectors.* --- ### Matrices A matrix is more than a vector outer-product (more precisely, a vector outer-product is merely one type of matrix). A matrix is a collection of vectors, arranged in a specific way, such that linear algebra operations can be carried out, systematically. For example, the arrangement of vectors in a matrix indicates how to find a series of inner-products and display their results. The most basic matrices for statistics are **data frames**. In data frames, column vectors are variables and row vectors are observations. We can demonstrate this easily in `R` with the `data.frame` and `matrix` functions. ``` ## y1 y2 ## obs.1 2 -1 ## obs.2 3 -2 ## obs.3 2 0 ## obs.4 5 0 ## obs.5 6 1 ## obs.6 8 -1 ``` --- ### Matrices (Cont.) ``` ## y1 y2 ## obs.1 2 -1 ## obs.2 3 -2 ## obs.3 2 0 ## obs.4 5 0 ## obs.5 6 1 ## obs.6 8 -1 ``` Which is R's way of saying, `$$\small{\mathbf{Y} = \begin{pmatrix} 2 & -1 \\ 3 & -2 \\ 2 & 0 \\ 5 & 0 \\ 6 & 1 \\ 8 & -1 \\ \end{pmatrix}}$$` `\(\mathbf{Y}\)` is a matrix and a data frame. Row vectors are observations for the variables represented as column vectors. The dimensions of `\(\mathbf{Y}\)` are `\(n \times p\)` for the `\(n\)` observations of subjects for `\(p\)` variables. --- ### Matrices (Cont.) #### This time with additional `R` code ``` r y1 <- c(2, 3, 2, 5, 6, 8) y2 <- c(-1, -2, 0, 0, 1, -1) *Y <- data.frame(y1 = y1, y2 = y2) rownames (Y) <- paste("obs", 1:6, sep = ".") Y ``` ``` ## y1 y2 ## obs.1 2 -1 ## obs.2 3 -2 ## obs.3 2 0 ## obs.4 5 0 ## obs.5 6 1 ## obs.6 8 -1 ``` --- ### Matrix Addition, Subtraction, Multiplication .pull-left[ #### Matrix addition and subtraction Matrix addition and subtraction is no different than vector addition and subtraction. Matrices have to have *commensurate* dimensions (same `\(n \times p\)`). #### Matrix multiplication Matrix multiplication is nothing more than systematic calculation of vector inner-products and arrangement of these into precise corresponding elements of a new matrix. Like vectors, inner dimensions must match and the product is defined by the outer dimensions. The simplest way to define this is as follows ] .pull-right[ Let `\(\mathbf{X}\)` be an `\(n \times k\)` matrix and let `\(\mathbf{Y}\)` by an `\(n \times p\)` matrix, such that `$$\small \mathbf{X} = \begin{pmatrix} \mathbf{x}_1 & \mathbf{x}_2 & \cdots & \mathbf{x}_k \end{pmatrix}$$` and `$$\small \mathbf{Y} = \begin{pmatrix} \mathbf{y}_1 & \mathbf{y}_2 & \cdots & \mathbf{y}_p \end{pmatrix}$$` Note that each `\(\mathbf{x}_i\)` column vector is `\(n \times 1\)` in dimension for `\(k\)` vectors and each `\(\mathbf{y}_i\)` is `\(n \times 1\)` in dimension for `\(p\)` vectors. Thus, ] --- ### Matrix Addition, Subtraction, Multiplication (Cont.) `$$\mathbf{X} = \begin{pmatrix} \mathbf{x}_1 & \mathbf{x}_2 & \cdots & \mathbf{x}_k \end{pmatrix}$$` and `$$\mathbf{Y} = \begin{pmatrix} \mathbf{y}_1 & \mathbf{y}_2 & \cdots & \mathbf{y}_p \end{pmatrix}$$` and `$$\mathbf{X}^T\mathbf{Y} = \begin{pmatrix} \mathbf{x}_1^T\mathbf{y}_1 & \mathbf{x}_1^T\mathbf{y}_2 & \cdots & \mathbf{x}_1^T\mathbf{y}_p\\ \mathbf{x}_2^T\mathbf{y}_1 & \mathbf{x}_2^T\mathbf{y}_2 & \cdots & \mathbf{x}_2^T\mathbf{y}_p\\ \vdots & \vdots & & \vdots\\ \mathbf{x}_k^T\mathbf{y}_1 & \mathbf{x}_k^T\mathbf{y}_2 & \cdots & \mathbf{x}_k^T\mathbf{y}_p\\ \end{pmatrix}$$` The matrix product is a matrix with `\(k \times p\)` inner-products. For the product to work, both matrices must have the same numer of rows, `\(n\)`. --- ### Matrix Multiplication (Cont.) #### Matrix multiplication example (using R script) ``` r n <- 50 p <- 3 y1 <- rnorm(50) y2 <- rnorm(50) y3 <- rnorm(50) Y <- as.matrix(data.frame(y1 = y1, y2 = y2, y3 = y3)) rownames (Y) <- paste("obs", 1:n, sep = ".") dim(Y) ``` ``` ## [1] 50 3 ``` --- ### Matrix Multiplication (Cont.) #### Example data frame, `\(\mathbf{Y}\)` <div class="datatables html-widget html-fill-item" id="htmlwidget-1a648e86760617d8a29a" style="width:100%;height:auto;"></div> <script type="application/json" data-for="htmlwidget-1a648e86760617d8a29a">{"x":{"filter":"none","vertical":false,"fillContainer":false,"data":[["obs.1","obs.2","obs.3","obs.4","obs.5","obs.6","obs.7","obs.8","obs.9","obs.10","obs.11","obs.12","obs.13","obs.14","obs.15","obs.16","obs.17","obs.18","obs.19","obs.20","obs.21","obs.22","obs.23","obs.24","obs.25","obs.26","obs.27","obs.28","obs.29","obs.30","obs.31","obs.32","obs.33","obs.34","obs.35","obs.36","obs.37","obs.38","obs.39","obs.40","obs.41","obs.42","obs.43","obs.44","obs.45","obs.46","obs.47","obs.48","obs.49","obs.50"],[1.4172,-0.7282,-0.2667,1.0735,0.723,-1.3032,0.762,1.3157,0.1616,0.5382,-1.125,-0.9247,-0.3069,0.8774,0.7183,-1.1153,-1.459,0.3562,0.1912,-0.3329,0.2793,-0.8332000000000001,2.4215,-0.6904,-0.4653,-0.0098,0.1828,0.9330000000000001,0.2655,-1.067,0.327,1.0057,0.8276,0.7262,1.0841,-1.3492,0.8317,1.4915,-0.3591,0.1575,0.0182,-1.1856,-1.2937,0.6163999999999999,0.8237,0.0788,1.17,0.0369,1.1106,-1.4118],[0.7752,0.7821,1.0205,-1.6889,0.0354,-0.8363,0.2619,1.4614,-1.0031,-0.5508999999999999,-0.1107,-1.1232,-0.5396,1.3207,1.732,0.0168,-0.3069,-2.4835,1.3488,-0.0394,0.7086,-0.2577,-1.0495,-1.2812,-0.4344,-0.1469,0.8479,0.6677999999999999,-1.2007,0.5141,0.4316,-0.4485,-0.5182,-1.6906,-0.1333,1.1896,-0.3844,0.3085,0.3534,-0.2103,0.1276,0.6246,2.2465,-0.9144,-0.1832,0.0494,-1.1033,0.5866,0.4888,1.4424],[-0.09130000000000001,-0.3958,0.8213,-0.5844,-1.1732,0.124,-0.1117,0.3667,1.8759,-1.0421,1.5979,0.8733,-1.3501,-0.6104000000000001,1.0627,0.1467,-0.7203000000000001,-0.5828,-0.2574,-0.8376,-0.8407,-0.1998,-0.2236,-0.0224,-0.2883,1.6402,-0.2349,-1.0472,0.5142,-0.9789,1.8749,0.7681,0.736,0.4188,1.14,-0.9332,-1.5754,-0.7969000000000001,-1.3347,-0.5171,-0.6405,-0.2055,1.197,0.8096,1.4685,0.0188,-0.8426,0.3668,-0.9631,1.9462]],"container":"<table class=\"display\">\n <thead>\n <tr>\n <th> <\/th>\n <th>y1<\/th>\n <th>y2<\/th>\n <th>y3<\/th>\n <\/tr>\n <\/thead>\n<\/table>","options":{"pageLength":8,"columnDefs":[{"className":"dt-right","targets":[1,2,3]},{"orderable":false,"targets":0},{"name":" ","targets":0},{"name":"y1","targets":1},{"name":"y2","targets":2},{"name":"y3","targets":3}],"order":[],"autoWidth":false,"orderClasses":false,"lengthMenu":[8,10,25,50,100]}},"evals":[],"jsHooks":[]}</script> --- ### Matrix Multiplication (Cont.) #### Matrix multiplication example (using R script) ``` r X <- cbind(Intercept = rep(1, n), group = rep(c(0, 1), n / 2)) rownames(X) <- rownames(Y) ``` <div class="datatables html-widget html-fill-item" id="htmlwidget-3f68752ae0db151dfa62" style="width:100%;height:auto;"></div> <script type="application/json" data-for="htmlwidget-3f68752ae0db151dfa62">{"x":{"filter":"none","vertical":false,"fillContainer":false,"data":[["obs.1","obs.2","obs.3","obs.4","obs.5","obs.6","obs.7","obs.8","obs.9","obs.10","obs.11","obs.12","obs.13","obs.14","obs.15","obs.16","obs.17","obs.18","obs.19","obs.20","obs.21","obs.22","obs.23","obs.24","obs.25","obs.26","obs.27","obs.28","obs.29","obs.30","obs.31","obs.32","obs.33","obs.34","obs.35","obs.36","obs.37","obs.38","obs.39","obs.40","obs.41","obs.42","obs.43","obs.44","obs.45","obs.46","obs.47","obs.48","obs.49","obs.50"],[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1],[0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1]],"container":"<table class=\"display\">\n <thead>\n <tr>\n <th> <\/th>\n <th>Intercept<\/th>\n <th>group<\/th>\n <\/tr>\n <\/thead>\n<\/table>","options":{"pageLength":6,"columnDefs":[{"className":"dt-right","targets":[1,2]},{"orderable":false,"targets":0},{"name":" ","targets":0},{"name":"Intercept","targets":1},{"name":"group","targets":2}],"order":[],"autoWidth":false,"orderClasses":false,"lengthMenu":[6,10,25,50,100]}},"evals":[],"jsHooks":[]}</script> --- ### Matrix Multiplication (Cont.) #### Matrix multiplication example (using R script) .med[ ``` r t(X) %*% Y ``` ``` ## y1 y2 y3 ## Intercept 6.295172 0.7029465 0.3656386 ## group -1.744153 -2.7077558 -1.2747868 ``` ``` r crossprod(X, Y) ``` ``` ## y1 y2 y3 ## Intercept 6.295172 0.7029465 0.3656386 ## group -1.744153 -2.7077558 -1.2747868 ``` ] + Notice that the row names are the column names of `\(\mathbf{X}\)` and the column names are the column names of `\(\mathbf{Y}\)`. + Also notice the `crossprod` function is short for *matrix of cross-products*, which we will refer to simply as the **matrix cross-product**. This is a bit unfortunate nomenclature, as the meaning is different than *vector cross-product*, a vector produced by multiplication in a certain way of two vectors, different from the *dot-product*, which is basically the same as the inner-product. (Try not to let the terminology get to you!) --- ### ~~Matrix Division~~ This is not something we can do. We can invert some matrices, but we will come back to this in a moment. --- ### Important Matrix and Vector Multiplication Properties - Any matrix or vector multiplied by a scalar multiplies every element by the scalar. The name, "scalar," means that every element is scaled. - Like vectors that have inner-products (dot products) or outer-products, sometimes matrices can be similar. The tendency though is to call these cross-product matrices. - A vector inner-product of the form, `\(\mathbf{x}^T\mathbf{x}\)` is a generalized method of squaring. It is precisely the sum of squared values, rather than just a squared value. - A matrix **.blue[cross-product]** of the form, `\(\mathbf{X}^T\mathbf{X}\)` is a generalized method of squaring. This matrix has summed squared values of column vectors of `\(\mathbf{X}\)` on a diagonal vector and summed cross-products between vectors in the off-diagonal elements. - A matrix **.blue[outer-product]** of the form, `\(\mathbf{X}\mathbf{X}^T\)` can also be considered a generalized method of squaring. This matrix has summed squared values of row vectors of `\(\mathbf{X}\)` on a diagonal vector and summed cross-products between vectors in the off-diagonal elements. - The **.blue[trace]** of a symmetric matrix is the sum of diagonal elements. One convenience is `\(trace(\mathbf{X}^T\mathbf{X}) = trace(\mathbf{X}\mathbf{X}^T)\)`. --- ### Important Matrix and Vector Multiplication Properties (cont.) - A matrix is **.blue[square]** when it has equal rows and columns, e.g., `\(n \times n\)`. - If `\(\mathbf{XY}\)` is possible and `\(\mathbf{YX}\)` is possible, they will not have the same product unless both matrices are square and `\(\mathbf{X} = \mathbf{Y}\)`. The dimensions of matrix products are also likely different. - The next slide will describe some special matrices. --- .center[ # Part II. Special Matrices ## And some properties that make them special ] --- ### Special Matrices | Matrix and description | Example | Comments | | ------ | ---- | ------ | |**Square matrix:** Any matrix with the same number of rows and columns, e.g., `\(n \times n\)` | `\(\small{\mathbf{S} = \begin{pmatrix} 2 & 3 & 1 \\ 4 & 7 & 5\\ 1 & 10 & 0 \\ \end{pmatrix}}\)` | These matrices, themselves, are generally not too compelling, unless there is a relationship between the row and column variables.| |**Symmetric matrix:** A square matrix that has symmetry above and below the diagonal| `\(\small{\mathbf{S} = \begin{pmatrix} 2 & 3 & 1 \\ 3 & 7 & 10\\ 1 & 10 & 0 \\ \end{pmatrix}}\)` | These matrices have the condition that `\(\mathbf{S} = \mathbf{S}^T\)`.| |**Diagonal matrix:** A symmetric matrix that has non-0 elements only along the diagonal; all others are 0.| `\(\small{\mathbf{D} = \begin{pmatrix} 2 & 0 & 0 \\ 0 & 7 & 0\\ 0 & 0 & 0 \\ \end{pmatrix}}\)` | Note that the diagonal elements can still be 0. If all diagonal elements are equal, we can say `\(\mathbf{D} = a\mathbf{I}\)`, for a scalar, `\(a\)`, and an identity matrix, `\(\mathbf{I}\)`.| |**Identity matrix:** A diagonal matrix that has diagonal elements equal to 1.| `\(\small{\mathbf{I} = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0\\ 0 & 0 & 1 \\ \end{pmatrix}}\)` | This is the matrix equivalent to a value of 1, such that `\(\mathbf{AI} = \mathbf{A}\)` --- ### Special Matrices (cont.) | Matrix and description | Example | Comments | | ------ | ---- | ------ | |**Orthogonal matrix:** A symmetric matrix that has the property that `\(\mathbf{AA}^T = \mathbf{A}^T\mathbf{A} = \mathbf{I}\)`.| `\(\small{\mathbf{A} = \begin{pmatrix} 0.7071 & 0.7071 \\ 0.7071 & -0.7071 \end{pmatrix}}\)` | `\(\small{\mathbf{AA}^T} = \begin{pmatrix} 0.7071 & 0.7071 \\ 0.7071 & -0.7071 \end{pmatrix} \begin{pmatrix} 0.7071 & 0.7071 \\ 0.7071 & -0.7071 \end{pmatrix} = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}\)` | |**Orthonormal matrix:** An orthogonal matrix with unit vectors (length of 1) `\(\mathbf{AA}^T = \mathbf{A}^T\mathbf{A} = \mathbf{I}\)`.| `\(\small{\mathbf{A} = \begin{pmatrix} 0.7071 & 0.7071 \\ 0.7071 & -0.7071 \end{pmatrix}}\)` | `\(\small{\lvert\lvert\mathbf{a}_1 \rvert\rvert =(\mathbf{a}_1^T\mathbf{a}_1)^{1/2} = (0.7071^2 +0.7071^2)^{1/2} =1}\)` `\(\small{\lvert\lvert\mathbf{a}_2\rvert\rvert=(\mathbf{a}_2^T\mathbf{a}_2)^{1/2} = (0.7071^2 + (-0.7071)^2)^{1/2} =1}\)`| |**Idempotent matrix:** A symmetric matrix with the property `\(\mathbf{HH} = \mathbf{H}\)` | `\(\small{\mathbf{H} = \begin{pmatrix} 0.5 & 0.5 \\ 0.5 & 0.5 \end{pmatrix}}\)` | `\(\begin{align*} \small{\mathbf{HH}} = \begin{pmatrix} 0.5 & 0.5 \\ 0.5 & 0.5 \end{pmatrix} \begin{pmatrix} 0.5 & 0.5 \\ 0.5 & 0.5 \end{pmatrix} \\ = \small{\begin{pmatrix} 0.5^2 + 0.5^2 & 0.5^2 + 0.5^2 \\ 0.5^2 + 0.5^2 & 0.5^2 + 0.5^2 \end{pmatrix} = \begin{pmatrix} 0.5 & 0.5 \\ 0.5 & 0.5 \end{pmatrix}} \end{align*}\)` ***Idempotency is a property that is important for orthogonal projection, covered next*** --- ### Special Matrices (cont.) and Special matrix properties | Matrix and description | Example | Comments | | ------ | ---- | ------ | |**Projection matrix:** Any matrix that can be multiplied with data that map the data to new vectors (transformation), taking the form, `\(\mathbf{PY}\)` or `\(\mathbf{YP}\)`.| `\(\small{\mathbf{P} = \begin{pmatrix} 0.5 & 0.5 \\ 0.5 & 0.5 \end{pmatrix}}\)` | Although projection (transformation) matrices can have many forms, typical ones are: "hat" matrices, `\(\mathbf{H}\)`, which if idempotent perform orthogonal projections, but if merely symmetric perform shears (more on this soon); and eigenvectors ( `\(\mathbf{V}\)` ). |**Positive-definite matrix:** If `\(\mathbf{X}\)` is symmetric, `\(\mathbf{a}^T \mathbf{Xa}\)` is a positive scalar if `\(\mathbf{a}\)` has non-zero values. Thus, `\(\mathbf{X}\)` is positive-definite.| `\(\small{\mathbf{X} = \begin{pmatrix} 1 & 0.5 \\ 0.5 & 2 \end{pmatrix}}\)` | This also means all eigenvalues of `\(\mathbf{X}\)` are non-negative (more later). ***Covariance matrices for shape data are often not positive-definite because (1) generalized Procrustes analysis creates redundancies in the data and (2) there are often more shape variables than observations. We will cover these issues in much detail, later*** It should be clear that matrices can have several properties. For example, a *projection matrix* can be symmetric and idempotent. --- ### Special Matrices (cont.) and Special matrix properties #### Some additional notes - Some sources refer to "triangular" matrices. These might be symmetric matrices, where the elements above a diagonal match the elements below. However, these might refer to a matrix that looks like a symmetric matrix but values below the diagonal are all 0. - A **matrix cross-product** involving a single matrix, e.g., `\(\mathbf{X}^T\mathbf{X}\)`, will always produce a symmetric matrix. --- .center[ # Part III. Matrix projection ## This needs further discussion because it is super important! ] --- Let's revisit this, focusing on the left sides: ### Linear Algebra: Single Goal (maybe not realized in this lecture, but in this week) Understand why these equations are foundational, and used in geometric morphometrics (GM) and multivariate statistics: `$$\mathbf{AZ}=\mathbf{UDV}^T$$` or `$$\mathbf{ZA}=\mathbf{VDU}^T$$` Let the matrix, `\(\mathbf{A}\)`, suggest ***association*** or ***alignmnet***. Let the matrix, `\(\mathbf{Z}\)`, be a manipulated form of the data matrix, `\(\mathbf{Y}\)`. (It might be unmanipulated) ** `\(\mathbf{AZ}\)` or `\(\mathbf{ZA}\)` could represent projection.** So, let `\(\mathbf{A} = \mathbf{H}\)`, a **projection matrix**. (A projection matrix is often called a "Hat" matrix.) --- ### Matrix Projection Matrix multiplication is often used for the purpose of projection of data, `\(\mathbf{Y}\)`. In a most technical and overly mathematical sense, projection is an alignment of data in a manifold of the data space (where Euclidean geometry is at least approximately appropriate). In a simpler sense, projection is a method of rotating, stretching, flipping, and/or scaling a data space. This is the essence of linear models - finding a constrained explanation of the data. Two ways to understand projection: `\(\mathbf{HY}\)`: where `\(\mathbf{H}\)` is an `\(n \times n\)` projection matrix `\(\mathbf{YH}\)`: where `\(\mathbf{H}\)` is a `\(p \times p\)` rotation/shear/scale matrix Both of these will make more sense with exposure to linear model uses, but for now, let's consider what the latter does. --- ## Matrix Projection: Geometric Interpretations To best understand attributes of high-dimensional data spaces and algebra applied to them, consider two-dimensional data spaces and realize the algebra can be generalized to higher dimensions. #### The Result of `\(\mathbf{YH}\)`, where `\(\mathbf{H}\)` is a `\(p \times p\)` matrix `\(\mathbf{H}\)` type | Description | Result/Comment :------- | :----------- | :-------------- `\(\mathbf{I}\)` | Identity matrix | No change `\(c\mathbf{I}\)` | Scaled identity matrix | Enlarge (expand) or reduce (shrink). `\(c\mathbf{I}\)` is not idempotent. `\(\mathbf{D}\)` | Diagonal matrix | Stretching of axes (variables). `\(\mathbf{D}\)` is not idempotent. `\(\mathbf{H}_{orthogonal}\)` | `\(\mathbf{H} = \mathbf{H}^T\)` | Rigid (orthogonal) rotation of `\(\mathbf{Y}\)` (see note). `\(\mathbf{H}\)` is idempotent. `\(\mathbf{H}_{oblique}\)` | `\(\mathbf{H} = \mathbf{H}^T\)` or `\(\mathbf{H} \neq \mathbf{H}^T\)` | Shear of `\(\mathbf{Y}\)`. `\(\mathbf{H}_{oblique}\)` is not orthogonal. `\(\mathbf{H}\)` is not idempotent. --- ### Matrix Projection: Geometric Interpretations To best understand attributes of high-dimensional data spaces and algebra applied to them, consider two-dimensional data spaces and realize the algebra can be generalized to higher dimensions. #### The Result of `\(\mathbf{YH}\)`, where `\(\mathbf{H}\)` is a `\(p \times p\)` matrix **Note** that `\(\mathbf{H}_{orthogonal}\)` can be viewed as a rotation matrix, which makes sense in low-dimensional data spaces. For example, in a two-dimension data space, it can be the same as: `$$\mathbf{H}_{orthogonal} = \begin{pmatrix} \cos\theta & -\sin\theta\\ \sin\theta & \cos\theta \end{pmatrix}$$` for a rotational angle, `\(\theta\)`, in a plane. --- ### Matrix Projection: Geometric Interpretations (Cont.) <img src="LectureData/02.matrices/projections.png" alt="" width="90%" /> --- .center[ # Part IV. Matrix inversion ## A rather challenging operation that thankfully computers can perform under certain circumstances. ] --- ### Matrix Inversion - Matrices cannot be divided. Much like division of a scalar is the same as multiplication by its inverse, some matrices can be inverted. - Matrix needs to be square, *non-singular*; symmetric, non-singular is optimal. (A positive-definite matrix is non-singular but a non-singular matrix can be either positive-definite or positive semidefinite.) - Matrix inversion fulfills the property that `\(\mathbf{AA}^{-1} = \mathbf{I}\)` - Illustration is simplest with a `\(2 \times 2\)` matrix Let `\(\mathbf{A}\)` be a `\(2 \times 2\)` matrix, `$$\mathbf{A} = \begin{pmatrix} a & b \\ c & d \end{pmatrix}$$` Then, `$$\mathbf{A}^{-1} = \left| \mathbf{A} \right|^{-1} \begin{pmatrix} d & -b\\ -c & a \end{pmatrix}$$` where `\(\left| \mathbf{A} \right| = ad-bc\)` If `\(|\mathbf{A}| = 0\)`, the matrix is *singular* and cannot be inverted (dividing by 0 from matrices). --- ### Matrix Inversion (Cont.) #### Example `$$\mathbf{A} = \begin{pmatrix} 3 & 4 \\ 4 & 6 \end{pmatrix}$$` `$$\left| \mathbf{A} \right| = 18-16$$` `$$\left| \mathbf{A} \right|^{-1} = 0.5$$` `$$\mathbf{A}^{-1} = \left| \mathbf{A} \right|^{-1} \begin{pmatrix} d & -b\\ -c & a \end{pmatrix} = 0.5 \begin{pmatrix} 6 & -4\\ -4 & 3 \end{pmatrix} = \begin{pmatrix} 3 & -2\\ -2 & 1.5 \end{pmatrix}$$` --- ### Matrix Inversion (Cont.) #### Example Confirm (using R script) .pull-left[ .med[ ``` r A = matrix(c(3, 4, 4, 6), nrow = 2, ncol = 2) A ``` ``` ## [,1] [,2] ## [1,] 3 4 ## [2,] 4 6 ``` ``` r solve(A) # solve function finds inverse ``` ``` ## [,1] [,2] ## [1,] 3 -2.0 ## [2,] -2 1.5 ``` ]] .pull-right[.med[ ``` r A %*% solve(A) ``` ``` ## [,1] [,2] ## [1,] 1 0 ## [2,] 0 1 ``` ]] --- ### Matrix Inversion (Cont.) - Note, inversion with 3 × 3 matrices or larger, not a trivial exercise .pull-left[ - Let `\(\mathbf{A}\)` be a `\(3 \times 3\)` matrix, `$$\mathbf{A} = \begin{pmatrix} a & b & c\\ d & e & f\\ g & h & i\\ \end{pmatrix}$$` ] .pull-right[ - Let's also assume before proof that `\(\mathbf{A}\)` is invertible. ] ***Algorithm:*** + Find the ***nine*** possible `\(2 \times 2\)` matrices within `\(\mathbf{A}\)`, e.g., .center[ `\(\begin{pmatrix} a & b\\ d & e\\ \end{pmatrix}\)`, `\(\begin{pmatrix} a & c\\ d & f\\ \end{pmatrix}\)`, `\(\begin{pmatrix} b & c\\ e & f\\ \end{pmatrix}\)`, etc. ] --- ### Matrix Inversion (Cont.) .pull-left[ - Let `\(\mathbf{A}\)` be a `\(3 \times 3\)` matrix, `$$\mathbf{A} = \begin{pmatrix} a & b & c\\ d & e & f\\ g & h & i\\ \end{pmatrix}$$` ] .pull-right[ - Let's also assume before proof that `\(\mathbf{A}\)` is invertible. ] ***Algorithm:*** + Find the determinant of each and give them names. Every other determinant is negative. The names correspond to the rows and columns not represented by each of the ***nine*** elements: .center[ `\(A = (ei - fh)\)`, `\(D = -(bi - ch)\)`, `\(G = (bf - ce)\)`, `\(B = -(di - fg)\)`, `\(E = (ai - cg)\)`, `\(H= -(af - cd)\)`, `\(C = (dh - eg)\)`, `\(F = -(ah - bg)\)`, `\(I = (ae - bd)\)` ] --- ### Matrix Inversion (Cont.) .pull-left[ - Let `\(\mathbf{A}\)` be a `\(3 \times 3\)` matrix, `$$\mathbf{A} = \begin{pmatrix} a & b & c\\ d & e & f\\ g & h & i\\ \end{pmatrix}$$` ***Algorithm:*** - Find the determinant of `\(\mathbf{A}\)`: `\(\left| \mathbf{A} \right| = aA + bB + cC\)` ] .pull-right[ - If the determinant is greater than 0, then `\(\mathbf{A}\)` is invertible. ] - Make a new matrix: `$$\mathbf{A}^{-1} = \frac{1}{\left| \mathbf{A} \right|} \begin{pmatrix} A & B & C\\ D & E & F\\ G & H & I\\ \end{pmatrix}^T = \frac{1}{\left| \mathbf{A} \right|} \begin{pmatrix} A & D & G\\ B & E & H\\ C & F & I\\ \end{pmatrix}$$` --- ### Matrix Inversion (Cont.) For the sake of sanity, we will not consider `\(4 \times 4\)` or larger matrices! Generally speaking, software packages also do not ***solve*** matrix inverses in a way that makes sense like the previous demonstration. (Also, `solve` is a function in `R` for finding a matrix inverse.) Rather, they use a library of numerical routines that find computational solutions that converge on algebraic solutions. For example, `R` relies on the **L**inear **A**lgebra **Pack**age (LAPACK), which relies on **B**asic **L**inear **A**lgebra **S**ub-programs (BLAS), the subroutines that provide building blocks for computational routines. `\(^*\)` **Basically, understanding the concept is sufficient.** One can trust that computers and their programs can be relied upon for efficient matrix inversion, within reason. Matrix inversion is, however, a ***Big Data*** concern, as computer memory is finite. .footnote[ `\(^*\)` LAPACK has roots in `Fortran`, which is used in `R`. `R` also relies on `C` and `C++`, which also use LAPACK plus other linear alegrba packages. ] --- .center[ # Part V. Practice using some matrix operations ## There is a point to this -- trust it will work out. ] --- ### Putting Our Current Knowledge to Good Use! - Recall that a matrix must be square in order to be inverted. - Recall that data frames are `\(n \times p\)` in dimension, and it might be rare that `\(n=p\)`. - Recall that matrix cross-products with a single matrix are essentially the squaring of matrices, and they are symmetric. - **Thus**, matrix cross-products are invertible (and frequently used as a step before inversion). #### An Example (recalling an earlier example) using R scripts .pull-left[.med[ ``` r y1 <- c(2, 3, 2, 5, 6, 8) y2 <- c(-1, -2, 0, 0, 1, -1) Y <- as.matrix(data.frame(y1 = y1, y2 = y2)) rownames (Y) <- paste("obs", 1:6, sep = ".") Y ``` ``` ## y1 y2 ## obs.1 2 -1 ## obs.2 3 -2 ## obs.3 2 0 ## obs.4 5 0 ## obs.5 6 1 ## obs.6 8 -1 ``` ]] .pull-right[.med[ ``` r x1 <- rep(1, 6) x2 <- c(rep(0, 3), rep(1,3)) X <- as.matrix(data.frame(x1 = x1, x2 = x2)) rownames(X) <- rownames(Y) X ``` ``` ## x1 x2 ## obs.1 1 0 ## obs.2 1 0 ## obs.3 1 0 ## obs.4 1 1 ## obs.5 1 1 ## obs.6 1 1 ``` ]] --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts We cannot invert `\(\mathbf{X}\)`, but let's invert the cross-product of `\(\mathbf{X}\)`. ``` r crossprod(X) ``` ``` ## x1 x2 ## x1 6 3 ## x2 3 3 ``` ``` r solve(crossprod(X)) ``` ``` ## x1 x2 ## x1 0.3333333 -0.3333333 ## x2 -0.3333333 0.6666667 ``` --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts While we are on the topic… Imagine that the second column is a variable, `\(X\)`, taking on the value either 0 or 1 The columns are “dummy” variables ``` r X ``` ``` ## x1 x2 ## obs.1 1 0 ## obs.2 1 0 ## obs.3 1 0 ## obs.4 1 1 ## obs.5 1 1 ## obs.6 1 1 ``` With a little bit of mind-bending, one might notice that `$$\mathbf{X}^T\mathbf{X} = \begin{pmatrix} n & \sum X\\ \sum X & \sum X^2 \end{pmatrix}$$` --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts Let's reduce the `\(\mathbf{X}\)` matrix, such that `$$\mathbf{X}_0 = \begin{pmatrix} 1\\ 1\\ 1\\ 1\\ 1\\ 1 \end{pmatrix}$$` --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts .pull-left[.med[ And now in R, let's find the cross-product of this, and its inverse ``` r X0 <- X[,1] # just the first column crossprod(X0) ``` ``` ## [,1] ## [1,] 6 ``` ]] .pull-right[.med[ ``` r solve(crossprod(X0)) ``` ``` ## [,1] ## [1,] 0.1666667 ``` ]] Thus, `\(\mathbf{X}_0^T\mathbf{X}_0 = n\)` and `\(\left( \mathbf{X}_0^T\mathbf{X}_0 \right)^{-1} = n^{-1}\)` --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts Now for fun, let's find the cross-product between `\(\mathbf{X}_0\)` and `\(\mathbf{Y}\)` .pull-left[.med[ ``` r crossprod(X0, Y) ``` ``` ## y1 y2 ## [1,] 26 -3 ``` ]] .pull-right[.med[ If we examine `\(\mathbf{Y}\)` again, you might notice something ``` ## y1 y2 ## obs.1 2 -1 ## obs.2 3 -2 ## obs.3 2 0 ## obs.4 5 0 ## obs.5 6 1 ## obs.6 8 -1 ``` ]] `\(\mathbf{X}_0^T\mathbf{Y}\)` produces vector sums, `\(\sum y_i\)` for each variable. --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts One more step! Let's find the result of two things: `\(\left( \mathbf{X}_0^T\mathbf{X}_0 \right)^{-1}\mathbf{X}_0^T\mathbf{Y}\)` and `\(\mathbf{X}_0 \left( \mathbf{X}_0^T\mathbf{X}_0 \right)^{-1}\mathbf{X}_0^T\mathbf{Y}\)` Note, the second one *.blue[can be written]* `\(\mathbf{H}_0\mathbf{Y}\)` --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts One more step! Let's find the result of two things: `\(\left( \mathbf{X}_0^T\mathbf{X}_0 \right)^{-1}\mathbf{X}_0^T\mathbf{Y}\)` and `\(\mathbf{X}_0 \left( \mathbf{X}_0^T\mathbf{X}_0 \right)^{-1}\mathbf{X}_0^T\mathbf{Y}\)` Note, the second one *.blue[can be written]* `\(\mathbf{H}_{0}\mathbf{Y}\)` .pull-left[.med[ ``` r solve(crossprod(X0)) %*% crossprod(X0, Y) ``` ``` ## y1 y2 ## [1,] 4.333333 -0.5 ``` ``` r X0 %*% solve(crossprod(X0)) %*% crossprod(X0, Y) ``` ``` ## y1 y2 ## [1,] 4.333333 -0.5 ## [2,] 4.333333 -0.5 ## [3,] 4.333333 -0.5 ## [4,] 4.333333 -0.5 ## [5,] 4.333333 -0.5 ## [6,] 4.333333 -0.5 ``` ]] .pull-right[ **The former finds a (row) vector of variable means. The latter creates a matrix of variable means for every observation in the original data matrix, `\(\mathbf{Y}\)`!** Although it might not yet be obvious, we just used the general linear model to find fitted (predicted) values for a model design that has a single mean (centroid). ] --- ### Putting Our Current Knowledge to Good Use! (Cont.) #### An Example (recalling an earlier example) using R scripts If we now substitute `\(\mathbf{X}\)` for `\(\mathbf{X}_0\)`, this is the result (`\(\mathbf{H}_A\mathbf{Y}\)`) .pull-left[ .med[ ``` r solve(crossprod(X)) %*% crossprod(X, Y) ``` ``` ## y1 y2 ## x1 2.333333 -1 ## x2 4.000000 1 ``` ``` r X %*% solve(crossprod(X)) %*% crossprod(X, Y) ``` ``` ## y1 y2 ## obs.1 2.333333 -1 ## obs.2 2.333333 -1 ## obs.3 2.333333 -1 ## obs.4 6.333333 0 ## obs.5 6.333333 0 ## obs.6 6.333333 0 ``` ]] Recall .pull.right[ ``` ## x1 x2 ## obs.1 1 0 ## obs.2 1 0 ## obs.3 1 0 ## obs.4 1 1 ## obs.5 1 1 ## obs.6 1 1 ``` ] One can see that the second variable is an indicator variable with 0 for one group and 1 for the other (two groups each with three observations). This approach allows one to estimate two group means and create a matrix of fitted values that has two means instead of one. --- ### Call back to Matrix Projections! Statistical analysis can largely be thought of as the comparison of alternate projections (transformations) of data, based on linear models. For example, we might compare the projections between two models: `$$\mathbf{H}_{A}\mathbf{Y} - \mathbf{H}_{0}\mathbf{Y} = \Delta\mathbf{H}\mathbf{Y}$$` using, perhaps, a statistic like `$$MS= trace \left((\Delta\mathbf{H}\mathbf{Y})^T \Delta\mathbf{H}\mathbf{Y}\right)/ (k_{alternative} - k_{null})$$` Note that `\(trace\)` means the sum of diagonal elements. We will learn more about statistical methods, but nearly each method will be either a comparison of projections or a decomposition of a matrix obtained from a projection. --- .center[ # Part VI. The General Linear Model ## Only briefly ] --- .center[ # Part VI. The General Linear Model ## Only briefly ] Let `\(\mathbf{Z}\)` be a transformation of the data matrix, `\(\mathbf{Y}\)`. (Could be mean-centering of the data, scaling of the data, or normalizing the data.) Let `\(\mathbf{H} = \mathbf{X} \left( \mathbf{X}^T\mathbf{X} \right)^{-1}\mathbf{X}^T\)` be a **projection** matrix for a linear model design, `\(\mathbf{X}\)`. This matrix projects (predicts) what values of `\(\mathbf{Z}\)` should be, based on `\(\mathbf{X}\)`. Then, `\(\mathbf{HZ}\)` are the estimates a linear model makes! --- .center[ # Part VII. Matrix Decompositions ## Also only briefly ] --- ### Matrix Decompositions + Also called *matrix factorization* in some treatments. + There are many (but we will not consider them all) + Operations to decompose a matrix into multiple matrices with unique properties + Decompositions might be used as steps to solve an algebraic problem, or might be the solution to an algebraic problem + **The following are a few important decompositions** --- ### Matrix Decompositions ### Eigen decomposition For an `\(p \times p\)` square matrix, `\(\mathbf{S}\)`, there exists a vector of `\(p\)` values, `\(\mathbf{\lambda}\)`, such that, `$$\mathbf{Sv} =\lambda\mathbf{v}$$` which can also be written (in a more computationally feasible way) as `$$\left( \mathbf{S}-\lambda\mathbf{I} \right) \mathbf{v} = 0$$` + Solving each *possible* eigenvector ( `\(\mathbf{v}\)` ) and eigenvalue ( `\(\lambda\)` ) via a characteristic polynomial is called, eigen decomposition or "eigenanalysis". + Systematically find the `\(p\)` solutions of `\(\lambda\)` for `\(\left| \mathbf{S}-\lambda\mathbf{I} \right|\)`. Up to `\(p\)` real numbers exist if `\(\mathbf{S}\)` is positive-definite. Some of these numbers might be 0 if `\(\mathbf{S}\)` is not full-rank. Some can be negative if `\(\mathbf{S}\)` is not symmetric. + By finding the values of `\(\lambda\)` that work, the vectors, `\(\mathbf{v}\)`, are estimated such that the solution above holds true. Computers can solve characteristic polynomials easily via numerical approximation algorithms. A popular choice is the QR algorithm, which uses QR decomposition until convergence. Do not try to perform the characteristic polynomial solution or QR algorithm by hand! --- ### Matrix Decompositions (cont.) ### QR Decomposition `$$\mathbf{X}_{n \times k} =\mathbf{Q}_{n \times k} \mathbf{R}_{k \times k}$$` `\(\mathbf{Q}\)` and `\(\mathbf{R}\)` are solved via numerical approximation algorithms `\(\mathbf{R}\)` is a triangular matrix (only diagonal and elements above diagonal can have values other than 0). If `\(\mathbf{R}\)` is not full-rank, some of the upper elements will be 0. The matrix, `\(\mathbf{X}\)` was used on purpose, as QR decomposition is often used on design matrices in linear models, for calculations of coefficients or `\(SSCP\)` matrices. --- ### Singular Value decomposition of a Rectangular Matrix `$$\mathbf{Z}_{n \times p} =\mathbf{U}_{n \times p'} \mathbf{D}_{p' \times p'} \mathbf{V}_{p' \times p'}^T$$` `\(\mathbf{U}\)` (left singular values) and `\(\mathbf{V}\)` (right singular values) are solved via numerical approximation algorithms. The number of singular values is `\(p'= \min(n - 1,p)\)`. If `\(\mathbf{Y}\)` is not full-rank, some of the singular values in `\(\mathbf{D}\)` will be 0. The singular values decrease from largest to smallest. --- ### Singular Value decomposition of a Symmetric Matrix `$$\mathbf{S}_{p \times p} =\mathbf{V}_{p \times p} \mathbf{\Lambda}_{p \times p} \mathbf{V}_{p \times p}^T$$` `\(\mathbf{V}\)` (left and right singular values) are solved via numerical approximation algorithms. If `\(\mathbf{S}\)` is symmetric, positive-definite, `\(\mathbf{V}\)` is the matrix of eigenvectors and `\(\mathbf{\Lambda}\)` the matrix of eigenvalues, equal to those found from eigenanalysis using the characteristic polynomial function. This is a "short-cut" method, often employed for principal component analysis when it is sure that `\(\mathbf{S}\)` is symmetric, positive-definite. **It is also a short-cut for eigen decomposition.** If `\(\mathbf{S}\)` is a positive-definite, symmetric matrix, the singular values are the same as eigenvalues, hence the use of `\(\mathbf{\Lambda}\)` as nomenclature, rather than `\(\mathbf{D}\)`. The matching left and right singular vectors, `\(\mathbf{V}\)` are the same as eigenvectors. --- ### Cholesky decomposition of a Symmetric Matrix `$$\mathbf{S}_{p \times p} =\mathbf{T}_{p \times p} \mathbf{T}_{p \times p}^T$$` where `\(\mathbf{T}\)` is a triangular matrix. No need for a lot of detail here but just realize that Cholesky decompositions are sometimes used as solutions to find the *square-root of a symmetric matrix*. *Matrix square-roots are more complex than just this. The can also be found with singular value decomposition, and other methods.* --- ### Matrix Decompositions (Cont.) + Matrix decompositions provide mechanisms for many multivariate methods, including ordinations, linear model coefficient calculations, and matrix covariation analyses. Those presented here will show up again later. One can return to this lecture to review. --- Let's revisit this, focusing on the right sides: ### Linear Algebra: Single Goal (maybe not realized in this lecture, but in this week) Understand why these equations are foundational, and used in geometric morphometrics (GM) and multivariate statistics: `$$\mathbf{AZ}=\mathbf{UDV}^T$$` or `$$\mathbf{ZA}=\mathbf{VDU}^T$$` Let the matrix, `\(\mathbf{U}\)`, be a set of orthogonal vectors, related to `\(\mathbf{A}\)`, found via SVD. Let the matrix, `\(\mathbf{V}\)`, be a set of orthogonal vectors, related to `\(\mathbf{Z}\)`, found via SVD. Let the matrix, `\(\mathbf{D}\)`, be a diagonal matrix with scalars describing the strength of association between `\(\mathbf{U}\)` and `\(\mathbf{V}\)` (and thus, `\(\mathbf{A}\)` and `\(\mathbf{Z}\)`), found via SVD. **So, for example, the left-hand side might be estimated values from a liner model and the right hand side might find the principal components of the estimated values.** --- Let's revisit this, focusing on the right sides: ### Linear Algebra: Single Goal (maybe not realized in this lecture, but in this week) Understand why these equations are foundational, and used in geometric morphometrics (GM) and multivariate statistics: `$$\mathbf{AZ}=\mathbf{UDV}^T$$` or `$$\mathbf{ZA}=\mathbf{VDU}^T$$` Let the matrix, `\(\mathbf{U}\)`, be a set of orthogonal vectors, related to `\(\mathbf{A}\)`, found via SVD. Let the matrix, `\(\mathbf{V}\)`, be a set of orthogonal vectors, related to `\(\mathbf{Z}\)`, found via SVD. Let the matrix, `\(\mathbf{D}\)`, be a diagonal matrix with scalars describing the strength of association between `\(\mathbf{U}\)` and `\(\mathbf{V}\)` (and thus, `\(\mathbf{A}\)` and `\(\mathbf{Z}\)`), found via SVD. **So, for example, the left-hand side might be alignment of two separate data matrices and the right hand side might describe the strength of association (via the singular values, `\(\mathbf{D}\)`).** --- ### Summary/Important Points + Vectors are important. They are the points in high-dimensional data spaces. They are the algebraic work horses. They have length and direction. They are the most important component of high-dimensional data analysis. + Matrices are not really more special than vectors, but they bundle vectors in important ways. + Many multivariate analyses involve inverting matrices, which means obtaining a symmetric matrix to invert. It is imperative to understand how and why these steps are important. + Once comfortable with the process of inverting matrices, recognizing that EVERY multivariate analysis is a permutation of linear model calculations. Linear model calculations make multivariate analysis easier. + Understanding projection is also important. Recognizing that different methods are just different ways of projecting data will help simplify some complex topics. --- ### Revisit Goal Understand why these equations are foundational, and used in geometric morphometrics (GM) and multivariate statistics: `$$\mathbf{AZ}=\mathbf{UDV}^T$$` or `$$\mathbf{ZA}=\mathbf{VDU}^T$$` The left describes either multivariate data (like GM data) or an alignment of the multivariate data, say to a linear model or other set of data. The right describes a process for finding major vectors in a multivariate space, plus interpreting how important the vectors are. **Everything else we do is some permutation of this!**