混和高斯模型的推导和实现
混合高斯模型的简要介绍

混合高斯模型跟高斯变量之和看起来有一点像, 注意不要把它们弄混淆了. 混合高斯模型给出的概率密度函数实际上是几个高斯概率密度函数的加权和:计算均值和方差的公式不仅适用于几个(多维)高斯分布混合的情况, 还适用于非高斯分布的情况.高斯变量之和就没什么好说的了, 几个高斯变量之和是一个新的高斯变量.原理: 高斯模型就是用高斯概率密度函数(正态分布曲线)精确地量化事物,将一个事物分解为若干的基于高斯概率密度函数(正态分布曲线)形成的模型。
对图像背景建立高斯模型的原理及过程:图像灰度直方图反映的是图像中某个灰度值出现的频次,也可以认为是图像灰度概率密度的估计。
如果图像所包含的目标区域和背景区域相比比较大,且背景区域和目标区域在灰度上有一定的差异,那么该图像的灰度直方图呈现双峰-谷形状,其中一个峰对应于目标,另一个峰对应于背景的中心灰度。
对于复杂的图像,尤其是医学图像,一般是多峰的。
通过将直方图的多峰特性看作是多个高斯分布的叠加,可以解决图像的分割问题。
在智能监控系统中,对于运动目标的检测是中心内容,而在运动目标检测提取中,背景目标对于目标的识别和跟踪至关重要。
而建模正是背景目标提取的一个重要环节。
我们首先要提起背景和前景的概念,前景是指在假设背景为静止的情况下,任何有意义的运动物体即为前景。
建模的基本思想是从当前帧中提取前景,其目的是使背景更接近当前视频帧的背景。
即利用当前帧和视频序列中的当前背景帧进行加权平均来更新背景,但是由于光照突变以及其他外界环境的影响,一般的建模后的背景并非十分干净清晰,而高斯混合模型是是建模最为成功的方法之一。
混合高斯模型使用K(基本为3到5个)个高斯模型来表征图像中各个像素点的特征,在新一帧图像获得后更新混合高斯模型, 用当前图像中的每个像素点与混合高斯模型匹配,如果成功则判定该点为背景点, 否则为前景点。
通观整个高斯模型,主要是有方差和均值两个参数决定,对均值和方差的学习,采取不同的学习机制,将直接影响到模型的稳定性、精确性和收敛性。
高斯混合模型详解

高斯混合模型详解【原创实用版】目录1.高斯混合模型的基本概念2.高斯混合模型的组成部分3.高斯混合模型的推导过程4.高斯混合模型的应用实例5.总结正文一、高斯混合模型的基本概念高斯混合模型(Gaussian Mixture Model,简称 GMM)是一种概率模型,用于对由多个高斯分布组成的数据集进行建模。
它是一个多元高斯分布,由多个一元高斯分布组合而成,每个一元高斯分布表示数据集中的一个子集。
通过学习这些子集的参数,即均值、协方差矩阵和权重,高斯混合模型可以对数据集进行有效的建模。
二、高斯混合模型的组成部分高斯混合模型由以下三个部分组成:1.均值(Means):每个高斯分布的均值向量,表示该分布的中心点。
2.协方差矩阵(Covariances):每个高斯分布的协方差矩阵,表示该分布的形状。
3.权重(Weights):每个高斯分布的权重,表示该分布在数据集中的重要性。
三、高斯混合模型的推导过程高斯混合模型的推导过程主要包括两个步骤:1.初始化:随机设置初始的均值、协方差矩阵和权重,这将影响优化过程的收敛速度。
2.优化:使用期望最大化(Expectation-Maximization,简称 EM)算法来优化模型参数。
EM 算法通过迭代更新均值、协方差矩阵和权重,使得模型对观测数据的似然函数最大化。
四、高斯混合模型的应用实例高斯混合模型在许多领域都有广泛的应用,例如:1.语音识别:高斯混合模型可以用来对语音信号进行建模,从而实现语音识别。
2.机器学习:高斯混合模型可以用来对数据集进行聚类,从而实现机器学习任务。
3.信号处理:高斯混合模型可以用来对信号进行建模,从而实现信号处理任务。
五、总结高斯混合模型是一种强大的概率模型,可以用来对复杂的数据集进行建模。
通过学习均值、协方差矩阵和权重,高斯混合模型可以有效地表示数据集中的潜在结构。
高斯混合模型 c语言算法

