Сгенерируйте данные из смеси двух двухмерных Гауссовых распределений.

mu1 = [1 2];
Sigma1 = [2 0; 0 0.5];
mu2 = [-3 -5];
Sigma2 = [1 0;0 1];
rng(1); % For reproducibility
X = [mvnrnd(mu1,Sigma1,1000); mvnrnd(mu2,Sigma2,1000)];

Подбор смешанной гауссовской модели. Укажите, что существует два компонента.

GMModel = fitgmdist(X,2);

Постройте график данных по установленным контурам смешанной гауссовской модели.

figure
y = [zeros(1000,1);ones(1000,1)];
h = gscatter(X(:,1),X(:,2),y);
hold on
gmPDF = @(x,y) arrayfun(@(x0,y0) pdf(GMModel,[x0 y0]),x,y);
g = gca;
fcontour(gmPDF,[g.XLim g.YLim])
title('{\bf Scatter Plot and Fitted Gaussian Mixture Contours}')
legend(h,'Model 0','Model1')
hold off

Сгенерируйте данные из смеси двух двухмерных Гауссовых распределений. Создайте третий предиктор, который является суммой первого и второго предикторов.

mu1 = [1 2];
Sigma1 = [1 0; 0 1];
mu2 = [3 4];
Sigma2 = [0.5 0; 0 0.5];
rng(3); % For reproducibility
X1 = [mvnrnd(mu1,Sigma1,100);mvnrnd(mu2,Sigma2,100)];
X = [X1,X1(:,1)+X1(:,2)];

Столбцы X являются линейно зависимыми. Это может вызвать плохо обусловленные ковариационные оценки.

Подбор Смешанной гауссовской модели к данным. Можно использовать try / catch операторы, которые помогают управлять сообщениями об ошибке.

rng(1); % Reset seed for common start values
try
    GMModel = fitgmdist(X,2)
catch exception
    disp('There was an error fitting the Gaussian mixture model')
    error = exception.message
end

There was an error fitting the Gaussian mixture model

error = 
'Ill-conditioned covariance created at iteration 3.'

Ковариационные оценки являются плохо обусловленными. Следовательно, оптимизация останавливается, и появляется ошибка.

Снова подгоняйте Смешанную гауссовскую модель, но используйте регуляризацию.

rng(3); % Reset seed for common start values
GMModel = fitgmdist(X,2,'RegularizationValue',0.1)

GMModel = 

Gaussian mixture distribution with 2 components in 3 dimensions
Component 1:
Mixing proportion: 0.536725
Mean:    2.8831    3.9506    6.8338

Component 2:
Mixing proportion: 0.463275
Mean:    0.8813    1.9758    2.8571

В этом случае алгоритм сходится к решению из-за регуляризации.

Смешанные гауссовские модели перед подгонкой к данным необходимо задать несколько компонентов. Для многих приложений может оказаться трудным узнать соответствующее количество компонентов. В этом примере показано, как исследовать данные и попытаться получить начальное предположение о количестве компонентов, используя анализ основных компонентов.

Загрузите набор данных радужки Фишера.

load fisheriris
classes = unique(species)

classes = 3x1 cell
    {'setosa'    }
    {'versicolor'}
    {'virginica' }

Набор данных содержит три класса видов радужки. Анализ выполняется так, как будто это неизвестно.

Используйте анализ основного компонента, чтобы уменьшить размерность данных до двух измерений для визуализации.

[~,score] = pca(meas,'NumComponents',2);

Подгонка трёх смешанных гауссовских моделей к данным путем определения 1, 2 и 3 компонентов. Увеличьте количество итераций оптимизации до 1000. Используйте запись через точку для хранения окончательных оценок параметров. По умолчанию программное обеспечение соответствует полным и различным ковариациям для каждого компонента.

GMModels = cell(3,1); % Preallocation
options = statset('MaxIter',1000);
rng(1); % For reproducibility

for j = 1:3
    GMModels{j} = fitgmdist(score,j,'Options',options);
    fprintf('\n GM Mean for %i Component(s)\n',j)
    Mu = GMModels{j}.mu
