This is my first post, please let me know if I am not being clear.

I am trying to understand the sklearn.mixture.GaussianMixture.score(X).

As I understand that the following is the equation for log likelihood:

In my case, I want to calculate log likelihood for a sklearn.mixture.GaussianMixture(covariances_type = 'diag'). Below is the code that is calculating .

def _compute_log_det_cholesky(matrix_chol, covariance_type, n_features):
    elif covariance_type == 'diag':
        log_det_chol = (np.sum(np.log(matrix_chol), axis=1))
    return log_det_chol

This post explains how to use Cholesky's precision matrix. My question is:

  • why are are not using sklearn.mixture.GaussianMixture.fit(X).covariances_.?

Since we are using sklearn.mixture.GaussianMixture(covariances_type = 'diag'), the output of sklearn.mixture.GaussianMixture.fit(X).covariances_ is k x m. Where k is number of components and m is the number of features.

Thanks in advance!


Working with Cholesky factors is cheaper because determinants of triangular matrices are given by the product of the elements on the diagonal. This massively reduces the computation time for the determinants because we effectiely get to skip the computations which are needed to compute the determinant in the general case.

Working with Cholesky factors also economizes the computations because we are also computing Cholesky factors when working with the $\Sigma^{-1}$ portion of the log-likelihood, so we're re-using the same intermediate steps instead of computing two different, expensive things (the inverse of a matrix and the determinant of a matrix).

The code is just collapsing some operations to avoid redundant computations. In the full code, the determinant for covariance_type=="diag" has the -0.5 factor combined in log_det.

def _estimate_log_gaussian_prob(X, means, precisions_chol, covariance_type):
    log_det = _compute_log_det_cholesky(
    precisions_chol, covariance_type, n_features)


    elif covariance_type == 'diag':
       precisions = precisions_chol ** 2
       log_prob = (np.sum((means ** 2 * precisions), 1) -
                  2. * np.dot(X, (means * precisions).T) +
                np.dot(X ** 2, precisions.T))

    return -.5 * (n_features * np.log(2 * np.pi) + log_prob) + log_det

Since the Gaussian log-likelihood doesn't ever use the covariance matrix $\Sigma$, but instead only uses its inverse $\Sigma^{-1}$ and a negative multiple of its determinant $-\frac{1}{2} \log | \Sigma|$, we actually only need to work with Cholesky factors of the precision matrix.

You can prove that this is correct using the properties of Cholesky factors, determinants, and matrix inverses. Our goal is to find a way to write $-\frac{1}{2} \ln |\Sigma|$ in terms of the Cholesky factor of the precision matrix: $\Sigma^{-1} = A^\top A$

$$\begin{align} \left| \Sigma^{-1} \right| &= \left| A^\top A \right| \\ &= \left| A \right| \left| A^{\top} \right| \\ &= \left| A \right|^2 \\ \ln \left| \Sigma^{-1} \right| &= \ln \left| A \right|^2 \\ \ln \left| \Sigma \right| &= \ln \left| A \right|^{-2} \\ &= -2 \ln \left| A \right| \\ -\frac{1}{2}\ln \left| \Sigma \right| &= \ln \left| A \right| \end{align}$$

This proves that the code is algebraically correct for all Cholesky factors of nonsingular $\Sigma$, including the special case when $\Sigma$ is diagonal.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.