R выделения 3 компонента линейной регрессии из смеси

Question

Я создал dataframe dfall, содержащую 3 линейных компоненты регрессии

x1=runif(50, min=0, max=100) 
e1=runif(50, min=0, max=10) 
y1 <- 0.2*x1+10+e1 
y1 


plot(x1,y1,col="blue") 


df1 <- data.frame(x=x1,y=y1) 
df1$ID <- 1 
df1$col <- "red" 


x2=runif(25, min=0, max=100) 
e2=runif(25, min=0, max=5) 
y2 <- 0.7*x2+15+e2 
y2 


plot(x2,y2,col="blue") 


df2 <- data.frame(x=x2,y=y2) 
df2$ID <- 2 
df2$col <- "green" 


x3=runif(35, min=0, max=100) 
e3=runif(35, min=0, max=15) 
y3 <- -0.5*x3+30+e3 
y3 


plot(x3,y3,col="blue") 


df3 <- data.frame(x=x3,y=y3) 
df3$ID <- 3 
df3$col <- "blue" 


dfall <-rbind(df1,df2,df3) 
dfall 
dfall <- dfall[sample(1:nrow(dfall)), ] 
dfall 

plot(dfall$x,dfall$y,col=dfall$col) 

Затем я попытался отделить линейные компоненты регрессии с помощью kmeans:

fitkm <- kmeans(dfall[,c(1:2)], 3) 

dfall <- data.frame(dfall, km=fitkm$cluster) 
dfall 

Однако Я получил очень плохие результаты классификации:

table(dfall$ID,dfall$km) 

Есть ли более способ точно отделить 3 линейные регрессионные компоненты? Благодарим вас за помощь.

Answer 1

Учитывая ваш пример, вы можете захотеть изучить модели конечных смесей, которые позволят вам восстановить параметры в ваших основных регрессиях и классификациях. Вот пример с вашими данными:

library(mixtools) 
mixmod <- regmixEM(dfall$y, dfall$x, k=3) 
summary(mixmod) 

Выход дает пропорции каждого вида наблюдений и коэффициенты - лямбда смесь пропорции, бета1 перехват и бета2 коэффициентов. Матч с вашей моделируемой данных довольно хорошо:

summary of regmixEM object: 
      comp 1 comp 2 comp 3 
lambda 0.315816 0.457191 0.226992 
sigma 3.758362 2.463029 1.259267 
beta1 36.675001 14.031268 17.338412 
beta2 -0.507215 0.213874 0.699148 
loglik at estimate: -357.4478 

Назначение наблюдения в данной категории сохраняется в mixmod объекта в матрице вероятностей в mixmod$posterior. Если мы извлечем назначенный класс и сравним его с истинным классом, то соответствие будет довольно хорошим (помните, что имена, присвоенные классам по модели смеси, являются произвольными, здесь comp 1 - ID 3 obvs и т. Д.):

predclass <- unlist(apply(mixmod$posterior, 1, function(x){names(which.max(x))})) 
table(dfall$ID, predclass) 
    predclass 
    comp.1 comp.2 comp.3 
    1  2  48  0 
    2  0  0  25 
    3  31  4  0 

Существует хороший обзор и описание моделей смеси и их реализация в R here.