高斯混合模型 c语言算法高斯混合模型 C 语言算法一、引言高斯混合模型(Gaussian Mixture Model,简称 GMM)是一种用于概率建模和数据聚类的统计模型。
它是由多个高斯分布组成的混合模型,每个高斯分布对应一个聚类簇。
C 语言是一种广泛应用于嵌入式系统和底层开发的编程语言。
本文将介绍如何使用 C 语言实现高斯混合模型算法。
二、高斯混合模型算法原理1. 高斯分布高斯分布是一种连续概率分布,也称为正态分布。
它的概率密度函数可以通过以下公式计算:```f(x) = (1 / (σ * √(2π))) * e^(-((x - μ)^2) / (2 * σ^2)) ```其中,μ 是分布的均值,σ 是分布的标准差。
2. 高斯混合模型高斯混合模型是由多个高斯分布组成的混合模型。
每个高斯分布都对应一个聚类簇,用来表示数据的不同类别或聚集程度。
高斯混合模型的概率密度函数可以表示为:```f(x) = Σ(w_i * f_i(x))```其中,w_i 是第 i 个高斯分布的权重,f_i(x) 是第 i 个高斯分布的概率密度函数。
3. 高斯混合模型的参数估计高斯混合模型的参数估计是通过最大似然估计方法来实现的。
具体步骤如下:- 初始化每个高斯分布的均值、标准差和权重;- 重复以下步骤直到收敛:- E 步:根据当前参数估计每个样本属于每个聚类的概率;- M 步:根据当前样本的权重更新每个聚类的参数估计;- 根据最终的参数估计得到高斯混合模型。
三、C 语言实现高斯混合模型算法1. 数据结构定义我们需要定义一些数据结构来表示高斯混合模型的参数和样本数据。
例如,可以定义一个结构体来表示每个高斯分布的参数:```ctypedef struct {double mean; // 均值double variance; // 方差double weight; // 权重} Gaussian;```2. 初始化参数在开始参数估计之前,我们需要初始化每个高斯分布的参数。
高斯混合模型(转)

