我在 R 中使用一個名為 BetaMixture 的包來擬合資料向量的 beta 分布的混合。輸出被提供給 a hist，它產生一個帶有混合模型組件的良好直方圖：
# Install and load the libraries
#install.packages("BetaModels")
library(BetaModels)
# Create a vector, fit mixture models and plot the histogram
vec <- c(rbeta(700, 5, 2), rbeta(300, 1, 10))
model <- BetaMixture(vec,2)
h <- hist(model, breaks = 35)
到現在為止還挺好。現在我如何在 ggplot 中得到這個？我檢查了h物件，但這與物件沒有什么不同model。它們完全一樣。我不知道這hist甚至對這門課是如何作業的。除了 ，它從 中提取什么model來生成這個圖@datavec？
uj5u.com熱心網友回復：
你可以得到hist的功能BetaMixed使用物件getMethod("hist", "BetaMixture")。
您可以在下面找到此功能到“ggplot2世界”的簡單翻譯。
myhist <- function (x, ...) {
.local <- function (x, mixcols = 1:7, breaks=25, ...)
{
df1 <- data.frame(x=x@datavec)
p <- ggplot(data=df1, aes(x=x))
geom_histogram(aes(y=..density..), bins=breaks, alpha=0.5, fill="gray50", color="black")
while (length(mixcols) < ncol(x@mle)) mixcols <- c(mixcols,
mixcols)
xv <- seq(0, 1, length = 502)[1:501]
for (J in 1:ncol(x@mle)) {
y <- x@phi[J] * dbeta(xv, x@mle[1, J], x@mle[2, J])
df2 <- data.frame(xv, y)
p <- p geom_line(data=df2, aes(xv, y), size=1, col=mixcols[J])
}
p <- p theme_bw()
invisible(p)
}
.local(x, ...)
}
library(ggplot2)
# Now p is a ggplot2 object.
p <- myhist(model, breaks=35)
print(p)
uj5u.com熱心網友回復：
回傳的物件BetaMixture是一個 S4 類物件，有 2 個感興趣的插槽。
SlotZ回傳屬于每個分布的每個資料點的概率矩陣。
所以在前 6 行中，所有點都屬于第二分布。
head(model@Z)
# [,1] [,2]
#[1,] 1.354527e-04 0.9998645
#[2,] 4.463074e-03 0.9955369
#[3,] 1.551999e-03 0.9984480
#[4,] 1.642579e-03 0.9983574
#[5,] 1.437047e-09 1.0000000
#[6,] 9.911427e-04 0.9990089
槽mle回傳引數的最大似然估計。
現在使用這些值來創建向量的 data.frame 和引數的 data.frame。
df1 <- data.frame(vec)
df1$component <- factor(apply(model@Z, 1, which.max))
colors <- as.integer(levels(df1$component))
params <- as.data.frame(t(model@mle))
names(params) <- c("shape1", "shape2")
繪制資料。
library(ggplot2)
g <- ggplot(df1, aes(x = vec, group = component))
geom_histogram(aes(y = ..density..),
bins = 35, fill = "grey", color = "grey40")
for(i in 1:nrow(params)){
sh1 <- params$shape1[i]
sh2 <- params$shape2[i]
g <- g stat_function(
fun = dbeta,
args = list(shape1 = sh1, shape2 = sh2),
color = colors[i]
)
}
suppressWarnings(print(g theme_bw()))
轉載請註明出處，本文鏈接：https://www.uj5u.com/net/373010.html
上一篇：R中的資料集操作