end

 GM Mean for 1 Component(s)

Mu = 1×2
10-14 ×

    0.1044   -0.0429

 GM Mean for 2 Component(s)

Mu = 2×2

    1.3212   -0.0954
   -2.6424    0.1909

 GM Mean for 3 Component(s)

Mu = 3×2

    0.4856   -0.1287
    1.4484   -0.0904
   -2.6424    0.1909

GMModels - массив ячеек, содержащий три, установленные gmdistribution модели. Средства в трёх моделях компонента различны, что позволяет предположить, что модель различает три вида радужки.

Постройте графики счетов по подобранным контурам смешанной гауссовской модели. Поскольку набор данных включает метки, используйте gscatter чтобы различить истинное количество компонентов.

figure
for j = 1:3
    subplot(2,2,j)
    h1 = gscatter(score(:,1),score(:,2),species);
    h = gca;
    hold on
    gmPDF = @(x,y) arrayfun(@(x0,y0) pdf(GMModels{j},[x0 y0]),x,y);
    fcontour(gmPDF,[h.XLim h.YLim],'MeshDensity',100)
    title(sprintf('GM Model - %i Component(s)',j));
    xlabel('1st principal component');
    ylabel('2nd principal component');
    if(j ~= 3)
        legend off;
    end
    hold off
end
g = legend(h1);
g.Position = [0.7 0.25 0.1 0.1];

Трехкомпонентная смешанная гауссовская модель в сочетании с РСА выглядит так, как будто она различает три вида радужной оболочки глаза.

Существуют другие опции, которые можно использовать, чтобы помочь выбрать соответствующее количество компонентов для смешанной гауссовской модели. Для примера,

Смешанные гауссовские модели перед подгонкой к данным необходимо задать несколько компонентов. Для многих приложений может оказаться трудным узнать соответствующее количество компонентов. Этот пример использует статистику АПК-подгонки, чтобы помочь вам выбрать наиболее подходящую Смешанную гауссовскую модель для меняющегося количества компонентов.

Сгенерируйте данные из смеси двух двухмерных Гауссовых распределений.

mu1 = [1 1];
Sigma1 = [0.5 0; 0 0.5];
mu2 = [2 4];
Sigma2 = [0.2 0; 0 0.2];
rng(1);
X = [mvnrnd(mu1,Sigma1,1000);mvnrnd(mu2,Sigma2,1000)];

plot(X(:,1),X(:,2),'ko')
title('Scatter Plot')
xlim([min(X(:)) max(X(:))]) % Make axes have the same scale
ylim([min(X(:)) max(X(:))])

Предположим, что вы не знаете базовых значений параметров, графики поля точек предлагают:

Подбор двухкомпонентной Смешанной гауссовской модели. Основываясь на проверке графика поля точек, задайте, что ковариационные матрицы являются диагональными. Печать окончательной статистики итерации и логарифмической правдоподобности в Командное окно путем передачи statset структура как значение Options аргумент пары "имя-значение".

options = statset('Display','final');
GMModel = fitgmdist(X,2,'CovarianceType','diagonal','Options',options);

11 iterations, log-likelihood = -4787.38

GMModel является установленным gmdistribution модель.

Исследуйте AIC по различному количеству компонентов.

AIC = zeros(1,4);
GMModels = cell(1,4);
options = statset('MaxIter',500);
for k = 1:4
    GMModels{k} = fitgmdist(X,k,'Options',options,'CovarianceType','diagonal');
    AIC(k)= GMModels{k}.AIC;
end

[minAIC,numComponents] = min(AIC);
numComponents

numComponents = 2

BestModel = GMModels{numComponents}

BestModel = 

Gaussian mixture distribution with 2 components in 2 dimensions
Component 1:
Mixing proportion: 0.501719
Mean:    1.9824    4.0013

Component 2:
Mixing proportion: 0.498281
Mean:    0.9880    1.0511

Наименьший AIC происходит, когда программное обеспечение подходит для двухкомпонентных Смешанных гауссовских моделей.

