The acobian matrix is 物-加 Breakdown occurs if p or s is 0 Example 7: Traceless Symmetric= Polar Decomposition(S=QAQ, trS=0) The reader will recall the usual definition of polar coordinates. If (p, s) are Cartesian coordinates, then the angle is 0= arctan(s/p) and the radius is r= vp2+s2. If we take a symmetric traceless 22 matrix and compute its eigenvalue and eigenvector decomposition, we find that the eigendecomposition is mathe- matically equivalent to the familiar transformation between Cartesian and polar coordinates. Indeed one of the eigenvectors of S is(cos o, sin o), where =0 /2. The Jacobian matrix is The Jacobian is the inverse of the radius. This corresponds to the familiar formula using the more usual notation drdy /r= drdo so that det =1/r. Breakdown occurs when r=0 Example 8: The Symmetric Eigenvalue Problem(S=QAQ) We compute the Jacobian for the general symmetric eigenvalue and eigenvector decomposition. Let ]-[=m1[2]如 We can compute the eigenvectors and eigenvalues of s directly in MATLAB and compute the Jacobian of the two eigenvalues and the eigenvector angle, but when we tried with the maple toolbox we found that the symbolic toolbox did not handle arctan"very well. Instead we found it easy to compute the Jacobian in the other direction We write S=QAQ, where Q is 2 x 2 orthogonal and A is diagonal. The Jacobian is 2 sin e cos 8(A1-A2) sin e cos 0(A1-A2) (cos20-sin20)(A1-A2) sin 0 cos 0 -sin 0 cos 0 det =A1-A2 Breakdown occurs if S is a multiple of the identity Example 9: Symmetric Congruence(Y=ASA
p " # " # The Jacobian matrix is 1 √p 0 0 1 √ s J = 0 0 1 2 √ r/p r/s 2 ps √ps √ − − ps so that 1 det J = 4ps Breakdown occurs if p or s is 0. Example 7: Traceless Symmetric = Polar Decomposition (S = QΛQT , tr S = 0) The reader will recall the usual definition of polar coordinates. If (p, s) are Cartesian coordinates, then the angle is θ = arctan(s/p) and the radius is r = p2 + s2 . If we take a symmetric traceless 2 × 2 matrix p s S = , s −p and compute its eigenvalue and eigenvector decomposition, we find that the eigendecomposition is mathematically equivalent to the familiar transformation between Cartesian and polar coordinates. Indeed one of the eigenvectors of S is (cos φ, sin φ), where φ = θ/2. The Jacobian matrix is √p2 p +s √p2 s +s 2 2 J = −s p p2+ 2 s2 p2+s The Jacobian is the inverse of the radius. This corresponds to the familiar formula using the more usual notation dxdy/r = drdθ so that det J = 1/r. Breakdown occurs when r = 0. Example 8: The Symmetric Eigenvalue Problem (S = QΛQT ) We compute the Jacobian for the general symmetric eigenvalue and eigenvector decomposition. Let p s cos θ − sin θ λ1 cos θ − sin θ T S = = . sin θ cos θ λ2 sin θ cos θ s r We can compute the eigenvectors and eigenvalues of S directly in MATLAB and compute the Jacobian of the two eigenvalues and the eigenvector angle, but when we tried with the Maple toolbox we found that the symbolic toolbox did not handle “arctan” very well. Instead we found it easy to compute the Jacobian in the other direction. We write S = QΛQT , where Q is 2 × 2 orthogonal and Λ is diagonal. The Jacobian is −2 sin θ cos θ (λ1 − λ2) cos2 θ sin 2 θ J = 2 sin θ cos θ (λ 1 − λ2) sin 2 θ cos2 θ (cos2 θ − sin 2 θ)(λ1 − λ2) sin θ cos θ − sin θ cos θ so that det J = λ1 − λ2 . Breakdown occurs if S is a multiple of the identity. Example 9: Symmetric Congruence (Y = ATSA)
Let Y= ASA, where Y and S are symmetric, but A b c d is arbitrary. The Jacobian matrix is J=b2 d22db cb+ and det =(det A)3 The cube on the determinant tends to surprise many people. Can you guess what it is for an n x n symmetric matrix(Y=ASA)? The answer(det =( det A)n+l)is in Example 3 of Section 3 7jacobian2by2m 7Code 8.1 of Random Eigenvalues by Alan Edelman EXperiment Compute the Jacobian of a 2x2 matrix function Symbolic tools are not perfect. The author exercised care in choosing the variables Lp ] A=La b;c d Compute Jacobians J=jacobian(Y(:),X(: )) JAC_square =factor(det(J)) Y=X-3; J=jacobian((: ), X( )) JAC_cube -factor(det(J)) Y=inv(X) J=jacobian(Y (:),X(: )) JAc_inv factor(det(s)) Y=A*X J=jacobian(Y(:), X(: )),JAC-linear = factor(det(J)) Y=Lp q;r/p det(X)/p]: J=jacobian(Y (:), X(: )) JAC-l1 -factor(det (j)) =Lp s r];y=Sqrt(p) sqrt(s) r/(sqrt(p)*sqrt(s))] J=jacobian(y, x) JAC DMD -factor(det(j)) x=Lp s]: y=[ sqrt(p"2+s"2) atan(s/p)] J=jacobian(y, x) JAC- notrace =factor(det()) Q=Icos(t) -sin(t): sin(t) cos(t)] D=[e10;0e2] y=[Y(1,1)Y(2,2)Y(1,2)];x[te1e2] J=jacobian(y, x) JAC-symeig =simplify(det (J)) X=[ Lp s; s r] Y=A,3*X*A y=[Y(1,1)Y(2,2)Y(1,2)];x=[prs] J=jacobian (y, x) JAC-symcong =factor(det()) Code 1
� � � � Let Y = ATSA, where Y and S are symmetric, but A = a b is arbitrary. The Jacobian matrix is c d 2 2 a c 2 ca J = b2 d2 2 db ab cd cb + ad and det J = (det A)3 . The cube on the determinant tends to surprise many people. Can you guess what it is for an n × n symmetric matrix (Y = ATSA)? The answer (det J = (det A)n+1)is in Example 3 of Section 3. %jacobian2by2.m %Code 8.1 of Random Eigenvalues by Alan Edelman %Experiment: Compute the Jacobian of a 2x2 matrix function %Comment: Symbolic tools are not perfect. The author % exercised care in choosing the variables. syms p q r s a b c d t e1 e2 X=[p q ; r s]; A=[a b;c d]; %% Compute Jacobians Y=X^2; J=jacobian(Y(:),X(:)), JAC_square =factor(det(J)) Y=X^3; J=jacobian(Y(:),X(:)), JAC_cube =factor(det(J)) Y=inv(X); J=jacobian(Y(:),X(:)), JAC_inv =factor(det(J)) Y=A*X; J=jacobian(Y(:),X(:)), JAC_linear =factor(det(J)) Y=[p q;r/p det(X)/p]; J=jacobian(Y(:),X(:)), JAC_lu =factor(det(J)) x=[p s r];y=[sqrt(p) sqrt(s) r/(sqrt(p)*sqrt(s))]; J=jacobian(y,x), JAC_DMD =factor(det(J)) x=[p s]; y=[ sqrt(p^2+s^2) atan(s/p)]; J=jacobian(y,x), JAC_notrace =factor(det(J)) Q=[cos(t) -sin(t); sin(t) cos(t)]; D=[e1 0;0 e2];Y=Q*D*Q.’; y=[Y(1,1) Y(2,2) Y(1,2)]; x=[t e1 e2]; J=jacobian(y,x), JAC_symeig =simplify(det(J)) X=[p s;s r]; Y=A.’*X*A; y=[Y(1,1) Y(2,2) Y(1,2)]; x=[p r s]; J=jacobian(y,x), JAC_symcong =factor(det(J)) Code 1
3 Y=BXA and the kronecker product 3.1 Jacobian of y= BXA(Kronecker Product Approach) There is a"nuts and bolts "approach to calculate some Jacobian determinants. A good example function the matrix inverse Y=X-. We recall from Example 4 of Section 1 that dY=-X-dXX-I In words, the perturbation dX is multiplied on the left and on the right by a fixed matrix. When this happens we are in a"Kronecker Product"situation, and can instantly write down the Jacobian We provide two definitions of the Kronecker Product for square matrices A E R", to BE R See 1 for a nice discussion of Kronecker products Operator Definition A⑧ B is the operator from X∈ Rm,n to y∈Rm, n where y=BXA. We write (A⑧B)X=BXA Matrix Definition(Tensor Product) A⑧ B is the matrix A8B= The following theorem is important for applications Theorem 1. det(A o B)=(det a)m(det B)" Application: If Y=X-I then dy=-(X-ToX-)dX so that I detJ=det XI Notational Note: The correspondence between the operator definition and matrix definition is worth spelling out. It corresponds to the following identity in MATLAB B Y(:)= kron(A, B)* X( i The second line does not change Y Here kron(A, B)is exactly the matrix in Equation(4), and X(: is the column vector consisting of the columns of X stacked on top of each other. (In computer science this is known as storing an array in "column major"order. Many authors write vec(X), where we use X(: ) Concretely, we have that vec(BXA)=(A@ B)vec(X) ⑧ B is as in(4) Proofs of theorem 8.1: Assume A and B are diagonalizable, with Au;= A;;(i=l,..., n) and Bu;=4;vi(i=l,..., m) Let Miy=v;u;. The mn matrices Mi, form a basis for Rmn and they are eigenmatrices of our map since BMij AT=HiA, Mig. The determinant is (det A)"(det B) The assumption of diagonalizability is not important
Y 3 Y = BXAT and the Kronecker Product 3.1 Jacobian of Y = BXAT (Kronecker Product Approach) There is a “nuts and bolts” approach to calculate some Jacobian determinants. A good example function is the matrix inverse Y = X−1 . We recall from Example 4 of Section 1 that dY = −X−1dXX−1 . In words, the perturbation dX is multiplied on the left and on the right by a fixed matrix. When this happens we are in a “Kronecker Product” situation, and can instantly write down the Jacobian. We provide two definitions of the Kronecker Product for square matrices A ∈ Rn,n to B ∈ Rm,m . See [1] for a nice discussion of Kronecker products. Operator Definition A ⊗ B is the operator from X ∈ Rm,n to Y ∈ Rm,n where Y = BXAT . We write (A ⊗ B)X = BXAT . Matrix Definition (Tensor Product) A ⊗ B is the matrix a11B . . . a1m2 B A ⊗ B = . . . . . . . (4) am11B . . . am1m2 B The following theorem is important for applications. Theorem 1. det(A ⊗ B) = (det A)m(det B)n Application: If Y = X−1 then dY = −(X−T ⊗ X−1) dX so that det J| = det X −2n | | | . Notational Note: The correspondence between the operator definition and matrix definition is worth spelling out. It corresponds to the following identity in MATLAB Y = B * X * A’ Y(:) = kron(A,B) * X(:) % The second line does not change Y Here kron(A,B) is exactly the matrix in Equation (4), and X(:) is the column vector consisting of the columns of X stacked on top of each other. (In computer science this is known as storing an array in “column major” order.) Many authors write vec(X), where we use X(:). Concretely, we have that vec(BXAT) = (A ⊗ B)vec(X) where A ⊗ B is as in (4). Proofs of Theorem 8.1: Assume A and B are diagonalizable, with Aui = λiui (i = 1, . . . , n) and Bvi = µivi (i = 1, . . . , m). Let Mij = T . The mn matrices Mij form a basis for Rmn v and they are eigenmatrices of our map since iuj BMijAT = µiλjMij . The determinant is µiλj or (det A) m(det B) n . (5) 1≤i≤n 1≤j≤m The assumption of diagonalizability is not important
Also one can directly manipulate the matrix since obviously A⑧B=(4⑧D)(I⑧B) as operators, and det I@B=(det B)" and det A o I which can be reshuffled into det I o A=(det A)m Other proofs may be obtained by working with the "LU"decomposition of A and B, the Svd of A and B, the QR decomposition, and many others Mathematical Note: That the operator Y= BXA can be expressed as a matrix is a co B(CIX,+C2X2)A= C1BX1A+C2BX2A ie.(A⑧B)(c1X1+c2X2)=c1(A⑧B)X1+c2(A⑧B)X2 It is important to realize that a linear transformation from Rm,n to R,n is defined by an element of Rmn, mn, i.e., by the m2n2 entries of an mn x mn matrix. The transformation defined by Kronecker products is an m-+ of this m2n2 Some Kronecker product properties 1.(A8B)=A8B 2.(A⑧B)-1=A-18B-1 3.det(A⑧B)=det(4)mdet(B)n,A∈Rn,n,B∈Rm,m 4.tr(A⑧B)=tr(Atr(B) 5. A@ B is orthogonal if A and B are orthogonal 6.(A⑧B)C⑧D)=(AC)⑧(BD 7. If Au= u, and Bu=Hu, then if X=vu, then BXA=AuX, and also AX"B=AuX.Therefore A8 B and B e A have the same eigenvalues, and transposed eigenvectors. It is easy to see that property 5 holds, since if A and B are orthogonal (AoB)=A8B=AB (A⑧B)-1 Linear Subspace Kronecker Products Some researchers consider the symmetric Kronecker product @ sym. In fact it has become clear that one might consider an anti-symmetric, upper-triangular or even a Toeplitz Kronecker product. We formulate Definition: Let S denote a linear subspace of Rmn and Ts a projection onto S. If AER, and B Rm then we define(A⑧B)sX=丌s(BXA)forX∈S Comments 1. If S is the set of symmetric matrices then m= n and (A 8 B)smx-BXA+AXBT 2. If S is the set of anti-sym matrices, then m= n (A)ant x= BXA+AXBT as well but this matrix is anti-sym
Also one can directly manipulate the matrix since obviously A ⊗ B = (A ⊗ I)(I ⊗ B) as operators, and det I ⊗ B = (det B)n and det A ⊗ I which can be reshuffled into det I ⊗ A = (det A)m . Other proofs may be obtained by working with the “LU” decomposition of A and B, the SVD of A and B, the QR decomposition, and many others. Mathematical Note: That the operator Y = BXAT can be expressed as a matrix is a consequence of linearity: B(c1X1 + c2X2)AT = c1BX1AT + c2BX2AT i.e. (A ⊗ B)(c1X1 + c2X2) = c1(A ⊗ B)X1 + c2(A ⊗ B)X2 It is important to realize that a linear transformation from Rm,n to Rm,n is defined by an element of Rmn,mn , i.e., by the m2n2 entries of an mn×mn matrix. The transformation defined by Kronecker products is an m2 + n2 subspace of this m2n2 dimensional space. Some Kronecker product properties: 1. (A ⊗ B)T = AT ⊗ BT 2. (A ⊗ B)−1 = A−1 ⊗ B−1 , A ∈ Rn,n , B ∈ Rm,m 3. det(A ⊗ B) = det(A)m det(B)n 4. tr(A ⊗ B) = tr(A)tr(B) 5. A ⊗ B is orthogonal if A and B are orthogonal 6. (A ⊗ B)(C ⊗ D) = (AC) ⊗ (BD) 7. If Au = λu, and Bv = µv, then if X = vuT , then BXAT = λµX, and also AXTBT = λµXT . Therefore A ⊗ B and B ⊗ A have the same eigenvalues, and transposed eigenvectors. It is easy to see that property 5 holds, since if A and B are orthogonal (A ⊗ B)T = AT ⊗ BT = A−1 ⊗ B−1 = (A ⊗ B)−1 . Linear Subspace Kronecker Products Some researchers consider the symmetric Kronecker product ⊗sym. In fact it has become clear that one might consider an anti-symmetric, upper-triangular or even a Toeplitz Kronecker product. We formulate a general notion: Definition: Let S denote a linear subspace of Rmn and πS a projection onto S. If A ∈ Rn,n and B ∈ Rm,m then we define (A ⊗ B)SX = πS (BXAT) for X ∈ S. Comments 1. If S is the set of symmetric matrices then m = n and BXAT + AXBT (A ⊗ B)symX = 2 2. If S is the set of anti-sym matrices, then m = n BXAT + AXBT (A ⊗ B)antiX = as well 2 but this matrix is anti-sym
3. If S is the set of upper triangular matrices, then m=n and(A8 B)upper upper(BXA) 4. We recall that Ts is a projection if Ts is linear and for all X E Rmn, TSXES (A⑧B)sym、BXA+AXB If b= a then by considering eigenmatrices E=uiu; for i < j 6. Jacobian of up Special case: A lower triangular, B upper triangular so that (A 8 B)upper X= bXA since BXA is upper triangular The eigenvalues of A and B are A;=Aii and u,=Bii respectively, where Au;=Aiu; and Bui=AiU (i=l,., n). The matrix Mii= viuf for i <j is upper triangular since v; and u, are zero below the ith and above the jth component respectively. (The eigenvectors of a triangular matrix are triangular. Since the Miy are a basis for upper triangular matrices BMiA=AiAj Mij. We then have detJ=Ⅱ4入=(A14…)吗-1-2…,pn) (A142413….Amn)(B1B21B332..Bn) Note that J is singular if an only if A or B is. 7. Jacobian of e Toeplitz Let X be a Toeplitz matrix. We can defin B)X=Toeplitz(BXA) where Toeplitz averages every diagonal 4 Jacobians of Linear Functions Powers and inverses The Jacobian of a linear map is just the determinant. This determinant is not always easily computed. The dimension of the underlying space of matrices plays a role. For example the Jacobian of Y= 2X is 2" for XERXn, 22 for upper triangular or symmetric X, 22 for antisymmetric X, and 2n for diagonal We will concentrate on general real matrices X and explore the symmetric case and case as well when appropriate
Y Y 3. If S is the set of upper triangular matrices, then m = n and (A ⊗ B)upper = upper(BXAT). 4. We recall that πS is a projection if π is linear and for all X ∈ Rmn S , πSX ∈ S. 5. Jacobian of ⊗sym: BXAT + AXBT (A ⊗ B)symX = 2 If B = A then det J = λiλj = (det A) n+1 i≤j T by considering eigenmatrices E = uiuj for i ≤ j. 6. Jacobian of ⊗upper: Special case: A lower triangular, B upper triangular so that (A ⊗ B)upperX = BXAT since BXAT is upper triangular. The eigenvalues of A and B are λi = Aii and µj = Bjj respectively, where Aui = λiui and Bvi = µivi T (i = 1, . . . , n). The matrix Mij = viuj for i ≤ j is upper triangular since vi and uj are zero below the ith and above the jth component respectively. (The eigenvectors of a triangular matrix are triangular.) j i M for i ≤ j . ij = Since the Mij are a basis for upper triangular matrices BMijAT = µiλjMij . We then have n n 2λ3 det J = µiλj = (λ1λ2 3 . . . λn)(µ1 µ n−1 µ n−2 . . . µn) 2 3 i≤j n = (A11A2 33 . . . Ann)(Bn 22A3 11Bn−1Bn−2 . . . Bnn). 22 33 Note that J is singular if an only if A or B is. 7. Jacobian of ⊗Toeplitz Let X be a Toeplitz matrix. We can define (A ⊗Toeplitz B)X = Toeplitz(BXAT) where Toeplitz averages every diagonal. 4 Jacobians of Linear Functions, Powers and Inverses The Jacobian of a linear map is just the determinant. This determinant is not always easily computed. The 2 dimension of the underlying space of matrices plays a role. For example the Jacobian of Y = 2X is 2n for n(n+1) n(n−1) X ∈ Rn×n, 2 for upper triangular or symmetric 2 X, 2 for antisymmetric 2 X, and 2n for diagonal X. We will concentrate on general real matrices X and explore the symmetric case and triangular case as well when appropriate