natural language toolkit expectation maximization clusterer c 20012023 nltk project trevor cohn tacohncs mu oz au url https www nltk org for license information see license txt the gaussian em clusterer models the vectors as being produced by a mixture of k gaussian sources the parameters of these sources prior probability mean and covariance matrix are then found to maximise the likelihood of the given data this is done with the expectation maximisation algorithm it starts with k arbitrarily chosen means priors and covariance matrices it then calculates the membership probabilities for each vector in each of the clusters this is the e step the cluster parameters are then updated in the m step using the maximum likelihood estimate from the cluster membership probabilities this process continues until the likelihood of the data does not significantly increase creates an em clusterer with the given starting parameters convergence threshold and vector mangling parameters param initialmeans the means of the gaussian cluster centers type initialmeans seq of numpy array or seq of sparsearray param priors the prior probability for each cluster type priors numpy array or seq of float param covariancematrices the covariance matrix for each cluster type covariancematrices seq of numpy array param convthreshold maximum change in likelihood before deemed convergent type convthreshold int or float param bias variance bias used to ensure nonsingular covariance matrices type bias float param normalise should vectors be normalised to length 1 type normalise boolean param svddimensions number of dimensions to use in reducing vector dimensionsionality with svd type svddimensions int set the parameters to initial values do the e and m steps until the likelihood plateaus estep calculate hidden variables hi j mstep update parameters cvm p mean bias term to stop covariance matrix being singular calculate likelihood fixme may be broken check for convergence happens when the exponent is negative infinity i e b 0 i e the inverse of cvm is huge cvm is almost zero noninteractive demonstration of the clusterers with simple 2d data example from figure 14 10 page 519 manning and schutze classify a new vector show the classification probabilities natural language toolkit expectation maximization clusterer c 2001 2023 nltk project trevor cohn tacohn cs mu oz au url https www nltk org for license information see license txt the gaussian em clusterer models the vectors as being produced by a mixture of k gaussian sources the parameters of these sources prior probability mean and covariance matrix are then found to maximise the likelihood of the given data this is done with the expectation maximisation algorithm it starts with k arbitrarily chosen means priors and covariance matrices it then calculates the membership probabilities for each vector in each of the clusters this is the e step the cluster parameters are then updated in the m step using the maximum likelihood estimate from the cluster membership probabilities this process continues until the likelihood of the data does not significantly increase creates an em clusterer with the given starting parameters convergence threshold and vector mangling parameters param initial_means the means of the gaussian cluster centers type initial_means seq of numpy array or seq of sparsearray param priors the prior probability for each cluster type priors numpy array or seq of float param covariance_matrices the covariance matrix for each cluster type covariance_matrices seq of numpy array param conv_threshold maximum change in likelihood before deemed convergent type conv_threshold int or float param bias variance bias used to ensure non singular covariance matrices type bias float param normalise should vectors be normalised to length 1 type normalise boolean param svd_dimensions number of dimensions to use in reducing vector dimensionsionality with svd type svd_dimensions int set the parameters to initial values do the e and m steps until the likelihood plateaus e step calculate hidden variables h i j m step update parameters cvm p mean bias term to stop covariance matrix being singular calculate likelihood fixme may be broken check for convergence happens when the exponent is negative infinity i e b 0 i e the inverse of cvm is huge cvm is almost zero non interactive demonstration of the clusterers with simple 2 d data example from figure 14 10 page 519 manning and schutze classify a new vector show the classification probabilities
try:
    import numpy
except ImportError:
    pass

from nltk.cluster.util import VectorSpaceClusterer


