close?all;clear;clc;
%???參考書籍Pattern.Recognition.and.Machine.Learning.pdf
%???http://www.pr-ml.cn
%???lwm@pr-ml.cn
%???2009/10/15
%%?
M=2;??????????%?number?of?Gaussian
N=200;????????%?total?number?of?data?samples
th=0.000001;??%?convergent?threshold
K=2;??????????%?demention?of?output?signal

%?待生成數(shù)據(jù)的參數(shù)
a_real?=[4/5;1/5];
mu_real=[3?4;
?????????5?3];
cov_real（::1）=[5?0;
?????????????????0?0.2];
cov_real（::2）=[0.1?0;
?????????????????0?0.1];
%?generate?the?data
x=[?mvnrnd（?mu_real（:1）??cov_real（::1）??round（N*a_real（1））?）‘????mvnrnd（mu_real（:2）cov_real（::2）N-round（N*a_real（1）））‘];

%?for?i=1:round（N*a_real（1））
%?????while?（~（（x（1i）>0）&&（x（2i）>0）&&（x（1i）<10）&&（x（2i）<10）））
%?????????x（:i）=mvnrnd（mu_real（:1）cov_real（::1）1）‘;
%?????end
%?end
%?
%?for?i=round（N*a_real（1））+1:N
%?????while?（~（（x（1i）>0）&&（x（2i）>0）&&（x（1i）<10）&&（x（2i）<10）））
%?????????x（:i）=mvnrnd（mu_real（:1）cov_real（::1）1）‘;
%?????end
%?end

figure（1）plot（x（1:）x（2:）‘.‘）

%這里生成的數(shù)據(jù)全部符合標(biāo)準(zhǔn)
%%?%%%%%%%%%%%%%%%%?參數(shù)初始化
a=[1/32/3];
mu=[1?2;2?1];%均值初始化完畢
cov（::1）=[1?0;
????????????0?1];
cov（::2）=[1?0;
????????????0?1];%協(xié)方差初始化
%%?EM?Algorothm
%?loop
count=0;
figure（2）hold?on
while?1
????a_old??=?a;
????mu_old?=?mu;
????cov_old=?cov;??????
????rznk_p=zeros（MN）;
????for?cm=1:M
????????mu_cm=mu（:cm）;
????????cov_cm=cov（::cm）;
????????for?cn=1:N
????????????p_cm=exp（-0.5*（x（:cn）-mu_cm）‘/cov_cm*（x（:cn）-mu_cm））;
????????????rznk_p（cmcn）=p_cm;
????????end
????????rznk_p（cm:）=rznk_p（cm:）/sqrt（det（cov_cm））;
????end
????rznk_p=rznk_p*（2*pi）^（-K/2）;
%E?step
????%開始求rznk
????rznk=zeros（MN）;%r（Z
????pikn=zeros（1M）;%r（Z
????pikn_sum=0;
????for?cn=1:N
????????for?cm=1:M
????????????pikn（1cm）=a（cm）*rznk_p（cmcn）;
%???????????pikn_sum=pikn_sum+pikn（1cm）;
????????end
????????for?cm=1:M
????????????rznk（cmcn）=pikn（1cm）/sum（pikn）;
????????end
????end
????????%求rank結(jié)束
%?M?step
????nk=zeros（1M）;
????for?cm=1:M
????????for?cn=1:N
????????????nk（1cm）=nk（1cm）+rznk（cmcn）;
????????end
????end
????a=nk/N;
????rznk_sum_mu=zeros（M1）;
????????
????%?求均值MU
????for?cm=1:M
????????rznk_sum_mu=0;%開始的時(shí)候就是錯(cuò)在這里，這里要置零。
????????for?cn=1:N
????????????rznk_sum_mu=rznk_sum_mu+rznk（cmcn）*x（:cn）;
????????end
????????mu（:cm）=rznk_sum_mu/nk（cm）;
????end
????
????%?求協(xié)方差COV???
????for?cm=1:M
????????rznk_sum_cov=zeros（KM）;
????????for?cn=1:N
????????????rznk_sum_cov=rznk_sum_cov+rznk（cmcn）*（x（:cn）-mu（:cm））*（x（:cn）-mu（:cm））‘;
????????end
????????cov（::cm）=rznk_sum_cov/nk（cm）;
????end
?????????
????t=max（[norm（a_old（:）-a（:））/norm（a_old（:））;norm（mu_old（:）-mu（:））/norm（mu_old（:））;norm（cov_old（:）-cov（:））/norm（cov_old（:））]）;
?
????temp_f=sum（log（sum（pikn）））;
????plot（counttemp_f‘r+‘）
????count=count+1;????????
????????
????if?t????????break;
????end
?
end??%while?1
???
hold?off
f=sum（log（sum（pikn）））;
????
a_best=a;
mu_best=mu;
cov_best=cov;
f_best=f;

?屬性????????????大小?????日期????時(shí)間???名稱
-----------?---------??----------?-----??----

?????文件????9028817??2009-12-08?18:43??PatternRecognitionAndMachineLearning.pdf

?????文件???????3765??2009-12-09?15:54??em_my3.m

?????文件???????3329??2009-12-08?22:43??em_my1.m

?????文件???????2927??2009-12-08?23:19??em_my2.m

-----------?---------??----------?-----??----

??????????????9038838????????????????????4