Графический выход плотности для функции gammamixEM (упаковка смесей)

Question

Я пользуюсь функцией gammamixEM из упаковки mixtools. Как вернуть графический выход плотности, как в функции normalmixEM (т. Е. Второй график в plot(...,which=2))?

Update:

Вот reproducible example для функции gammamixEM:

x <- c(rgamma(200, shape = 0.2, scale = 14), rgamma(200, 
    shape = 32, scale = 10), rgamma(200, shape = 5, scale = 6)) 
out <- gammamixEM(x, lambda = c(1, 1, 1)/3, verb = TRUE) 

Вот reproducible example для функции normalmixEM:

data(faithful) 
attach(faithful) 
out <- normalmixEM(waiting, arbvar = FALSE, epsilon = 1e-03) 
plot(out, which=2) 

Я хотел бы получить этот графический выход плотности из функции gammamixEM.

Answer 1

Здесь вы идете.

out <- normalmixEM(waiting, arbvar = FALSE, epsilon = 1e-03) 

x   <- out 
whichplots <- 2 
density = 2 %in% whichplots 
loglik = 1 %in% whichplots 

def.par <- par(ask=(loglik + density > 1), "mar") # only ask and mar are changed 
mix.object <- x 


k <- ncol(mix.object$posterior) 
x <- sort(mix.object$x) 
a <- hist(x, plot = FALSE) 
maxy <- max(max(a$density), .3989*mix.object$lambda/mix.object$sigma) 

Я просто должен был вырыть в исходный код plot.mixEM

Итак, теперь, чтобы сделать это с gammamixEM:

x <- c(rgamma(200, shape = 0.2, scale = 14), rgamma(200, 
                shape = 32, scale = 10), rgamma(200, shape = 5, scale = 6)) 
gammamixEM.out <- gammamixEM(x, lambda = c(1, 1, 1)/3, verb = TRUE) 



mix.object <- gammamixEM.out 

k <- ncol(mix.object$posterior) 
x <- sort(mix.object$x) 
a <- hist(x, plot = FALSE) 
maxy <- max(max(a$density), .3989*mix.object$lambda/mix.object$sigma) 

main2 <- "Density Curves" 
xlab2 <- "Data" 
col2 <- 2:(k+1) 

hist(x, prob = TRUE, main = main2, xlab = xlab2, 
    ylim = c(0,maxy)) 


for (i in 1:k) { 
    lines(x, mix.object$lambda[i] * 
      dnorm(x, 
       sd = sd(x))) 
} 

Я считаю, что это должно быть довольно просто вперед, чтобы продолжить этот пример немного, если вы хотите добавить ярлыки, плавные линии и т. д. Вот sourceplot.mixEM функция.