Смешанные гауссовские модели оценки параметров могут варьироваться при различных начальных значениях. В этом примере показов, как управлять начальными значениями, когда вы подгоняете Смешанные гауссовские модели используя fitgmdist.

Загрузите набор данных радужки Фишера. Используйте длины и ширины лепестков в качестве предикторов.

load fisheriris
X = meas(:,3:4);

Подбор Смешанной гауссовской модели к данным с помощью начальных значений по умолчанию. Существует три вида радужки, поэтому задайте k = 3 компонента.

rng(10); % For reproducibility
GMModel1 = fitgmdist(X,3);

По умолчанию программное обеспечение:

Подгонка Смешанной гауссовской модели путем соединения каждого наблюдения с его меткой.

y = ones(size(X,1),1);
y(strcmp(species,'setosa')) = 2;
y(strcmp(species,'virginica')) = 3;

GMModel2 = fitgmdist(X,3,'Start',y);

Подгонка смешанной гауссовской модели путем явного определения начальных средств, ковариационных матриц и пропорций смешения.

Mu = [1 1; 2 2; 3 3];
Sigma(:,:,1) = [1 1; 1 2];
Sigma(:,:,2) = 2*[1 1; 1 2];
Sigma(:,:,3) = 3*[1 1; 1 2];
PComponents = [1/2,1/4,1/4];
S = struct('mu',Mu,'Sigma',Sigma,'ComponentProportion',PComponents);

GMModel3 = fitgmdist(X,3,'Start',S);

Использование gscatter чтобы построить диаграмму рассеяния, которая различает виды радужки. Для каждой модели постройте график подобранных контуров смешанной гауссовской модели.

figure
subplot(2,2,1)
h = gscatter(X(:,1),X(:,2),species,[],'o',4);
haxis = gca;
xlim = haxis.XLim;
ylim = haxis.YLim;
d = (max([xlim ylim])-min([xlim ylim]))/1000;
[X1Grid,X2Grid] = meshgrid(xlim(1):d:xlim(2),ylim(1):d:ylim(2));
hold on
contour(X1Grid,X2Grid,reshape(pdf(GMModel1,[X1Grid(:) X2Grid(:)]),...
    size(X1Grid,1),size(X1Grid,2)),20)
uistack(h,'top')
title('{\bf Random Initial Values}');
xlabel('Sepal length');
ylabel('Sepal width');
legend off;
hold off
subplot(2,2,2)
h = gscatter(X(:,1),X(:,2),species,[],'o',4);
hold on
contour(X1Grid,X2Grid,reshape(pdf(GMModel2,[X1Grid(:) X2Grid(:)]),...
    size(X1Grid,1),size(X1Grid,2)),20)
uistack(h,'top')
title('{\bf Initial Values from Labels}');
xlabel('Sepal length');
ylabel('Sepal width');
legend off
hold off
subplot(2,2,3)
h = gscatter(X(:,1),X(:,2),species,[],'o',4);
hold on
contour(X1Grid,X2Grid,reshape(pdf(GMModel3,[X1Grid(:) X2Grid(:)]),...
    size(X1Grid,1),size(X1Grid,2)),20)
uistack(h,'top')
title('{\bf Initial Values from the Structure}');
xlabel('Sepal length');
ylabel('Sepal width');
legend('Location',[0.7,0.25,0.1,0.1]);
hold off

По контурам GMModel2 по-видимому, предполагает небольшую тримодальность, в то время как другие предполагают бимодальные распределения.

Отобразите предполагаемые средства компонента.

table(GMModel1.mu,GMModel2.mu,GMModel3.mu,'VariableNames',...
    {'Model1','Model2','Model3'})

ans=3×3 table
         Model1               Model2              Model3     
    _________________    ________________    ________________

    5.2115     2.0119    4.2857    1.3339    1.4604    0.2429
     1.461    0.24423     1.462     0.246    4.7509    1.4629
    4.6829     1.4429    5.5507    2.0316    5.0158    1.8592

GMModel2 кажется, различают виды радужки лучшими.