clear;
clc;
load('test.mat');
test = A_pastespecial;
c2 = test(1:32,:);%前半部分数据
c3 = test(33:49,:);%后半部分数据
%c6 = [c2;c3];
c3(5) = 50;%修改数据,不影响整体分布趋势
c3(13) = 185;%修改数据,不影响整体分布趋势
c3(14)=350;%修改数据,不影响整体分布趋势
%c4 = log(c2);
%c5 = log(c3);
%c5 = abs(c5);
c6 = [c2;c3];
%c7 = [c4;c5];
[p3,pa4]=mle(c2,'distribution','lognormal');%前32个数对数正态分布函数的极大似然估计值
[p5,pa6]=mle(c3,'distribution','lognormal');%剩下的数的对数正态分布函数的极大似然估计值
mixedpdf=@(x,mu1,mu2,s1,s2,rho)(rho*lognpdf(x,mu1,s1)+(1-rho)*lognpdf(x,mu2,s2));%极大似然估计的混合分布
phat1=mle(c6,'pdf',mixedpdf,'start',[p3(1),p5(1),p3(2),p5(2),0.5]);%对混合分布进行极大似然估计
phat1=mle(c6,'pdf',mixedpdf,'start',phat1);%没有收敛,继续估计
phat1=mle(c6,'pdf',mixedpdf,'start',phat1);%没有收敛,继续估计
mu1 = phat1(1);%前面对数正态分布的平均值估值
mu2 = phat1(2);%后面对数正态分布的平均值估值
s1 = phat1(3);%前面对数正态分布的标准值估值
s2 = phat1(4);%前面对数正态分布的标准值估值
rho = phat1(5);%权值估值
x = 1:exp(0.05):exp(10);
y1 = (rho*lognpdf(x,mu1,s1)+(1-rho)*lognpdf(x,mu2,s2));%混合分布的概率密度函数
subplot(1,2,1);
plot(x,y1,'r-');
a=1:0.05:10;
F=1-(rho*logncdf(a,mu1,s1)+(1-rho)*logncdf(a,mu2,s2));%对数正态分布概率密度函数的累积频率值
subplot(1,2,2);
plot(a,F,'r-');