⾼斯混合模型(转)混合⾼斯模型的原理说⽩了就⼀句话:任意形状的概率分布都可以⽤多个⾼斯分布函数去近似。
换句话说,⾼斯混合模型(GMM)可以平滑任意形状的概率分布。
其参数求解⽅法⼀般使⽤极⼤似然估计求解,但使⽤极⼤似然估计法往往不能获得完整数据(⽐如样本已知,但样本类别(属于哪个⾼斯分布)未知),于是出现了EM(最⼤期望值)求解⽅法。
虽然上⾯说的简单,但是混合⾼斯模型和EM求解的理论还是⽐较复杂的,我把我所找到的我认为能够快速掌握⾼斯混合模型的资料打包到了附件中,⼤家可以去下载,了解混合⾼斯模型以及EM的完整推导过程。
附件下载地址:⼤致抽取下⾼斯混合模型的重要概念:1)任意数据分布可⽤⾼斯混合模型(M个单⾼斯)表⽰((1)式)(1)其中:(2)(3)表⽰第j个⾼斯混合模型2)N个样本集X的log似然函数如下:(4)其中:参数:,(5)下⾯具体讲讲基于EM迭代的混合⾼斯模型参数求解算法流程:1)初始参数由k-means(其实也是⼀种特殊的⾼斯混合模型)决定:function [Priors, Mu, Sigma] = EM_init_kmeans(Data, nbStates)% Inputs -----------------------------------------------------------------% o Data: D x N array representing N datapoints of D dimensions.% o nbStates: Number K of GMM components.% Outputs ----------------------------------------------------------------% o Priors: 1 x K array representing the prior probabilities of the% K GMM components.% o Mu: D x K array representing the centers of the K GMM components.% o Sigma: D x D x K array representing the covariance matrices of the% K GMM components.% Comments ---------------------------------------------------------------% o This function uses the 'kmeans' function from the MATLAB Statistics% toolbox. If you are using a version of the 'netlab' toolbox that also% uses a function named 'kmeans', please rename the netlab function to% 'kmeans_netlab.m' to avoid conflicts.[nbVar, nbData] = size(Data);%Use of the 'kmeans' function from the MATLAB Statistics toolbox[Data_id, Centers] = kmeans(Data', nbStates);Mu = Centers';for i=1:nbStatesidtmp = find(Data_id==i);Priors(i) = length(idtmp);Sigma(:,:,i) = cov([Data(:,idtmp) Data(:,idtmp)]');%Add a tiny variance to avoid numerical instabilitySigma(:,:,i) = Sigma(:,:,i) + 1E-5.*diag(ones(nbVar,1));endPriors = Priors ./ sum(Priors);2)E步(求期望)求取:(7)其实,这⾥最贴切的式⼦应该是log似然函数的期望表达式事实上,3)中参数的求取也是⽤log似然函数的期望表达式求偏导等于0解得的简化后的式⼦,⽽这些式⼦⾄于(7)式有关,于是E步可以只求(7)式,以此简化计算,不需要每次都求偏导了。
高斯混合聚类及EM实现

⾼斯混合聚类及EM实现⼀、引⾔ 我们谈到了⽤ k-means 进⾏聚类的⽅法,这次我们来说⼀下另⼀个很流⾏的算法:Gaussian Mixture Model (GMM)。
事实上,GMM 和 k-means 很像,不过 GMM 是学习出⼀些概率密度函数来(所以 GMM 除了⽤在 clustering 上之外,还经常被⽤于 density estimation ),简单地说,k-means 的结果是每个数据点被 assign 到其中某⼀个 cluster 了,⽽ GMM 则给出这些数据点被 assign 到每个 cluster 的概率,⼜称作soft assignment。
给出⼀个概率有很多好处,因为它的信息量⽐简单的⼀个结果要多,⽐如,我可以把这个概率转换为⼀个 score ,表⽰算法对⾃⼰得出的这个结果的把握。
也许我可以对同⼀个任务,⽤多个⽅法得到结果,最后选取“把握”最⼤的那个结果;另⼀个很常见的⽅法是在诸如疾病诊断之类的场所,机器对于那些很容易分辨的情况(患病或者不患病的概率很⾼)可以⾃动区分,⽽对于那种很难分辨的情况,⽐如,49% 的概率患病,51% 的概率正常,如果仅仅简单地使⽤ 50% 的阈值将患者诊断为“正常”的话,风险是⾮常⼤的,因此,在机器对⾃⼰的结果把握很⼩的情况下,会“拒绝发表评论”,⽽把这个任务留给有经验的医⽣去解决。
⼆、归纳法 ⼀系列数据N要求对他拟合,如果不作要求的话,可以⽤⼀个N-1次多项式来完美拟合这N个点,⽐如拉格朗⽇插值,⽜顿插值等,或者如果不限制次数的话可以找到⽆穷个完美函数,但是往往咱们要求指数型或者线性,就是需要对信息有机结合,GMM就是这样,要求⽤⾼斯模型来拟合,当然了可以构造任意的混合模型,不过根据中⼼极限定理等⾼斯模型⽐较合适。
另外,Mixture Model 本⾝其实也是可以变得任意复杂的,通过增加 Model 的个数,我们可以任意地逼近任何连续的概率密分布。
高斯混合模型详细推导

⾼斯混合模型详细推导⼀、⾼斯混合模型定义 ⾼斯混合模型具有如下概率分布形式:P(y|θ)=K∑k=1αkϕ(y|θk) (9.24)其中,α是系数,αk⩾,\sum\limits_{k=1}^{K}\alpha_k=1;\phi(y|\theta_k)是⾼斯分布,\theta_k=(\mu_k,\sigma_{k}^{2}),\phi(y|\theta_k)=\frac{1}{\sqrt{2\pi}\sigma_k}exp\left(-\frac{(y-\mu_k)^2}{2\sigma_{k}^{2}}\right) (9.25)为第k个⾼斯分布。
⼆、确定模型的对数似然函数 设想观测数据y_j,j=1,2,\cdots,N是这样产⽣的:⾸先依概率\alpha_k选择第k个⾼斯分布模型\phi(y|\theta_k);然后依第k个分模型的概率分布\phi(y|\theta_k)⽣成观测数据y_j。
反映观测数据y_j来⾃第k个分模型的数据是未知的,k=1,2,\cdots,K,以隐变量\gamma_{jk}表⽰,其定义如下:\gamma_{jk}=1,\ if\ y_j\ from\ \phi_k\\ \gamma_{jk}=0,\ else (9.27)\gamma_{jk}是0-1随机变量。
模型的对数似然函数为logP(y|\theta)=log\prod\limits_{j=1}^{N}[\sum\limits_{k=1}^{K}\alpha_k\phi(y_j|\theta_k)]=\sum\limits_{j=1}^{N}log\sum\limits_{k=1}^{K}\alpha_k\phi(y_j|\theta_k) (9.28) 由于该式⼦包含累和的对数的形式,直接⽤极⼤似然法处理很困难,在这⾥采⽤EM算法进⾏参数求解。
三、EM算法的E步:确定Q函数Q(\theta,\theta^{(i)})=\sum\limits_{Z}P(Z|Y,\theta^{(i)})logP(Y,Z|\theta)=\sum\limits_{k=1}^{K}\sum\limits_{j=1}^{N}P(\gamma_{jk}|y_j,\theta^{(i)})log\alpha_k\phi(y_j|\theta_k) 由于P(\gamma_{jk}|y_j,\theta^{(i)})表⽰在观测数据和该次迭代参数的条件下,数据y_j来⾃⾼斯分布k的概率,因此易得P(\gamma_{jk}|y_j,\theta^{(i)})=\frac{\alpha_k^{(i)}\phi(y_j|\theta_k^{(i)})}{\sum\limits_{k=1}^{K}\alpha_k^{(i)}\phi(y_j|\theta_k^{(i)})}=\hat{\gamma_{jk}}\hat{\gamma_{jk}}表⽰分模型k对观测数据j的响应度, 所以Q(\theta,\theta^{(i)})=\sum\limits_{k=1}^{K}\left\{\sum\limits_{j=1}^{N}\hat{\gamma_{jk}}log\alpha_k+\sum\limits_{j=1}^{N}\hat{\gamma_{jk}}\left[log\left(\frac{1}{\sqrt(2\pi)}\right)-log\sigma_k-\frac{1}{2\sigma_{k}^{2}}(y_j-\mu_k)^2\right]\right\}\\=\sum\limits_{k=1}^{K}n_klog\alpha_k+\sum\limits_{k=1}^{K}\sum\limits_{j=1}^{N}\hat{\gamma_{jk}}\left[log\left(\frac{1}{\sqrt{2\pi}}\right)-log\sigma_k-\frac{1}{2\sigma_{k}^{2}}(y_j-\mu_k)^2\right] (9.29)其中n_k=\sum_{j=1}^{N}\hat{\gamma_{jk}}。
聚类之高斯混合模型与EM算法

聚类之⾼斯混合模型与EM算法⼀、⾼斯混合模型概述1、公式⾼斯混合模型是指具有如下形式的概率分布模型:其中,αk≥0,且∑αk=1,是每⼀个⾼斯分布的权重。
Ø(y|θk)是第k个⾼斯分布的概率密度,被称为第k个分模型,参数为θk=(µk, αk2),概率密度的表达式为:⾼斯混合模型就是K个⾼斯分布的线性组合,它假设所有的样本可以分为K类,每⼀类的样本服从⼀个⾼斯分布,那么⾼斯混合模型的学习过程就是去估计K个⾼斯分布的概率密度Ø(y|θk),以及每个⾼斯分布的权重αk。
每个观测样本出现的概率就表⽰为K个⾼斯分布概率的加权。
所谓聚类,就是对于某个样本y j,把该样本代⼊到K个⾼斯分布中求出属于每个类别的概率:然后选择概率值最⾼的那个类别作为它最终的归属。
把所有的样本分别归⼊K个类,也就完成了聚类的过程。
2、案例假设有 20 个⾝⾼样本数据,并不知道每个样本数据是来⾃男⽣还是⼥⽣。
在这种情况下,如何将这 20 个⾝⾼数据聚成男⼥⽣两⼤类呢?⽤⾼斯混合模型来聚类,那么假设男⼥⽣⾝⾼分别服从两个不同的⾼斯分布,⾼斯混合模型就是由男⽣⾝⾼和⼥⽣⾝⾼这两个⾼斯分布混合⽽成。
在⾼斯混合模型中,样本点属于某⼀类的概率不是⾮0即 1 的,⽽是属于不同类有不同的概率值。
如下图,有两个⾼斯分布,均值分别为µ1和µ2,⽽⾼斯混合模型就是⼜这两个⾼斯分布的概率密度线性组合⽽成。
⼆、⾼斯混合模型参数估计的EM算法假设观测数据y1, y2, ...y N由⾼斯混合模型⽣成:其中,要估计的参数θ=(α1, α2, ...αK; θ1, θ2, ..., θK),θk=(µk, αk2),k=1,2,...,K。
因此如果⾼斯混合模型由K个⾼斯分布混合⽽成,那么就有3K个参数需要估计。
我们⽤极⼤似然估计法来估计参数θ,也就是求参数θ,使得观测数据y的对数似然函数L(θ)=logP(y|θ)的极⼤化:由于对数似然函数L(θ)中包含了和的对数,⽐较难以求解,因此考虑⽤EM算法。
EM算法详细例子及推导

EM算法详细例子及推导EM算法(Expectation-Maximization Algorithm)是一种用于求解含有隐变量(latent variable)的概率模型的参数估计方法。
其基本思想是通过迭代的方式,通过观测数据得到对隐变量的估计,然后再基于该估计对模型参数进行优化。
下面我们以一个简单的高斯混合模型为例,详细介绍EM算法的推导和实例。
1. 高斯混合模型(Gaussian Mixture Model, GMM)高斯混合模型是一种概率模型,由多个高斯分布组合而成。
假设我们观测到的数据由K个高斯分布组成,每个高斯分布对应一个参数向量:均值miu和方差sigma^2、同时,我们还有一个隐变量Z,表示观测数据属于哪个高斯分布,取值范围为{1,2,...,K}。
2.EM算法EM算法的核心思想是通过交替进行两个步骤:E步(Expectation)和M步(Maximization)。
在E步中,我们对当前模型参数下的隐变量进行估计,得到对隐变量的最大似然估计。
在M步中,我们利用得到的隐变量估计更新模型参数,使模型对观测数据的似然函数最大化。
不断重复这两步直至模型收敛。
下面我们通过具体的例子来推导EM算法。
假设我们观测到了一个数据集X = {x1, x2, ..., xn},我们希望通过EM算法对其进行建模。
Step1: 初始化模型参数首先,我们需要初始化模型参数。
选择K个高斯分布的参数miu和sigma^2,并假设所有的高斯分布对应的隐变量Z服从均匀分布。
这时,我们得到了初始模型参数Theta = {miu1, sigma^21, ..., miuK,sigma^K, pi1, pi2, ..., piK}。
Step2: E步,计算隐变量的后验分布在E步中,我们计算隐变量的后验分布。
对于每个观测样本xi,我们计算其属于每个高斯分布的概率,即:gamma(k,i) = P(Zi=k,xi, Theta) = P(Zi=k,xi, miu_k,sigma_k^2) = pi_k * N(xi,miu_k, sigma_k^2) / sum(pi_j * N(xi,miu_j, sigma_j^2), j=1 to K其中N(xi,miu_k, sigma_k^2)表示xi在第k个高斯分布下服从的概率密度函数。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
基于GMM 的运动目标检测方法研究一、GMM 数学公式推导1、预备知识:(1)设离散型随机变量X 的分布率为: {} 2,1,P ===k p a X k k 则称()∑=kk kp aX E 为X 的数学期望或均值(2)设连续型随机变量X 的概率密度函数(PDF )为f(x) 其数学期望定义为:()()dx x xf X E ⎰+∞∞-=(3)()()()[]2X E X E X D -=称为随机变量x 的方差,()X D 称为X的标准差(4)正态分布:()2,~σμN X 概率密度函数为:()()⎥⎥⎦⎤⎢⎢⎣⎡--=22221σμσπx e x p(5)设(x,y)为二维随机变量,()[]()[]{}Y E Y X E X E --若存在,则 称其为X 和Y 的协方差,记为cov(x,y)()()[]()[]{}()XY E Y E Y X E X E Y X =--=,cov 2、单高斯模型:SGM (也就是多维正态分布) 其概率密度函数PDF 定义如下: ()()()()μμπμ----=x C x nT eCC x N 12121,;其中,x 是维数为n 的样本向量(列向量),μ是期望,C 是协方差矩阵,|C|表示C 的行列式,1-C 表示C 的逆矩阵,()Tx μ-表示()μ-x 的转置。
3、混合高斯模型:GMM设想有 m 个类:m 321ϖϖϖϖ,,,, ,每类均服从正态分布。
各分布的中心点(均值)分别为:m 321μμμμ,,,,方差分别为:m 321σσσσ,,,,每一类在所有的类中所占的比例为 ()()()()m P P P P ϖϖϖϖ,,,,321 其中()11=∑=mi i P ϖ。
同时,已知个观察点:。
其中,用大写P 表示概率,用小写p 表示概率密度。
则依此构想,可得概率密度函数为:()()()()()()()()()()()μμπϖϖσμϖσμϖσμ---=-∑=⋅++⋅+⋅=x C x mi d i m m m T eCP P N P N P N x p 12112221112,,,其中d 是维数,|·|是行列式但是在利用GMM 进行目标检测时,这些模型的参数可能已知,也可能不知道,当参数已知时,可以直接利用GMM 进行目标检测,在未知的情况下,需要对参数进行估计。
对参数估计时,还要考虑样本分类是否已知。
(1)样本已知: 最大似然估计:可以直接采用MLE (最大似然估计)进行参数估计: 未知量为集合:()()()m P P C C ϖϖμμλ,,1m 1m 1 ,,,,,,= 将衡量概率密度函数优劣的标准写出:()()∏==nk k x P x p 1||λλ即为:()()()()()i k T i k x C x n k mi di eC P x p μμπϖλ---==-∏∑=12111||2|只要定出该标准的最大值位置,就可以求出最优的待定参数。
为了 求出这个最大值的位置,就需用导数求极点,具体求解过程于下:()()()()()∑∑==∑==∏===n k P x N nk x P x P x p mi i i k k nk k 1,1|||11ln ln ln ln ϖλλλλ求导:()()()()()()()()()()()()}||2{1},{1ln ln 1-112111|11,1,1|i k T i k k mi i i k mi i i k x C x mi di mk x p mi i i k mk P x N nk P x N nk x p eC P P x N μμλϖλϖλλπϖλϖλλλλ---======∑∂∂=∑∂∂∑=∑∂∂=∂∂∑∑∑∑==然后再分别对各个参数求导:①求参数iμ :②对 感兴趣,求偏导数有:③对感兴趣,接下来的求导比较复杂,在此就没有继续推导。
(2)样本未知: EM 估计,算法流程: ①初始化:方案1:协方差矩阵0j C 设为单位矩阵,每个模型比例的先验概率设为Mj 10=α,均值0j μ为随机数。
方案2:有K 均值(K-means)聚类算法对样本进行聚类,利用各类的均值作为0j μ,并计算0j C ,0j α去各类样本占总数的比例。
②估计步骤(E-step ): 令j α的后验概率为: ()()M j n i x N x N Mk i k ki j j ij ≤≤≤≤=∑=1,1,||1φαφαβ③最大化步骤(M-step ):更新权值:nni ijj ∑==1βα更新均值:∑∑===n i ni iji j ijx 11ββμ更新方差矩阵:()()∑∑==--=ni ijni Tiiiiijj x x C 11βμμβ④收敛条件:不断地迭代步骤②和③,重复更新上面的三个值,直到()()εφφ<-||'|X p X p ,其中为更新参数后计算的值,即前后两次迭代得到的结果变化小于一定程度则终止迭代,通常-510=ε 二、GMM 发展历史及现状背景建模方法有很多种,如中值法、均值法、卡尔曼滤波器模型、码本背景模型等,其中混合高斯模型是最经典的算法。
GMM 最早是由CHris Stauffer 等在[1]中提出的,该方法是按照高斯分布对每个像素建立模型, 并通过基于回归滤波的在线 EM 近似方法对模型参数进行更新,它能鲁棒地克服光照变化、 树枝摇动等造成的影响,但该方法也存在一些问题:1)该方法对运动物体在场景中停止不动或者长时间停止时检测失效,而且带有初始学习速度慢,在线更新费时、计算量大;2)无法完整准确地检测大并且运动缓慢的运动目标,运动目标的像素点不集中,只能检测到运动目标的部分轮廓,无法提取出目标对象的完整区域;3)无法将背景显露区域与运动目标区域很好地区分开;4)当运动目标由静止缓慢转化为运动时,易将背景显露区检测为前景,出现“影子”现象。
三、GMM 缺点及改进方法针对上述问题,一些科学研究者又在GMM 算法的基础上做了很多的改进:张、白等人[2]引入分块思想,把图像分为L*L 块;黄、胡等人[3]也引入了分块的思想,但是他们的分块理念是以当前像素点的8邻域作为一块;华、刘[4]把GMM 与改进的帧差法(相邻两帧图像对应像素点8邻域像素值相减之和)相结合,提高了计算效率;Suo 等人[5]是将混合高斯模型中的模型个数采改进为自适应的;刘等人[6]融合帧间差分法,检测背景显露区域和运动区域,很好的解决了问题4。
除此之外,还有基于纹理的混合高斯模型。
四、GMM 算法流程(1)用第一帧图像对高斯混合模型进行初始化 ()()0,,,0y x I y x =μ ① ()init std y x _,0=σ ②()init std init std y x __,20⨯=σ ③Mw 10=④ 一般模型的个数M 为3-6个,其中std_init 设置为20(2)对于t 时刻的像素()y x I t ,,分别与已经存在的M 个高斯模型依次进行匹配:()1,1,5.2|),(,--<-t i t i t y x y x I σμ ⑤(3)如果满足匹配条件,则该像素值与高斯模型匹配成功。
如果匹配不成功: a :当k<K 时,增加新的高斯模型; b :当k=K 时,用新高斯模型代替优先级σϖi最小的模型。
新的高斯模型,用当前像素值作为新模型的均值,即()y x I i ,=μ,协方差为init std i _=σ,权重为α=i w ,其中α为学习速率。
(4)未匹配模式的均值和方差不变,对匹配模式的第i 个高斯模型参数进行更新:()()y x I t t i t i ,11,,αμαμ+-=- ⑥ ()()()21,21,2,,1---+-=t i t t i t i y x I μασασ ⑦()αα+-=-1,,1t i t i w w ⑧(5)高斯模型参数更新完毕后,对每个像素点的K 歌高斯模型按优先级σϖi降序排序。
取前B 个高斯模型作为背景像素的最佳描述:15.0;min arg 1<<⎪⎭⎫⎝⎛>=∑=T T w B M k i k⑨(6)继续对()y x I t ,与上述B 个高斯模型进行匹配检验,如果()y x I t ,与前B 个高斯模型的任意一个匹配,则该像素点为背景点;否则为前景点。
(7)重复步骤(2)-(6),直到视频结束。
五、GMM 代码实现 #include<opencv.hpp> #include<highgui.h> #include<cv.h>using namespace cv; using namespace std;#define COMPONET 5 //混合高斯模型个数 #define ALPHA 0.03 //学习率 #define SD_INIT 6 //方差初始值 #define THRESHOLD 0.25 //前景所占比例 #define D 2.5int main() {CvCapture*capture=cvCreateFileCapture("E:\\project2\\videos\\video.avi");IplImage *frame, *grayFrame, *foreground, *background;int *foreg, *backg, *rank_index;double *weight, *mean, *sigma, *u_diff, *rank;double p = ALPHA / (1 / (double)COMPONET);double rank_temp = 0;int rank_index_temp = 0;CvRNG state; //随机生成状态器int match, height, width;frame = cvQueryFrame(capture);grayFrame = cvCreateImage(CvSize(frame->width, frame->height), IPL_DEPTH_8U, 1);foreground = cvCreateImage(CvSize(frame->width, frame->height), IPL_DEPTH_8U, 1);background = cvCreateImage(CvSize(frame->width, frame->height), IPL_DEPTH_8U, 1);height = grayFrame->height;width = grayFrame->widthStep;foreg = (int*)malloc(sizeof(int)*width*height);backg = (int*)malloc(sizeof(int)*width*height);rank = (double*)malloc(sizeof(double) * 1 * COMPONET); //优先级weight = (double*)malloc(sizeof(double)*width*height*COMPONET); //权重mean = (double *)malloc(sizeof(double)*width*height*COMPONET);//pixel meanssigma = (double *)malloc(sizeof(double)*width*height*COMPONET);//pixel standard deviationsu_diff = (double *)malloc(sizeof(double)*width*height*COMPONET);//difference of each pixel from mean//初始化均值、方差、权重for (int i = 0; i < height; i++){for (int j = 0; j < width; j++){for (int k = 0; k < COMPONET; k++){mean[i*width*COMPONET + j*COMPONET + k] = cvRandReal(&state) * 255;sigma[i*width*COMPONET + j*COMPONET + k] = SD_INIT;weight[i*width*COMPONET + j*COMPONET + k] = (double)1 / COMPONET;}}}while (1){rank_index = (int *)malloc(sizeof(int)*COMPONET);cvCvtColor(frame, grayFrame, CV_BGR2GRAY);// calculate difference of pixel values from meanfor (int i = 0; i < height; i++){for (int j = 0; j < width; j++){for (int k = 0; k < COMPONET; k++){u_diff[i*width*COMPONET + j*COMPONET + k] =abs((uchar)grayFrame->imageData[i*width + j] -mean[i*width*COMPONET + j*COMPONET + k]);}}}//update gaussian components for each pixelfor (int i = 0; i < height; i++){for (int j = 0; j < width; j++){match = 0;double sum_weight = 0;for (int k = 0; k < COMPONET; k++){if (u_diff[i*width*COMPONET + j*COMPONET + k]<= D*sigma[i*width*COMPONET + j*COMPONET + k]) //pixelmatches component{match = 1;// variable to signal component match//update weights, mean, sd, pweight[i*width*COMPONET + j*COMPONET + k] = (1 - ALPHA)*weight[i*width*COMPONET + j*COMPONET + k] + ALPHA;/*p = ALPHA / weight[i*width*COMPONET + j*COMPONET + k];mean[i*width*COMPONET + j*COMPONET + k] = (1 - p)*mean[i*width*COMPONET + j*COMPONET + k] + p*(uchar)grayFrame->imageData[i*width + j];sigma[i*width*COMPONET + j*COMPONET + k] = sqrt((1 - p)*(sigma[i*width*COMPONET + j*COMPONET + k] * sigma[i*width*COMPONET + j*COMPONET + k]) + p*(pow((uchar)grayFrame->imageData[i*width + j] - mean[i*width*COMPONET + j*COMPONET + k], 2)));*/mean[i*width*COMPONET + j*COMPONET + k] = (1 - ALPHA)*mean[i*width*COMPONET + j*COMPONET + k] + ALPHA*(uchar)grayFrame->imageData[i*width + j];sigma[i*width*COMPONET + j*COMPONET + k] = sqrt((1 - ALPHA)*(sigma[i*width*COMPONET + j*COMPONET + k] * sigma[i*width*COMPONET + j*COMPONET + k]) + ALPHA*(pow((uchar)grayFrame->imageData[i*width + j] -mean[i*width*COMPONET + j*COMPONET + k], 2)));}//else{// weight[i*width*COMPONET + j*COMPONET + k] =(1 - ALPHA)*weight[i*width*COMPONET + j*COMPONET + k]; // weight slighly decreases//}sum_weight += weight[i*width*COMPONET +j*COMPONET + k];}//权重归一化for (int k = 0; k < COMPONET; k++){weight[i*width*COMPONET + j*COMPONET + k] =weight[i*width*COMPONET + j*COMPONET + k] / sum_weight;}//获取权重最小下标double temp = weight[i*width*COMPONET +j*COMPONET];int min_index = 0;backg[i*width + j] = 0;for (int k = 0; k < COMPONET; k++){backg[i*width + j] = backg[i*width + j] + mean[i*width*COMPONET + j*COMPONET + k] * weight[i*width*COMPONET + j*COMPONET + k];if (weight[i*width*COMPONET + j*COMPONET + k] < temp){min_index = k;temp = weight[i*width*COMPONET + j*COMPONET + k];}rank_index[k] = k;}background->imageData[i*width + j] = (uchar)backg[i*width + j];//if no components match, create new componentif (match == 0){mean[i*width*COMPONET + j*COMPONET + min_index]= (uchar)grayFrame->imageData[i*width + j];sigma[i*width*COMPONET + j*COMPONET + min_index] = SD_INIT;weight[i*width*COMPONET + j*COMPONET + min_index] = 1 / COMPONET;}//计算优先级for (int k = 0; k < COMPONET; k++){rank[k] = weight[i*width*COMPONET + j*COMPONET + k] / sigma[i*width*COMPONET + j*COMPONET + k];}//sort rank valuesfor (int k = 1; k < COMPONET; k++){for (int m = 0; m < k; m++){if (rank[k] > rank[m]){//swap max valuesrank_temp = rank[m];rank[m] = rank[k];rank[k] = rank_temp;//swap max index valuesrank_index_temp = rank_index[m];rank_index[m] = rank_index[k];rank_index[k] = rank_index_temp;}}}//calculate foregroundmatch = 0;int b = 0;while ((match == 0) && (b < COMPONET)){if (weight[i*width*COMPONET + j*COMPONET + rank_index[b]] >= THRESHOLD){if (abs(u_diff[i*width*COMPONET + j*COMPONET + rank_index[b]]) <= D*sigma[i*width*COMPONET + j*COMPONET + rank_index[b]]){foreground->imageData[i*width + j] = 0;match = 1;}else{foreground->imageData[i*width + j] = (uchar)grayFrame->imageData[i*width + j];}}b++;}}}frame = cvQueryFrame(capture);cvShowImage("fore", foreground);cvShowImage("back", background);cvShowImage("frame", frame);char s = cvWaitKey(33);if (s == 27) break;free(rank_index);}return 0;}六、参考文献[1]Chris Stauffer,W.E.L Grimson.Adaptive background mixture models for real-time tracking[2]张燕平、白云球.应用改进混合高斯模型的运动目标检测[3]黄大卫、胡文翔。