gaussian_mixed.Rmd

---
title: "gaussian mixture model"
author: "liuc"
date: "11/11/2021"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

## 高斯混合模型

> https://towardsdatascience.com/mixture-modelling-from-scratch-in-r-5ab7bfc83eef


高斯混合模型是非监督的聚类方法，其和kmeans都是利用EM(Expectation-Maximization) algorithm 进行实现的。

高斯混合模型 (Gaussian Mixture Model，GMM) 是一种概率模型，用于对由多个高斯分布组成的数据进行建模。它假设数据是由 K 个高斯分布混合而成的，每个高斯分布被称为一个组件，每个组件都有自己的均值和方差。GMM 的目标是通过最大化似然函数来估计模型参数，即均值、方差和混合系数，使得模型能够最好地拟合数据。

GMM 可以用于多种任务，例如聚类、异常检测、密度估计等。在聚类任务中，GMM 可以将数据点分配到不同的组件中，从而实现聚类。在异常检测任务中，如果某个数据点的概率密度低于一个阈值，则认为该数据点是异常的。在密度估计任务中，GMM 可以估计数据的概率密度函数。

GMM 的主要缺点是它对初始值敏感，如果初始值选的不好，则有可能会收敛到局部最优解。为了克服这个问题，可以使用多个初始值运行算法，并选择具有最大似然的模型。

```{r}
library('easystats')
library(mclust)
# library('mixtools')
# library('ClusterR')

```

### 首先利用stats包自带的kmeans & mclust::Mclust 进行聚类

```{r}
X <- iris[,1:4]
y <- iris$Species

set.seed(42)
km_res <- kmeans(X, 3)
km_res
gmm.mclust <- Mclust(X, 3)
summary(gmm.mclust)

table(y, km_res$cluster)
table(y, gmm.mclust$classification) # 结果可以看到GMM 对iris数据集的聚类效果更好一些
```

### 如何确定最佳的K值

对于partitioning clustering，其不同于层次聚类，往往需要一个K值。

```{r}
# mlust BIC 确认最佳的K值
BIC <- mclust::mclustBIC(X)
summary(BIC) # 结果显示K=2 为最佳，和实际出入较大
plot(BIC)

mod1 <- Mclust(X, x = BIC)
summary(mod1, parameters = TRUE)

plot(mod1, what = "classification")


```

```{r}
# ICL 推测最佳K值
ICL <- mclustICL(X)
summary(ICL)

plot(ICL)

#
d_clust <- Mclust(as.matrix(X), G=1:15, 
                  modelNames = mclust.options("emModelNames"))
d_clust$BIC
plot(d_clust)
```


继续最佳K值的确认，采用elbow法

```{r}
# the elbow method
k.max <- 15
wss <- sapply(1:k.max, 
              function(k){kmeans(X, k)$tot.withinss})
wss
plot(1:k.max, wss,
     type="b", pch = 19, frame = FALSE, 
     xlab="Number of clusters K",
     ylab="Total within-clusters sum of squares") # 也显示K为2为好。。这个。。。

# The Silhouette Method

```


### use clustR package for GMM

```{r}
gmm = clustR::GMM(X, 2, "maha_dist", "random_subset", 10, 10)

summary(gmm)

opt_gmm = Optimal_Clusters_GMM(X, max_clusters = 10, criterion = "BIC", 
                               dist_mode = "maha_dist", seed_mode = "random_subset",
                               km_iter = 10, em_iter = 10, var_floor = 1e-10,
                               plot_data = T)
```