Problèmes de singularité dans le modèle de mélange gaussien

15

Dans le chapitre 9 du livre Reconnaissance de formes et apprentissage automatique, il y a cette partie sur le modèle de mélange gaussien:

entrez la description de l'image ici entrez la description de l'image ici Pour être honnête, je ne comprends pas vraiment pourquoi cela créerait une singularité. Quelqu'un peut-il m'expliquer cela? Je suis désolé mais je suis juste un étudiant de premier cycle et un novice en apprentissage automatique, donc ma question peut sembler un peu idiote, mais aidez-moi s'il vous plaît. Merci beaucoup

Dang Manh Truong
la source
On dirait que c'est aussi facile à réparer, reparamétrer à σk2=τ2γk et ensuite pénaliser γk pour être trop proche de zéro lors de l'optimisation.
probabilitéislogic
1
@probabilityislogic Je ne sais pas si je suis ici :(
Dang Manh Truong

Réponses:

11

Si nous voulons adapter un gaussien à un seul point de données en utilisant le maximum de vraisemblance, nous obtiendrons un gaussien très hérissé qui "s'effondrera" à ce point. La variance est nulle lorsqu'il n'y a qu'un seul point, ce qui dans le cas gaussien multivarié, conduit à une matrice de covariance singulière, donc on l'appelle le problème de singularité.

Lorsque la variance atteint zéro, la probabilité de la composante gaussienne (formule 9.15) passe à l'infini et le modèle devient surajusté. Cela ne se produit pas lorsque nous ajustons un seul gaussien à un certain nombre de points, car la variance ne peut pas être nulle. Mais cela peut arriver quand nous avons un mélange de gaussiens, comme illustré sur la même page de PRML.

entrez la description de l'image ici

Mise à jour :
le livre propose deux méthodes pour résoudre le problème de la singularité, qui sont

1) réinitialiser la moyenne et la variance en cas de singularité entrez la description de l'image ici

2) en utilisant MAP au lieu de MLE en ajoutant un prior. entrez la description de l'image ici

