R語言代寫實現有限混合模型建模分析

原文連接：http://tecdat.cn/?p=6129

 

介紹

有限混合模型在應用於數據時很是有用，其中觀察來自不一樣的羣體，而且羣體隸屬關係未知。  

模擬數據

首先，咱們將模擬一些數據。讓咱們模擬兩個正態分佈 - 一個平均值爲0，另外一個平均值爲50，二者的標準差爲5。

m1 <- 0
m2 <- 50
sd1 <- sd2 <- 5
N1 <- 100
N2 <- 10

a <- rnorm(n=N1, mean=m1, sd=sd1)
b <- rnorm(n=N2, mean=m2, sd=sd2)

如今讓咱們將數據「混合」在一塊兒......

​  

 

 

 

print(table(clusters(flexfit), data$class))
##    
##       1   2
##   1 100   0
##   2   0  10

參數怎麼樣？

cat('pred:', c1[1], '\n')
cat('true:', m1, '\n\n')
cat('pred:', c1[2], '\n')
cat('true:', sd1, '\n\n')

cat('pred:', c2[1], '\n')
cat('true:', m2, '\n\n')
cat('pred:', c2[2], '\n')
cat('true:', sd2, '\n\n')
## pred: -0.5613484 
## true: 0 
## 
## pred: 4.799484 
## true: 5 
## 
## pred: 52.86911 
## true: 50 
## 
## pred: 6.89413 
## true: 5

讓咱們可視化真實數據和咱們擬合的混合模型。

ggplot(data) +
geom_histogram(aes(x, ..density..), binwidth = 1, colour = "black", fill = "white") +
stat_function(geom = "line", fun = plot_mix_comps,
                args = list(c1[1], c1[2], lam[1]/sum(lam)),
 stat_function(geom = "line", fun = plot_mix_comps,
                args = list(c2[1], c2[2], lam[2]/sum(lam)),
                colour = "blue", lwd = 1.5) +
ylab("Density")

​

看起來咱們作得很好！ 

​

 ​

 

​

例子

如今，讓咱們考慮一個花瓣寬度爲鳶尾花的真實例子。 

p <- ggplot(iris, aes(x = Petal.Width)) + 
  geom_histogram(aes(x = Petal.Width, ..density..), binwidth = 0.1, colour = "black", fill = "white")
p

 

​ 

 

flexfit <- flexmix(Petal.Width ~ 1, data = iris, k = 3, model = list(mo1, mo2, mo3))

print(table(clusters(flexfit), iris$Species))
##    
##     setosa versicolor virginica
##   1      0          2        46
##   2      0         48         4
##   3     50          0         0

geom_histogram(aes(x = Petal.Width, ..density..), binwidth = 0.1, colour = "black", fill = "white") +
                 args = list(c1[1], c1[2], lam[1]/sum(lam)),
                colour = "red", lwd = 1.5) +
stat_function(geom = "line", fun = plot_mix_comps,
                args = list(c2[1], c2[2], lam[2]/sum(lam)),
 stat_function(geom = "line", fun = plot_mix_comps,
                args = list(c3[1], c3[2], lam[3]/sum(lam)),
                colour = "green", lwd = 1.5) +
ylab("Density")

 

​

即便咱們不知道潛在的物種分配，咱們也可以對花瓣寬度的基本分佈作出某些陳述 。

 

 

若是您有任何疑問，請在下面發表評論。   

大數據部落 -中國專業的第三方數據服務提供商，提供定製化的一站式數據挖掘和統計分析諮詢服務

統計分析和數據挖掘諮詢服務：y0.cn/teradat（諮詢服務請聯繫官網客服）

​QQ：3025393450

 

​QQ交流羣：186388004 

【服務場景】  

科研項目; 公司項目外包;線上線下一對一培訓;數據爬蟲採集;學術研究;報告撰寫;市場調查。

【大數據部落】提供定製化的一站式數據挖掘和統計分析諮詢

歡迎選修咱們的R語言數據分析挖掘必知必會課程！

 

 

歡迎關注 微信公衆號，瞭解更多數據乾貨資訊！