xponential said:
Homework Statement
if M is dxd positive semi-definite matrix, then:
M = L'(or L transpose) * L
where L is a matrix of dimensions rxd
How can I get L given M and the number of rows, r?
for example,
if M is 800x800 p.s.d, I want L such that it is a 30x800 matrix and so
M (800x800) = L' (800x30) * L (30x800)
If you know a MATLAB function that can do just that, it would be perfect!
If M is an nxn positive semi-definite matrix, then one can show that there exist an upper-triangular matrix U such that M = U
T U. One obtains U by the algorithm for Cholesky factorization. Although many discussions of that algorithm assume M is positive-definite, if one is careful to put the argument in the correct form one can prove that it also applies to the positive-semidefnite case, with very slight modifications. I don't have access to Matlab, but the latest versions of Maple can do Cholesky factorization for positive-semidefinite matrices. Basically, the algorithm for a psd matrix M = (m
ij) is:
(1) If m
11 = 0, all elements in row 1 and column 1 must also be zero (otherwise, M could not be psd). If that holds, ignore row and column 1 and start again. Therefore, we may as well assume m
11 > 0. Set
u_{11}= \sqrt{m_{11}} \text{ and } u_{1j} = \frac{m_{1j}}{u_{11}},\; j=2,\ldots, n
(2) For i ≥ 2: let
u_{ii}^2 = m_{ii} - \sum_{p=1}^{i-1} u_{pi}^2. If ##u_{ii}^2 = 0## we must have ##m_{ij} - \sum_{p=1}^{i-1} u_{pi} u_{pj} = 0 \; \forall j > i.## In that case, set u
ij = 0 for all j ≥ i. Otherwise, set
u_{ii} = \sqrt{u_{ii}^2}, \; u_{ij} = \frac{m_{ij} - <br />
\sum_{p=1}^{i-1} u_{pi} u_{pj}}{u_{ii}} \; \forall j > i.
So, in the psd case some rows of U will be all zeros (unlike in the pd case); we can omit the zero rows of U to get a new matrix V which is no longer upper-triangular, but still gives
M = V
T V.
This algorithm is so simple that I used to have it programmed into an HP hand-held calculator; it could decompose matrices up to about 10x10.
RGV
