2、
（1）
clear;
load 2_data1;
hist(Y);
title(‘data 1’);
clear;
load 2_data2;
hist(Y);
title(‘data 2’);
（2）（3）
Data1
clear;
load 2_data1;
hist(Y);
title(‘data 1’);
alpha1=0.32;
alpha2=0.69;
mean1=0;
mean2=4.50;
var1=0.85;
var2=0.27; % initial guess
n=size(Y,2);
alpha1(1) = alpha1;
alpha2(1) = alpha2;
mean1(1)= mean1;
mean2(1)= mean2;
std1(1)= sqrt(var1);
std2(1)= sqrt(var2);
f1=exp(-(Y-mean1).^2/(2*var1))./sqrt(2*pi*var1);
f2=exp(-(Y-mean2).^2/(2*var2))./sqrt(2*pi*var2);
P1=(alpha1*f1)./(alpha1*f1+alpha2*f2);
P2=(alpha2*f2)./(alpha1*f1+alpha2*f2);
alpha1nth = sum(P1)/n;
alpha2nth = sum(P2)/n;
mean1nth = sum(Y.*P1)/sum(P1);
mean2nth = sum(Y.*P2)/sum(P2);
std1nth = sqrt(sum((Y-mean1nth).^2 .* P1) ./ sum(P1));
std2nth = sqrt(sum((Y-mean2nth).^2 .* P2) ./ sum(P2));
i = 1;
while (abs(alpha1(i)-alpha1nth) > 1e-6 || abs(std1(i)-std1nth) > 1e-6 || abs(mean1(i)-mean1nth) > 1e-6 || abs(mean2(i)-mean2nth) > 1e-6 || abs(alpha2(i)-alpha2nth) > 1e-6 || abs(std2(i)-std2nth) > 1e-6)
alpha1nth = alpha1(i);
mean1nth = mean1(i);
std1nth = std1(i);
alpha2nth = alpha2(i);
mean2nth = mean2(i);
std2nth = std2(i);
var1(i)=(std1(i)).^2;
var2(i)=(std2(i)).^2;
f1=exp(-(Y-mean1(i)).^2/(2*var1(i)))./sqrt(2*pi*var1(i));
f2=exp(-(Y-mean2(i)).^2/(2*var2(i)))./sqrt(2*pi*var2(i));
ll(i)=sum(log(alpha1(i)*f1+alpha2(i)*f2));
P1=(alpha1(i)*f1)./(alpha1(i)*f1+alpha2(i)*f2);
P2=(alpha2(i)*f2)./(alpha1(i)*f1+alpha2(i)*f2);
alpha1(i+1) = sum(P1)/n;
alpha2(i+1) = sum(P2)/n;
mean1(i+1) = sum(Y.*P1)/sum(P1);
mean2(i+1) = sum(Y.*P2)/sum(P2);
std1(i+1) = sqrt(sum((Y-mean1(i+1)).^2 .* P1) / sum(P1));
std2(i+1) = sqrt(sum((Y-mean2(i+1)).^2 .* P2) / sum(P2));
i = i+1;
end
plot(ll,’ro-‘);
title(‘data 1’);
alpha1nth
mean1nth
std1nth
alpha2nth
mean2nth
std2nth
>>
alpha1nth =
0.3241
mean1nth =
0.0269
std1nth =
0.9195
alpha2nth =
0.6759
mean2nth =
4.9724
std2nth =
0.5178
Data2
clear;
load 2_data2;
hist(Y);
title(‘data 2’);
alpha1=0.5;
alpha2=0.5;
mean1=1.8;
mean2=3.2;
var1=0.9;
var2=1.0; % initial guess
n=size(Y,2);
alpha1(1) = alpha1;
alpha2(1) = alpha2;
mean1(1)= mean1;
mean2(1)= mean2;
std1(1)= sqrt(var1);
std2(1)= sqrt(var2);
f1=exp(-(Y-mean1).^2/(2*var1))./sqrt(2*pi*var1);
f2=exp(-(Y-mean2).^2/(2*var2))./sqrt(2*pi*var2);
P1=(alpha1*f1)./(alpha1*f1+alpha2*f2);
P2=(alpha2*f2)./(alpha1*f1+alpha2*f2);
alpha1nth = sum(P1)/n;
alpha2nth = sum(P2)/n;
mean1nth = sum(Y.*P1)/sum(P1);
mean2nth = sum(Y.*P2)/sum(P2);
std1nth = sqrt(sum((Y-mean1nth).^2 .* P1) ./ sum(P1));
std2nth = sqrt(sum((Y-mean2nth).^2 .* P2) ./ sum(P2));
i = 1;
while (abs(alpha1(i)-alpha1nth) > 1e-4 || abs(std1(i)-std1nth) > 1e-4 || abs(mean1(i)-mean1nth) > 1e-4 || abs(mean2(i)-mean2nth) > 1e-4 || abs(alpha2(i)-alpha2nth) > 1e-4 || abs(std2(i)-std2nth) > 1e-4)
alpha1nth = alpha1(i);
mean1nth = mean1(i);
std1nth = std1(i);
alpha2nth = alpha2(i);
mean2nth = mean2(i);
std2nth = std2(i);
var1(i)=(std1(i)).^2;
var2(i)=(std2(i)).^2;
f1=exp(-(Y-mean1(i)).^2/(2*var1(i)))./sqrt(2*pi*var1(i));
f2=exp(-(Y-mean2(i)).^2/(2*var2(i)))./sqrt(2*pi*var2(i));
ll(i)=sum(log(alpha1(i)*f1+alpha2(i)*f2));
P1=(alpha1(i)*f1)./(alpha1(i)*f1+alpha2(i)*f2);
P2=(alpha2(i)*f2)./(alpha1(i)*f1+alpha2(i)*f2);
alpha1(i+1) = sum(P1)/n;
alpha2(i+1) = sum(P2)/n;
mean1(i+1) = sum(Y.*P1)/sum(P1);
mean2(i+1) = sum(Y.*P2)/sum(P2);
std1(i+1) = sqrt(sum((Y-mean1(i+1)).^2 .* P1) / sum(P1));
std2(i+1) = sqrt(sum((Y-mean2(i+1)).^2 .* P2) / sum(P2));
i = i+1;
end
plot(ll,’ro-‘)
title(‘data 2’)
alpha1nth
mean1nth
std1nth
alpha2nth
mean2nth
std2nth
>>
alpha1nth =
0.2266
mean1nth =
1.8647
std1nth =
0.9309
alpha2nth =
0.7734
mean2nth =
3.0173
std2nth =
0.5058