# Closed-form expression for differential of matrix function

Let $X$ be a real $n\times n$ positive semidefinite matrix of rank $m\le n$ and let $Y\in\mathbb{R}^{m\times n}$ be the unique matrix satisfying (i) $X=Y^\top Y$, and (ii) $Y\, [I\, |\, 0]^\top = L$ with $L\in\mathbb{R}^{m\times m}$ being upper triangular with positive diagonal entries. (Notice that when $n=m$, $Y$ coincides with the (unique) Cholesky factor of $X$.) Let $A\in\mathbb{R}^{m\times m}$ be a positive definite matrix and $B\in\mathbb{R}^{n\times m}$ be such that $YB$ is nonsingular. Consider the following matrix-valued function $$f(X)=(YB)^{-1}A (YB)^{-\top}.$$

I'm looking for a closed-form expression of the matrix differential $\mathrm{d}\, f(X)$.

My attempt. A simple observation is that, in case $A=\alpha I$, $\alpha\in \mathbb{R}$, the sought differential reduces to $$\mathrm{d}\, f(X) = -\alpha(B^\top XB)^{-1} {B^{\top} \mathrm{d}\, X\, B} (B^\top XB)^{-1}.$$

For the case of general positive definite $A$, I describe below my attempt that deals with vectorized ($\mathrm{vec}$ operation) differentials. First, using chain rule, $$\tag{1}\label{eq:1} \mathrm{vec}(\mathrm{d}\, f(X)) = \frac{\mathrm{vec}(\mathrm{d}\, f(X))}{\mathrm{vec}(\mathrm{d}\, Y)}\frac{\mathrm{vec}(\mathrm{d}\, Y)}{\mathrm{vec}(\mathrm{d}\, X)} \mathrm{vec}(\mathrm{d}\, X).$$ Concerning the first term, we have $$\frac{\mathrm{vec}(\mathrm{d}\, f(X))}{\mathrm{vec}(\mathrm{d}\, Y)}=-\left((YB)^{-1}A(YB)^{-\top}B^\top \otimes (YB)^{-1}\right)-\left((YB)^{-1}\otimes (YB)^{-1}A(YB)^{-\top}B^\top \right)K^{n,m},$$ where $K^{n,m}$ is an $nm\times nm$ commutation matrix and $\otimes$ denotes Kronecker product. For the second term, using the same argument of the accepted answer to this MO question, we have $$\frac{\mathrm{vec}(\mathrm{d}\, Y)}{\mathrm{vec}(\mathrm{d}\, X)} = \left((Y^\top\otimes I)K^{n,m}+(I\otimes Y^\top)\right)^{-L},$$ where $(\cdot)^{-L}$ denotes left-inversion. Plugging the above-derived expressions into \eqref{eq:1}, yields a (vectorized) formula for the sought differential.

At this point, provided that my calculations are correct, I wonder whether Eq. \eqref{eq:1} can be further simplified and written in matrix form (if possible). I conjecture (or, more frankly, I hope) that the final expression is similar to the scalar ($A=\alpha I$) case, i.e., something of the form $$\mathrm{d}\, f(X) = -(YB)^{-1}A (YB)^{-\top} {B^{\top} \mathrm{d}\, X\, B} (YB)^{-1}A (YB)^{-\top}.$$

Edit. I believe that the hard part of my question regards the computation of the differential of $Y$ w.r.t. $X$. Of course, this can be accomplished via vectorization, as described above, and subsequent "matricization". Nevertheless, I wonder whether a more "genuine" matrix expression for $\mathrm{d}Y$ exists.

Any comment/suggestion is very appreciated.

