R наложение двумерной нормальной плотности (эллипсы) на диаграмму рассеяния

Подобные вопросы есть на сайте, но я не смог найти ответ на эту, казалось бы, очень простую проблему. Я подгоняю смесь двух гауссов к набору данных Old Faithful:

if(!require("mixtools")) { install.packages("mixtools");  require("mixtools") }
data_f <- faithful
plot(data_f$waiting, data_f$eruptions)
data_f.k2 = mvnormalmixEM(as.matrix(data_f), k=2, maxit=100, epsilon=0.01) 
data_f.k2$mu # estimated mean coordinates for the 2 multivariate Gaussians
data_f.k2$sigma # estimated covariance matrix

Я просто хочу наложить друг на друга два эллипса для двух гауссовских компонентов модели, описываемой средними векторами data_f.k2$mu и ковариационными матрицами data_f.k2$sigma. Чтобы получить что-то вроде:

введите здесь описание изображения

Кому интересно, здесь — это решение MatLab, которое создало приведенный выше график.


r mean plot covariance normal-distribution
person Zhubarb    schedule 04.02.2015    source источник
comment
Я еще не использовал его, но, возможно, пакет ellipse cran.r- project.org/web/packages/ellipse/index.html ? различные процедуры для рисования эллипсов и эллипсоподобных доверительных областей   -  person Julian    schedule 04.02.2015

Ответы (3)


arrow_upward
2
arrow_downward

Вы можете использовать функцию ellipse из пакета mixtools. Первоначальная проблема заключалась в том, что эта функция меняет местами x и y с вашего графика. Я постараюсь понять это и обновить ответ. (Я оставлю цвета кому-то другому...)

plot( data_f$eruptions,data_f$waiting)
for (i in 1: length(data_f.k2$mu))  ellipse(data_f.k2$mu[[i]],data_f.k2$sigma[[i]])
person Julian    schedule 04.02.2015
comment
извините, код работает, но я не вижу многоточия на графике. Вы не против опубликовать получившуюся цифру? - person Zhubarb; 04.02.2015

arrow_upward
3
arrow_downward

Если вас также интересуют цвета, вы можете использовать апостериор, чтобы получить соответствующие группы. Я сделал это с ggplot2, но сначала покажу цветное решение, используя код @Julian.

# group data for coloring
data_f$group <- factor(apply(data_f.k2$posterior, 1, which.max))
# plotting
plot(data_f$eruptions, data_f$waiting, col = data_f$group)
for (i in 1: length(data_f.k2$mu))  ellipse(data_f.k2$mu[[i]],data_f.k2$sigma[[i]], col=i)

И для моей версии с использованием ggplot2.

# needs ggplot2 package
require("ggplot2")
# ellipsis data 
ell <- cbind(data.frame(group=factor(rep(1:length(data_f.k2$mu), each=250))), 
             do.call(rbind, mapply(ellipse, data_f.k2$mu, data_f.k2$sigma, 
                                   npoints=250, SIMPLIFY=FALSE)))

# plotting command
p <- ggplot(data_f, aes(color=group)) + 
  geom_point(aes(waiting, eruptions)) +
  geom_path(data=ell, aes(x=`2`, y=`1`)) +
  theme_bw(base_size=16)
print(p)
person shadow    schedule 04.02.2015

arrow_upward
0
arrow_downward

Использование внутренней функции построения графиков mixtools:

plot.mixEM(data_f.k2, whichplots=2)

person Luca Pandolfini    schedule 26.11.2018
comment
Не могли бы вы добавить какое-либо объяснение, пожалуйста? - person rudolf_franek; 26.11.2018
comment
Добро пожаловать в stackoverflow.com. При ответе на вопрос, пожалуйста, дайте полное объяснение. Ваш нынешний пост слишком короток, чтобы быть полезным. - person Simon.S.A.; 26.11.2018