[SciPy-Dev] matrix_sqrt for singular symmetric square matrices

josef.pktd at gmail.com josef.pktd at gmail.com
Sun Oct 28 10:02:07 EDT 2018


On Sun, Oct 28, 2018 at 9:11 AM Ilhan Polat <ilhanpolat at gmail.com> wrote:

> This is covered by LDLt decomposition with an extra step of taking the
> square root of each block in D. The economy mode of this would be removing
> rows/cols 0 blocks from D and L. For the second case I think a polar
> decomposition would be a better approach.
>
> Calling these factors a square root might take you out of the common
> terminology though
>

It's common terminology in statistics and econometrics with notation like
R = S^{1/2} or Q = S^{-1/2}

(In a related case where statsmodels uses cholesky for regression we have
boring but descriptive names like
cholsigma and cholsigmainv )

It would be fine if there is a "linalg" name that is discoverable.
The point would be to have a function for when cholesky doesn't work and
that I don't have to figure out each time I need it.
I had also written a similar function as `matrix_half` in the past and
inlined it several times, but it's the first time that I tried to get a
reduced number of rows or columns for it.

My usecase is similar to np.linalg.pinv that I use very often because it is
convenient and simple even if working with the SVD directly would be
computationally more efficient for getting additional results.

Josef


>
> On Sun, Oct 28, 2018 at 9:00 AM Gael Varoquaux <
> gael.varoquaux at normalesup.org> wrote:
>
>> On Sun, Oct 28, 2018 at 08:56:37AM +0100, Gael Varoquaux wrote:
>> > using '@' to denote the matrix product, and S = np.diag(s) is the
>> > diagonal matrix of eigenvalues. The matrix square root is then given by:
>>
>> >    sqrt(M) = U' @ np.diag(np.sqrt(s)) @ U
>>
>> I forgot to say: this is the definition used by Joseph, in his original
>> post (so basically, I am backing his choice), with the only difference
>> that I would not use an SVD, but and "eigh", which should be faster and
>> more stable for SPD matrices (and non SPD matrices do not have a square
>> root).
>>
>> G
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at python.org
>> https://mail.python.org/mailman/listinfo/scipy-dev
>>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at python.org
> https://mail.python.org/mailman/listinfo/scipy-dev
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20181028/8b5d4813/attachment.html>


More information about the SciPy-Dev mailing list