Mixed Gaussian model (MATLAB)

Recommended Blog: http://blog.csdn.net/crzy_sparrow/article/details/7413019

There are many kinds of background models, many of which have insufficient adaptability to the mutation of illumination and other factors, and Gaussian mixture model is one of the most successful background modeling methods.

Gaussian background model is a classical adaptive mixed Gaussian background extraction method proposed by Stauffer and other people, and is a background modeling method, which is based on the distribution of each pixel in the time domain to construct the color distribution model of each pixel, in order to achieve the purpose of background modeling. The mixed Gaussian background model is a weighted sum of the finite Gaussian functions, which can describe the multi-peak state of pixels, and is suitable for accurate modeling of the complex backgrounds such as light gradients and tree sway. After the continuous improvement of many researchers, the method has become a more commonly used background extraction method.

The principle of mixed Gaussian model in a nutshell: The probability distributions of arbitrary shapes can be approximated by multiple Gaussian distribution functions. In other words, the Gaussian mixture model (GMM) can smooth the probability distribution of arbitrary shapes. The method of parametric solution is usually solved by maximum likelihood estimation, but the maximum likelihood estimation method is often not used to obtain complete data (such as the sample is known, but the sample class (which belongs to which Gaussian distribution) is unknown), so the EM (maximal expected value) solution is presented.

Although the above-mentioned simple, but mixed Gaussian model and EM solution theory is still relatively complex, I found I think I can quickly grasp the Gaussian mixture model of the data package to the attachment, you can download, understand the mixed Gaussian model and the full derivation of EM.

Attachment Download Address:

http://download.csdn.net/detail/crzy_sparrow/4187994

1) Arbitrary data distribution can be represented by a Gaussian mixture model (M-Gaussian) ((1))

which

2) The log likelihood function for n sample set X is as follows:

(1) Initial parameters

(2) E-step (request for expectation)

Fetch:

In fact, the most appropriate formula here should be the expected expression of the log-likelihood function

In fact, 3) the calculation of the parameters of the log likelihood function is also used to find the desired expression of the derivative is equal to 0 of the simplified formula, and these formulas for (7), so the e-step can only ask for (7) formula, in order to simplify the calculations, do not need to seek partial guidance each time.

(3) m step (Maximize Step)

(4) Closing conditions

Iterative em step, update parameters, until the likelihood function before and after the difference is less than a threshold, or before and after the difference between the parameters (generally choose Euclidean distance) less than a certain threshold, terminate the iteration, here to choose the former.

The matlab corresponding program is as follows:

(1) Initial parameters

    function [Pmiu pPi Psigma] = init_params ()% initialization parameter
        Pmiu = centroids;% k-by-d matrix
        pPi = zeros (1, K);%1-by-k matrix
        Psigma = zeros (d, D, k);% hard
 
        assign × to each centroids
        Distmat = Repmat (sum (x.*x, 2), 1, K) + ...% x is A n-by-d data matrix.
            Repmat (SUM (Pmiu.*pmiu, 2) ', N, 1)-...% x->k column u->n row xu^t
            is n-by-k 2*x*pmiu ';% calculates the distance from each point to the K center
        [~, labels] = m In (Distmat, [], 2);% found closest to X Pmiu,[c,i] labels represents this minimum value is selected from that column for
 
        k=1:k
            Xk = X (Labels = = k,:);% Xk is all X vectors classified into K class The Matrix
            pPi (k) = Size (Xk, 1)/n;% count a few of the
            Psigma (:,:, k) = CoV (Xk) of the K class; Calculate covariance matrix, d-by-d matrix, unbiased estimation of minimum variance
        end
    End

Where Pmiu is the class center of the K class, K*d,d is the number of columns in the sample x, and a column in the sample x is very much a sample, (for example, the dots and z in space)

PPi for each class of proportional number, 1*k, such as a total of sample points, which belong to the first category accounted for 0.2, belonging to the low two class accounted for 0.1, belongs to the category K ...

       psigma is a covariance matrix for each class, D*d*k, representing the covariance matrix of class K as Psigma (:,:, K)

</pre></p><p><span style= "FONT-SIZE:16PX; line-height:26px; font-family:arial; Background-color:rgb (255, 255, 255); " ><span style= "color: #ff0000;" ><span style= "Color:rgb (255, 0, 0); font-family:arial; font-size:16px; line-height:26px; " ><span style= "Color:rgb (51, 51, 51); font-family:arial; font-size:16px; line-height:26px; " >       </span><span style= "font-family:arial; font-size:16px; line-height:26px; Color:rgb (255, 0, 0); " > (</span><span style= "font-family:arial; line-height:26px; Color:rgb (255, 0, 0); font-size:18px; " &GT;2) calculates N</span></span></span></span><span style= "Color:rgb (255, 0, 0); font-family:arial; font-size:18px; line-height:26px; " >         </span></p><p><span style= "font-family:arial; font-size:18px; line-height:26px; " ></span><pre name= "code" class= "CPP" > Function Px = calc_Prob () Px = Zeros (N, K);
            For k = 1:k Xshift = X-repmat (Pmiu (k,:), N, 1);%x-u lemda=1e-5; Conv=psigma (:,:, K) +lemda*diag (Diag (Ones (D))); The% handles the singular problem here, plus a very small lemda*i Inv_psigma = INV (CONV) for the covariance matrix;% covariance The inverse tmp = SUM ((xshift*inv_psigma). * Xshift, 2);% ( X-u_k) sigma.* (x-u_k), TMP is a n*1 vector coef = (2*pi) ^ (-D/2) * sqrt (det (inv_psigma));% Front parameters Px (:, k) = Co EF * EXP ( -0.5*tmp);% takes data point x into the Gaussian model to get the value end end

N in this n-Extreme equation (1)

(3) Calculate and update Ppi,pmiu, Psigma;

Lprev =-inf;
    While true
        Px = Calc_prob (),% calculates N (x|mu,sigma)
 
        % new value for Pgamma
        pgamma = Px. * Repmat (PPi, N, 1);% estimated Gamma is a n*k matrix
        pgamma = Pgamma./Repmat (SUM (Pgamma, 2), 1, K);% of each sample again K clustering class, also known as componen generation probability
 
        % new value for parameters O F each Component
        Nk = SUM (pgamma, 1);%n_k
        Pmiu = diag (1./nk) * pgamma ' * X;% number * (K-BY-N * n-by-d) plus parentheses to help understand
        p Pi = nk/n;
        For KK = 1:k
            xshift = X-repmat (Pmiu (KK,:), N, 1);%x-u
            Psigma (:,:, KK) = (xshift ' * ...
                ) (Diag (Pgamma (:, KK)) * Xshift)) /Nk (KK);% update sigma
        end
 
        % check for convergence
        L = SUM (log (Px*ppi '));
        If L-lprev < threshold break
            ;
        End
        Lprev = L;
    End

The PX in the equation is the formula (1)

The Pgamma in the above equation is the formula (7)

Ppi,pmiu in the above style, Psigma is updated as per m step.

The stop criterion in the above equation is calculated according to the formula (4)