Square Root of a Matrix

Discussion in 'MATLAB' started by Greg Heath, Jul 14, 2005.

  1. Greg Heath

    Greg Heath Guest

    To find the square root of the n X n matrix C, I use the eigenvalue
    decomposition

    [E L] = eig(C);

    to obtain

    sqrtC = E*sqrt(L)*E'.

    Now, I would like an analytic expression for sqrtC(n,c) when C
    is a correlation coefficient matrix with equal nondiagonal
    components c:

    For n = 2, C = [1 c; c 1]
    For n = 3, C = [1 c c; c 1 c; c c 1]

    etc.

    Any ideas?

    Greg
     
    Greg Heath, Jul 14, 2005
    #1
    1. Advertisements

  2. Greg Heath

    Daniel Ennis Guest

    This isn't really the answer to your question, but just in case you
    didn't come across this you can also use,

    sqrtm

    -Dan-
     
    Daniel Ennis, Jul 14, 2005
    #2
    1. Advertisements

  3. ----------------------
    Hello Greg,

    The formula below should give you a square root of your matrix but not
    always the one given by 'sqrtm'. I have a strong suspicion 'sqrtm' is not
    consistent in its choice of square root, at least with respect to your
    kind of problem. Your eigenvalues can easily be negative so their square
    roots would then be imaginary. If your matrix is n x n, it has in general
    2^n square roots, obtained by all possible combinations of signs on the
    square roots of your n eigenvalues, but I don't know formulas like the one
    below for them all, only that one. You can get one more by reversing all
    the signs, but I'd have to put more work into the problem to find formulas
    for them all. It took some sweat just to find this one!

    x = repmat(sqrt((n-1)*c+1)/n-sqrt(c-1)/n*i,n,n);
    x = x + diag(repmat(sqrt(c-1)*i,n,1));

    (Remove "xyzzy" and ".invalid" to send me email.)
    Roger Stafford
     
    Roger Stafford, Jul 14, 2005
    #3
  4. ---------------------
    Hello Greg,

    I forgot that you stated c was a correlation coefficient, which
    therefore must necessarily be less than or equal to 1, so the supposedly
    imaginary part in the sqrt(c-1)*i term becomes real. The formula should
    then read:

    x = repmat(sqrt((n-1)*c+1)/n-sqrt(1-c)/n,n,n) + ...
    diag(repmat(sqrt(1-c),n,1));

    where x is a square root of your matrix. This agrees with 'sqrtm' as far
    as I can observe. With c <= 1 both are choosing the positive square roots
    of all the original eigenvalues. There is one eigenvalue in them equal to
    +sqrt((n-1)*c+1) and n-1 eigenvalues all equal to +sqrt(1-c).

    (Remove "xyzzy" and ".invalid" to send me email.)
    Roger Stafford
     
    Roger Stafford, Jul 14, 2005
    #4
  5. Greg Heath

    Greg Heath Guest

    Thanks.

    How did you figure it out? I quit at the cubic eigenvalue equation
    for n = 3. I would have tried MAPLE except I forgot how.

    Greg
     
    Greg Heath, Jul 14, 2005
    #5
  6. -----------------
    Hello Greg,

    I have to confess it was a case of blind luck. As I mentioned earlier,
    I had forgotten that your c <= 1 and evaluated the square root matrices
    for some particular integers c > 1 and for several values of n. I
    discovered that the real and imaginary parts of their elements were each
    square roots of something rational and thereby soon discovered the
    formulas in terms of c and n (i.e., sqrt((n-1)*c+1) and sqrt(1-c).) These
    formulas carry over to the case with c <= 1 except that there is no
    imaginary part since sqrt(1-c) becomes real - the two square roots are
    combined in real sums and differences. If I had tried things originally
    with only c <= 1, it is doubtful I could have disentangled the two square
    roots to find their formulas without a lot more work. There's fortuity
    for you!

    (Remove "xyzzy" and ".invalid" to send me email.)
    Roger Stafford
     
    Roger Stafford, Jul 14, 2005
    #6
  7. Greg Heath

    Greg Heath Guest

    I'll say. Did it ever occur to you to try MAPLE?

    Greg
     
    Greg Heath, Jul 15, 2005
    #7
  8. Hello Greg,

    I thought I'd let you know that, thanks to your matrix example, I
    managed to correct a mistaken notion I had previously held about the
    number of possible square roots an n x n matrix can have. When an
    eigenvalue is repeated, as occurs for your matrices for n greater than 2,
    there are infinitely many possible sets of corresponding eigenvectors.
    You can get them by performing arbitrary rotations (unitary
    transformations) on the space of such eigenvectors. Hence, by assigning
    both possible signs to the square roots of these repeated eigenvalues,
    arbitrarily rotating the corresponding eigenvector space, and then working
    your way back to a square root matrix, you can get infinitely many
    different square roots! That was a surprise. Your matrices have
    infinitely many square roots. I used TMW's Symbolic Math Toolbox (Maple)
    to verify that this is true with your 3 x 3 matrix. I do think the 2^n
    limit still applies to any matrix with non-repeating eigenvalues, however.

    (Remove "xyzzy" and ".invalid" to send me email.)
    Roger Stafford
     
    Roger Stafford, Jul 15, 2005
    #8
    1. Advertisements

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments (here). After that, you can post your question and our members will help you out.