class EMClusterer(VectorSpaceClusterer):
    

    def __init__(
        self,
        initial_means,
        priors=None,
        covariance_matrices=None,
        conv_threshold=1e-6,
        bias=0.1,
        normalise=False,
        svd_dimensions=None,
    ):
        
        VectorSpaceClusterer.__init__(self, normalise, svd_dimensions)
        self._means = numpy.array(initial_means, numpy.float64)
        self._num_clusters = len(initial_means)
        self._conv_threshold = conv_threshold
        self._covariance_matrices = covariance_matrices
        self._priors = priors
        self._bias = bias

    def num_clusters(self):
        return self._num_clusters

    def cluster_vectorspace(self, vectors, trace=False):
        assert len(vectors) > 0

        
        dimensions = len(vectors[0])
        means = self._means
        priors = self._priors
        if not priors:
            priors = self._priors = (
                numpy.ones(self._num_clusters, numpy.float64) / self._num_clusters
            )
        covariances = self._covariance_matrices
        if not covariances:
            covariances = self._covariance_matrices = [
                numpy.identity(dimensions, numpy.float64)
                for i in range(self._num_clusters)
            ]

        
        lastl = self._loglikelihood(vectors, priors, means, covariances)
        converged = False

        while not converged:
            if trace:
                print("iteration; loglikelihood", lastl)
            
            h = numpy.zeros((len(vectors), self._num_clusters), numpy.float64)
            for i in range(len(vectors)):
                for j in range(self._num_clusters):
                    h[i, j] = priors[j] * self._gaussian(
                        means[j], covariances[j], vectors[i]
                    )
                h[i, :] /= sum(h[i, :])

            
            for j in range(self._num_clusters):
                covariance_before = covariances[j]
                new_covariance = numpy.zeros((dimensions, dimensions), numpy.float64)
                new_mean = numpy.zeros(dimensions, numpy.float64)
                sum_hj = 0.0
                for i in range(len(vectors)):
                    delta = vectors[i] - means[j]
                    new_covariance += h[i, j] * numpy.multiply.outer(delta, delta)
                    sum_hj += h[i, j]
                    new_mean += h[i, j] * vectors[i]
                covariances[j] = new_covariance / sum_hj
                means[j] = new_mean / sum_hj
                priors[j] = sum_hj / len(vectors)

                
                covariances[j] += self._bias * numpy.identity(dimensions, numpy.float64)

            
            l = self._loglikelihood(vectors, priors, means, covariances)

            
            if abs(lastl - l) < self._conv_threshold:
                converged = True
            lastl = l

    def classify_vectorspace(self, vector):
        best = None
        for j in range(self._num_clusters):
            p = self._priors[j] * self._gaussian(
                self._means[j], self._covariance_matrices[j], vector
            )
            if not best or p > best[0]:
                best = (p, j)
        return best[1]

    def likelihood_vectorspace(self, vector, cluster):
        cid = self.cluster_names().index(cluster)
        return self._priors[cluster] * self._gaussian(
            self._means[cluster], self._covariance_matrices[cluster], vector
        )

    def _gaussian(self, mean, cvm, x):
        m = len(mean)
        assert cvm.shape == (m, m), "bad sized covariance matrix, %s" % str(cvm.shape)
        try:
            det = numpy.linalg.det(cvm)
            inv = numpy.linalg.inv(cvm)
            a = det**-0.5 * (2 * numpy.pi) ** (-m / 2.0)
            dx = x - mean
            print(dx, inv)
            b = -0.5 * numpy.dot(numpy.dot(dx, inv), dx)
            return a * numpy.exp(b)
        except OverflowError:
            
            
            return 0

    def _loglikelihood(self, vectors, priors, means, covariances):
        llh = 0.0
        for vector in vectors:
            p = 0
            for j in range(len(priors)):
                p += priors[j] * self._gaussian(means[j], covariances[j], vector)
            llh += numpy.log(p)
        return llh

    def __repr__(self):
        return "<EMClusterer means=%s>" % list(self._means)


def demo():
    

    from nltk import cluster

    

    vectors = [numpy.array(f) for f in [[0.5, 0.5], [1.5, 0.5], [1, 3]]]
    means = [[4, 2], [4, 2.01]]

    clusterer = cluster.EMClusterer(means, bias=0.1)
    clusters = clusterer.cluster(vectors, True, trace=True)

    print("Clustered:", vectors)
    print("As:       ", clusters)
    print()

    for c in range(2):
        print("Cluster:", c)
        print("Prior:  ", clusterer._priors[c])
        print("Mean:   ", clusterer._means[c])
        print("Covar:  ", clusterer._covariance_matrices[c])
        print()

    
    vector = numpy.array([2, 2])
    print("classify(%s):" % vector, end=" ")
    print(clusterer.classify(vector))

    
    vector = numpy.array([2, 2])
    print("classification_probdist(%s):" % vector)
    pdist = clusterer.classification_probdist(vector)
    for sample in pdist.samples():
        print(f"{sample} => {pdist.prob(sample) * 100:.0f}%")


if __name__ == "__main__":
    demo()