J. Estermont
J. Estermont

Reputation: 31

R: plot of normalmixEM shows truncated density plots (mixtools)

I'm currently trying to plot the components found via EM algorithm. However, the estimated densities do not extend fully to the end. It looks like this:enter image description here

My code is:

plot(EM_data, which=2, xlim= c(0, 80), xlab2= "", yaxt= "n", main2 ="", lwd2=0.8, border = "azure3")
lines(density(EM_data), lty=2, lwd=0.8)

The plot is truncated wether I specify xlim or not. xlim2 is not defined for this type of plot. Where am I going wrong?

Upvotes: 1

Views: 984

Answers (1)

The method to plot mixEM only draws within the range of the data, if you want to extend the densities you must build your own function.

Use something like this:

Example data:

library(mixtools)

data(faithful)
attach(faithful)
set.seed(100)
EM_data<-normalmixEM(waiting, arbvar = FALSE, epsilon = 1e-03)

mixtools plot:

plot(EM_data, which=2, xlim= c(30, 110), xlab2= "", yaxt= "n", main2 ="", 
     lwd2=0.8, border = "azure3")
lines(density(EM_data$x), lty=2, lwd=0.8)

enter image description here

Adaptation by extending densities:

a <- hist(EM_data$x, plot = FALSE)
maxy <- max(max(a$density), 0.3989 * EM_data$lambda/EM_data$sigma)
hist(EM_data$x, prob = TRUE, main = "", xlab = "", xlim= c(30, 110), 
     ylim = c(0, maxy), yaxt= "n", border = "azure3")
for (i in 1:ncol(EM_data$posterior)) {
  curve(EM_data$lambda[i] * dnorm(x, mean = EM_data$mu[i], sd = EM_data$sigma[i]), 
        col = 1 + i, lwd = 0.8, add = TRUE)
}
lines(density(EM_data$x), lty=2, lwd=0.8)

enter image description here

Upvotes: 2

Related Questions