Hello, In using kmeans2 with the "random" initialization method, I often get: numpy.linalg.linalg.LinAlgError: Matrix is not positive definite - Cholesky decomposition cannot be computed I haven't dug too deeply into the implementation; I'm wondering if there's a gotcha someone can point out quickly here. What exactly is the "random" initialization method doing with the data set, and what constraints must the data have so as to keep whatever intermediate matrix it computes to be positive definite? Thanks.
On Mon, Mar 31, 2008 at 1:19 PM, Alec Solway <asolway@seas.upenn.edu> wrote:
Hello,
In using kmeans2 with the "random" initialization method, I often get:
numpy.linalg.linalg.LinAlgError: Matrix is not positive definite - Cholesky decomposition cannot be computed
I haven't dug too deeply into the implementation; I'm wondering if there's a gotcha someone can point out quickly here. What exactly is the "random" initialization method doing with the data set, and what constraints must the data have so as to keep whatever intermediate matrix it computes to be positive definite?
The relevant function is scipy/cluster/vq.py:_krandinit(). It is finding the covariance matrix and manually doing a multivariate normal sampling. Your data is most likely degenerate and not of full rank. It's arguable whether or not this should fail, but numpy.random.multivariate_normal() uses the SVD instead of a Cholesky decomposition to find the matrix square root, so it sort of ignores non-positive definiteness. Try replacing the last 2 lines of _krandinit() with x = N.random.multivariate_normal(mu, cov, k) and see if that helps you. -- Robert Kern "I have come to believe that the whole world is an enigma, a harmless enigma that is made terrible by our own mad attempt to interpret it as though it had an underlying truth." -- Umberto Eco
On 31-Mar-08, at 11:17 PM, Robert Kern wrote:
The relevant function is scipy/cluster/vq.py:_krandinit(). It is finding the covariance matrix and manually doing a multivariate normal sampling. Your data is most likely degenerate and not of full rank. It's arguable whether or not this should fail, but numpy.random.multivariate_normal() uses the SVD instead of a Cholesky decomposition to find the matrix square root, so it sort of ignores non-positive definiteness.
This might not be relevant, depending on how the covariance is computed, but one 'gotcha' I've seen with numerical algorithms that assume positive-definiteness is that occasionally floating point oddities will induce (very slight) non-symmetry of the input matrix, and thus the algorithm will choke; it's easily solved by averaging the matrix with it's transpose (though there are probably more efficient ways). David
On Tue, Apr 1, 2008 at 3:41 AM, David Warde-Farley <dwf@cs.toronto.edu> wrote:
This might not be relevant, depending on how the covariance is computed, but one 'gotcha' I've seen with numerical algorithms that assume positive-definiteness is that occasionally floating point oddities will induce (very slight) non-symmetry of the input matrix, and thus the algorithm will choke; it's easily solved by averaging the matrix with it's transpose (though there are probably more efficient ways).
Ultimately, the covariance is computed like so: dot(X, X.T.conj()) Ultimately, it is up to the implementation of dot() as to whether or not it will ensure exact symmetry. At least on my MacBook Pro with the Accelerate.framework providing the BLAS (an ATLAS derivative), I could not generate a matrix that was not exactly symmetrical over many random inputs. I think that less-than-full-rank data is the more likely problem. -- Robert Kern "I have come to believe that the whole world is an enigma, a harmless enigma that is made terrible by our own mad attempt to interpret it as though it had an underlying truth." -- Umberto Eco
participants (3)
-
Alec Solway
-
David Warde-Farley
-
Robert Kern