dontloo
la source
Concernant le cas gaussien unique, pourquoi la variance ne peut-elle pas être nulle? Le manuel dit: "Rappelez-vous que ce problème ne s'est pas posé dans le cas d'une distribution gaussienne unique. Pour comprendre la différence, notez que si un seul gaussien s'effondre sur un point de données, il contribuera à multiplier les facteurs pour la fonction de vraisemblance résultant de l'autre les points de données et ces facteurs iront à zéro exponentiellement rapidement, donnant une probabilité globale qui va à zéro plutôt qu'à l'infini. "mais je ne le comprends pas beaucoup :(
Dang Manh Truong
@DangManhTruong c'est parce que selon la définition de la variance, , à moins que tous les points soient de la même valeur, nous avons toujours une variance non nulle. var(x)=E[(xμ)2]
dontloo
Je vois! Merci: D Alors en pratique, que devons-nous faire pour l'éviter? Le livre n'explique pas cela.
Dang Manh Truong
@DangManhTruong salut, je l'ai ajouté à la réponse, veuillez jeter un oeil :)
dontloo
@DangManhTruong vous êtes les bienvenus
dontloo
3

Rappelons que ce problème ne s'est pas posé dans le cas d'une seule distribution gaussienne. Pour comprendre la différence, notez que si un seul gaussien s'effondre sur un point de données, il contribuera à des facteurs multiplicatifs à la fonction de vraisemblance résultant des autres points de données et ces facteurs iront à zéro exponentiellement rapidement, donnant une probabilité globale qui va à zéro plutôt que l'infini.

Je suis également un peu confus par cette partie, et voici mon interprétation. Prenez le boîtier 1D pour plus de simplicité.

Lorsqu'un seul gaussien "s'effondre" sur un point de données , c'est-à-dire μ = x i , la probabilité globale devient:xiμ=xi

p(x)=p(xi)p(xi)=(12πσ)(niN12πσe(xnμ)22σ2)

Vous voyez comme , le terme à gauche p ( x i ) , qui est comme le cas pathologique dans GMM, mais le terme à droite, qui est la probabilité d'autres points de données p ( xi ) , contient toujours des termes comme e - ( x n - μ ) 2σ0p(xi)p(xi) qui0exponentiellement rapide commeσ0, donc l'effet global sur la probabilité est qu'il aille à zéro.e(xnμ)22σ20σ0

Le point principal ici est que lors de l'ajustement d'un seul gaussien, tous les points de données doivent partager un ensemble de paramètres , contrairement au cas du mélange où un composant peut "se concentrer" sur un point de données sans pénaliser la probabilité globale des données .μ,σ

Yibo Yang
la source
2

Cette réponse donnera un aperçu de ce qui se passe qui conduit à une matrice de covariance singulière lors de l'ajustement d'un GMM à un ensemble de données, pourquoi cela se produit ainsi que ce que nous pouvons faire pour éviter cela.

Par conséquent, il est préférable de commencer par récapituler les étapes lors de l'ajustement d'un modèle de mélange gaussien à un ensemble de données.


0. Décidez du nombre de sources / clusters (c) que vous souhaitez adapter à vos données
1. Initialisez les paramètres moyenne , covariance Σ c et fraction_per_class π c par cluster c μcΣcπc

EStep_

  1. Calculez pour chaque point de données la probabilité r i c que le point de données x i appartient au cluster c avec: xiricxi

    N(x|μ
    ric=πcN(xi | μc,Σc)Σk=1KπkN(xi | μk,Σk)
    décrit la gaussienne multivariée avec: N ( x i , μ c , Σ c ) = 1N(x | μ,Σ)

    ricnous donne pour chaque point de donnéexila mesure de:Probabilitythatxibelongstoclas
    N(xi,μc,Σc) = 1(2π)n2|Σc|12exp(12(xiμc)TΣc1(xiμc))


    ricxi doncsixiest très proche deune gaussienne c, on obtiendra une hautericvaleur pour cette gaussienne et des valeurs relativement faibles sinon. M-Step_ Pour chaque cluster c: calculer le poids totalmcProbability that xi belongs to class cProbability of xi over all classesxiric

    MStep_

    mc(en gros, la fraction des points alloués au cluster c) et mettre à jour , μ c et Σ c en utilisant r i c avec: m c = Σ i r i c π c = m cπcμcΣcric

    mc = Σiric

    μc=1
    πc = mcm

    Σc=1
    μc = 1mcΣiricxi

    N'oubliez pas que vous devez utiliser les moyens mis à jour dans cette dernière formule. Répétez itérativement les étapes E et M jusqu'à ce que la fonction log-vraisemblance de notre modèle converge où la log-vraisemblance est calculée avec: lnp(X|π,μ,Σ)=Σ N i = 1 ln(Σ K
    Σc = 1mcΣiric(xiμc)T(xiμc)





    ln p(X | π,μ,Σ) = Σi=1N ln(Σk=1KπkN(xi | μk,Σk))



XAX=XA=I

[0000]


UNEXjeΣc-1
0matrice de covariance ci-dessus si la gaussienne multivariée tombe en un point pendant l'itération entre le pas E et M. Cela pourrait se produire si nous avons par exemple un ensemble de données auquel nous voulons adapter 3 gaussiens mais qui ne se compose en fait que de deux classes (clusters) de sorte que, grosso modo, deux de ces trois gaussiens capturent leur propre cluster tandis que le dernier gaussien ne le gère que pour attraper un seul point sur lequel il se trouve. Nous verrons à quoi cela ressemble ci-dessous. Mais étape par étape: supposons que vous ayez un ensemble de données bidimensionnel composé de deux grappes mais que vous ne le savez pas et que vous souhaitez y adapter trois modèles gaussiens, c'est-à-dire c = 3. Vous initialisez vos paramètres dans l'étape E et tracez les gaussiens au-dessus de vos données qui semblent smth. comme (peut-être vous pouvez voir les deux clusters relativement dispersés en bas à gauche et en haut à droite): enter image description hereμcπcenter image description here

rjeccovrjec
rjec=πcN(Xje | μc,Σc)Σk=1KπkN(Xje | μk,Σk)
vous voyez que là rjecLes valeurs auraient de grandes valeurs si elles sont très probablement sous le groupe c et des valeurs faibles sinon. Pour rendre cela plus apparent, considérons le cas où nous avons deux gaussiens relativement répartis et un gaussien très serré et nous calculons lerjec pour chaque point de donnée Xjecomme illustré dans la figure: enter image description here Parcourez donc les points de données de gauche à droite et imaginez que vous écririez la probabilité pour chaqueXjequ'il appartient au gaussien rouge, bleu et jaune. Ce que vous voyez, c'est que pour la plupart desXjela probabilité d'appartenir au gaussien jaune est très faible. Dans le cas ci-dessus où le troisième gaussien se trouve sur un seul point de données,rjec est seulement supérieur à zéro pour ce point de données, alors qu'il est nul pour tous les autres Xje. (s'effondre sur ce point de données -> Cela se produit si tous les autres points font plus probablement partie d'un ou deux gaussiens et c'est donc le seul point qui reste pour les trois gaussiens -> La raison pour laquelle cela se produit peut être trouvée dans l'interaction entre l'ensemble de données lui-même dans l'initialisation des gaussiens. Autrement dit, si nous avions choisi d'autres valeurs initiales pour les gaussiens, nous aurions vu une autre image et le troisième gaussien ne s'effondrerait peut-être pas). C'est suffisant si vous enfoncez de plus en plus ce gaussien. lerjecla table ressemble alors à smth. comme: enter image description here Comme vous pouvez le voir,rjecde la troisième colonne, c'est-à-dire pour le troisième gaussien sont nuls au lieu de cette seule ligne. Si nous recherchons quel point de données est représenté ici, nous obtenons le point de données: [23.38566343 8.07067598]. Ok, mais pourquoi obtient-on une matrice de singularité dans ce cas? Eh bien, et c'est notre dernière étape, nous devons donc à nouveau considérer le calcul de la matrice de covariance qui est:
Σc = Σjerjec(Xje-μc)T(Xje-μc)
nous avons vu que tout rjec sont zéro à la place pour celui Xjeavec [23.38566343 8.07067598]. Maintenant, la formule veut que nous calculions(Xje-μc). Si nous regardons leμcpour ce troisième gaussien on obtient [23.38566343 8.07067598]. Oh, mais attendez, c'est exactement la même chose queXje et c'est ce que Bishop a écrit avec: "Supposons que l'un des composants du modèle de mélange, disons le j th component, has its mean μj exactly equal to one of the data points so that μj=xn for some value of n" (Bishop, 2006, p.434). So what will happen? Well, this term will be zero and hence this datapoint was the only chance for the covariance-matrix not to get zero (since this datapoint was the only one where ric>0), it now gets zero and looks like:

[0000]


Consequently as said above, this is a singular matrix and will lead to an error during the calculations of the multivariate gaussian. So how can we prevent such a situation. Well, we have seen that the covariance matrix is singular if it is the 0 matrix. Hence to prevent singularity we simply have to prevent that the covariance matrix becomes a 0 matrix. This is done by adding a very little value (in sklearn's GaussianMixture this value is set to 1e-6) to the digonal of the covariance matrix. There are also other ways to prevent singularity such as noticing when a gaussian collapses and setting its mean and/or covariance matrix to a new, arbitrarily high value(s). This covariance regularization is also implemented in the code below with which you get the described results. Maybe you have to run the code several times to get a singular covariance matrix since, as said. this must not happen each time but also depends on the initial set up of the gaussians.
import matplotlib.pyplot as plt
from matplotlib import style
style.use('fivethirtyeight')
from sklearn.datasets.samples_generator import make_blobs
import numpy as np
from scipy.stats import multivariate_normal


# 0. Create dataset
X,Y = make_blobs(cluster_std=2.5,random_state=20,n_samples=500,centers=3)

# Stratch dataset to get ellipsoid data
X = np.dot(X,np.random.RandomState(0).randn(2,2))


class EMM:

    def __init__(self,X,number_of_sources,iterations):
        self.iterations = iterations
        self.number_of_sources = number_of_sources
        self.X = X
        self.mu = None
        self.pi = None
        self.cov = None
        self.XY = None



    # Define a function which runs for i iterations:
    def run(self):
        self.reg_cov = 1e-6*np.identity(len(self.X[0]))
        x,y = np.meshgrid(np.sort(self.X[:,0]),np.sort(self.X[:,1]))
        self.XY = np.array([x.flatten(),y.flatten()]).T


        # 1. Set the initial mu, covariance and pi values
        self.mu = np.random.randint(min(self.X[:,0]),max(self.X[:,0]),size=(self.number_of_sources,len(self.X[0]))) # This is a nxm matrix since we assume n sources (n Gaussians) where each has m dimensions
        self.cov = np.zeros((self.number_of_sources,len(X[0]),len(X[0]))) # We need a nxmxm covariance matrix for each source since we have m features --> We create symmetric covariance matrices with ones on the digonal
        for dim in range(len(self.cov)):
            np.fill_diagonal(self.cov[dim],5)


        self.pi = np.ones(self.number_of_sources)/self.number_of_sources # Are "Fractions"
        log_likelihoods = [] # In this list we store the log likehoods per iteration and plot them in the end to check if
                             # if we have converged

        # Plot the initial state    
        fig = plt.figure(figsize=(10,10))
        ax0 = fig.add_subplot(111)
        ax0.scatter(self.X[:,0],self.X[:,1])
        for m,c in zip(self.mu,self.cov):
            c += self.reg_cov
            multi_normal = multivariate_normal(mean=m,cov=c)
            ax0.contour(np.sort(self.X[:,0]),np.sort(self.X[:,1]),multi_normal.pdf(self.XY).reshape(len(self.X),len(self.X)),colors='black',alpha=0.3)
            ax0.scatter(m[0],m[1],c='grey',zorder=10,s=100)


        mu = []
        cov = []
        R = []


        for i in range(self.iterations):               

            mu.append(self.mu)
            cov.append(self.cov)


            # E Step
            r_ic = np.zeros((len(self.X),len(self.cov)))

            for m,co,p,r in zip(self.mu,self.cov,self.pi,range(len(r_ic[0]))):
                co+=self.reg_cov
                mn = multivariate_normal(mean=m,cov=co)
                r_ic[:,r] = p*mn.pdf(self.X)/np.sum([pi_c*multivariate_normal(mean=mu_c,cov=cov_c).pdf(X) for pi_c,mu_c,cov_c in zip(self.pi,self.mu,self.cov+self.reg_cov)],axis=0)
            R.append(r_ic)

            # M Step

            # Calculate the new mean vector and new covariance matrices, based on the probable membership of the single x_i to classes c --> r_ic
            self.mu = []
            self.cov = []
            self.pi = []
            log_likelihood = []

            for c in range(len(r_ic[0])):
                m_c = np.sum(r_ic[:,c],axis=0)
                mu_c = (1/m_c)*np.sum(self.X*r_ic[:,c].reshape(len(self.X),1),axis=0)
                self.mu.append(mu_c)

                # Calculate the covariance matrix per source based on the new mean
                self.cov.append(((1/m_c)*np.dot((np.array(r_ic[:,c]).reshape(len(self.X),1)*(self.X-mu_c)).T,(self.X-mu_c)))+self.reg_cov)
                # Calculate pi_new which is the "fraction of points" respectively the fraction of the probability assigned to each source 
                self.pi.append(m_c/np.sum(r_ic)) 



            # Log likelihood
            log_likelihoods.append(np.log(np.sum([k*multivariate_normal(self.mu[i],self.cov[j]).pdf(X) for k,i,j in zip(self.pi,range(len(self.mu)),range(len(self.cov)))])))



        fig2 = plt.figure(figsize=(10,10))
        ax1 = fig2.add_subplot(111) 
        ax1.plot(range(0,self.iterations,1),log_likelihoods)
        #plt.show()
        print(mu[-1])
        print(cov[-1])
        for r in np.array(R[-1]):
            print(r)
        print(X)

    def predict(self):
        # PLot the point onto the fittet gaussians
        fig3 = plt.figure(figsize=(10,10))
        ax2 = fig3.add_subplot(111)
        ax2.scatter(self.X[:,0],self.X[:,1])
        for m,c in zip(self.mu,self.cov):
            multi_normal = multivariate_normal(mean=m,cov=c)
            ax2.contour(np.sort(self.X[:,0]),np.sort(self.X[:,1]),multi_normal.pdf(self.XY).reshape(len(self.X),len(self.X)),colors='black',alpha=0.3)




EMM = EMM(X,3,100)     
EMM.run()
EMM.predict()
2Obe
la source
0

À mon humble avis, toutes les réponses manquent un fait fondamental. Si l'on regarde l'espace des paramètres d'un modèle de mélange gaussien, cet espace est singulier le long du sous-espace où il y a moins que le nombre total de composants dans le mélange. Cela signifie que les dérivées sont automatiquement nulles et généralement tout le sous-espace apparaîtra comme un mle. Plus philosophiquement, le sous-espace des covariances de rang inférieur au rang complet est la limite de l'espace des paramètres et il faut toujours se méfier lorsque le mle se produit sur la frontière - cela indique généralement qu'il y a un plus grand espace de paramètres qui se cache dans lequel on peut trouver le «vrai» mle. Il existe un livre intitulé "Statistiques algébriques" de Drton, Sturmfeld et Sullivant. Cette question est discutée dans ce livre en détail. Si vous êtes vraiment curieux, vous devriez regarder cela.

meh
la source
-2

For a single Gaussian, the mean may possibly equal one of the data points (xn for example) and then there is the following term in the likelihood function:

N(xn|xn,σj11)limσjxn1(2π)1/2σjexp(1σj|xnσj|2)=1(2π)1/2σj
The limit σj0 is now clearly divergent since the argument of the exponential vanishes.

However for a data point xm different from the mean σj, we will have

N(xm|xm,σj11)=1(2π)1/2σjexp(1σj|xmσj|2)
and now the argument of the exponential diverges (and is negative) in the limit σj0. As a result the product of these two terms in the likelihood function will vanish.
Nick
la source
This answer is incorrect as there is no reason to identify mean μj and standard deviation σj.
Xi'an