• Since $Y$ is not well-defined, what exactly do you mean by its differential? – Igor Rivin Dec 28 '17 at 18:35
• @IgorRivin: According to the accepted answer to this MO question mathoverflow.net/questions/150427/…, I think that the (classical notion of) differential of $Y$ w.r.t. $X$ is well-defined under certain assumption concerning invertibility of $(I\otimes Y^\top)$. See also my edit. – Ludwig Dec 30 '17 at 17:35
• Since multiplying $Y$ on the left by an orthogonal matrix gives you another $Y,$ I still don't know what you mean... – Igor Rivin Dec 30 '17 at 17:54
• @IgorRivin: I edited the OP in order to deal with a well-defined $Y$. – Ludwig Dec 31 '17 at 9:04

It suffices to consider the case $n=2$, $m=1$. Namely, write $Y=[Y_1|Y_2]$ etc, then $L=Y_1$ is upper triangular with positives on the diagonal, and \begin{align*} X &= Y^\top Y = \begin{pmatrix} Y_1^\top Y_1 & Y_1^\top Y_2 \\ Y_2^\top Y_1 & Y_2^\top Y_2\end{pmatrix}\,, \\ dX &= \begin{pmatrix} (dY_1)^\top. Y_1 + Y_1^\top.dY_1 & (dY_1)^\top. Y_2 + Y_1^\top.dY_2 \\ (dY_2)^\top. Y_1 + Y_2^\top. dY_1 & (dY_2)^\top. Y_2 + Y_2^\top.dY_2 \end{pmatrix} \\ YB &= Y_1B_1 + Y_2B_2 \end{align*} Moreover, \begin{align*} df(X) =& -(YB)^{-1}.dY.B.(YB)^{-\top} - (YB)^{-1}.A.(YB)^{-\top}.(dY.B)^\top.(YB)^{-\top} \\=& -(YB)^{-1}.dY.B.f(X) - f(X)((YB)^{-1}.dY.B)^\top \\ dY.B =& dY_1.B_1 + dY_2.B_2 \end{align*}
Now you can compute $dY_1$ from $dX_{1,1}$ as in the invertible situation. Let us drop subindices and go to the invertible situation. $C\mapsto C^{-\top}$ is the Cartan involution on the reductive Lie group $GL^+(m)$. Consider the Iwasawa decomposition $GL^+(m) = SO(m).A.N$, $A$ the diagonal matrices with positive entries, and $N$ the upper unipotent matrices (which equals here the Gram-Schmidt orthonormalisation procedure with the coefficients arranged in $A.N$). First note that $Y:S_+(m)\to AN$ is a smooth map into a Lie group, so $dY$, better $TY: TS_+(m)\to T(AN)$, and the right logarithmic derivative $\delta Y:= TY.Y^{-1}:TS_+(m)\to \mathfrak{an}$ is a Lie algebra valued 1-form, $\delta Y\in \Omega^1(S_+(m);\mathfrak{an})$. You can get back $Y$ from the 1-form $\delta Y$ by Cartan development. Namely, $\delta Y$ describes a flat principal connection on the trivial principal $AN$-bundle $S_+(m)\times AN \to S_+(m)$, and any horizontal leaf of it is a right translate of the mapping $Y$. Moreover $Z\mapsto Z^\top$ restricts to the mapping $\mathfrak{an}\to \mathfrak{an}^*$ corresponding to the inner product $\operatorname{Trace}(U^\top.V)$ on $\mathfrak{gl}(m)$. We have \begin{gather*} (dY)^\top.Y + Y^{\top}.dY = dX = dX^{\top} \\ \delta Y + (\delta Y)^\top = dY.Y^{-1} + Y^{-\top}.(dY)^\top = Y^{-\top}.dX.Y^{-1} = (dX.Y^{-1})^{\top}.Y^{-1} \end{gather*} Now, $\delta Y + (\delta Y)^\top$ allows in a simple way to compute $\delta Y$ (take the upper triangular part and 1/2 of the diagonal entries).
• Thanks for the answer! Could you please elaborate a little more how $\mathrm{d} Y_1$ can be computed from $\mathrm{d} X_{1,1}$ (without using vectorization, if possible)? – Ludwig Jan 1 '18 at 9:40