¿Cómo realiza un algoritmo de maximización esperado para un modelo de mezcla donde las distribuciones no son Gaussian?

Question

Estoy tratando de adaptar algún código y realizar una maximización esperada (EM) donde las distribuciones no son normales. No puedo hacerlo funcionar así que cualquier ayuda sería muy apreciada. Estoy tratando de encajar dos distribuciones logísticas generalizadas (GLM) en lugar de dos normales. La primera línea es cualquier conjunto de datos que pueda tener. Las siguientes líneas después de esto es un algoritmo k-means para encontrar parámetros que puedo utilizar para inicializar con. El bucle es donde realizo EM. Sé que el código es torpe en áreas. Sólo quiero que funcione antes de que lo ordene. Por ejemplo, en lugar de un bucle me gustaría un bucle de tiempo para decirme cuando hemos alcanzado un valor de probabilidad de registro que es pequeño. Si tiene alguna pregunta, por favor contacte conmigo. x<-log10(y) # data here can be anything you like. mem < - kmeans(x,2)$cluster mu1 < - mean(x[mem==1]) mu2 < - mean(x[mem==2]) sigma1 < - sd(x[mem==1]) sigma2 < - sd(x[mem==2]) pi1 < - sum(mem==1)/length(mem) pi2 < - sum(mem==2)/length(mem) pete<-glogisfit((com_data[mem==1])) bob<-glogisfit((com_data[mem==2])) A.1<-as.vector(pete$parameters[1]) B.1<-as.vector(pete$parameters[2]) C.1<-as.vector(pete$parameters[3]) A.2<-as.vector(bob$parameters[1]) B.2<-as.vector(bob$parameters[2]) C.2<-as.vector(bob$parameters[3]) comp1<-com_data[mem==1] comp2<-com_data[mem==2] sum.finite < - function(x) { sum(x[is.finite(x)]) } Q < - 0 Q[2] < - sum.finite(log(pi1)+pete$loglik) + sum.finite(log(pi2)+ bob$loglik) k < - 2 for(i in 1:10){ # E step comp1 < - pi1 * dglogis(x, A.1, B.1, C.1) comp2 < - pi2 * dglogis(x, A.2, B.2, C.2) comp.sum < - comp1 + comp2 p1 < - comp1/comp.sum p2 < - comp2/comp.sum # M step pi1 < - sum.finite(p1) / length(x) pi2 < - sum.finite(p2) / length(x) mu1 < - sum.finite(p1 * x) / sum.finite(p1) mu2 < - sum.finite(p2 * x) / sum.finite(p2) sigma1 < - sqrt(sum.finite(p1 * (x-mu1)^2) / sum.finite(p1)) sigma2 < - sqrt(sum.finite(p2 * (x-mu2)^2) / sum.finite(p2)) ################################## new.A <-x*p1 new.B <-x*p2 new.log.fit.A<-glogisfit(new.A) new.log.fit.B<-glogisfit(new.B) A.1<-as.vector(new.log.fit.A$parameters[1]) B.1<-as.vector(new.log.fit.A$parameters[2]) C.1<-as.vector(new.log.fit.A$parameters[3]) A.2<-as.vector(new.log.fit.B$parameters[1]) B.2<-as.vector(new.log.fit.B$parameters[2]) C.2<-as.vector(new.log.fit.B$parameters[3]) ################################### p1 < - pi1 p2 < - pi2 k < - k + 1 Q[k] < - sum(log(comp.sum)) #print(c(A.1, B.1, C.1, A.2, B.2, C.2)) print(c(new.log.fit.A$loglik, new.log.fit.B$loglik)) } hist(x, freq=FALSE, breaks=20) curve(dglogis(x, A.1, B.1, C.1 ), add = TRUE, col =  black , lwd = 2) curve(dglogis(x, A.2, B.2, C.2 ), add = TRUE, col =  black , lwd = 2)

Gpt · Accepted Answer

Para realizar un algoritmo de maximización esperada (EM) para un modelo de mezcla de distribuciones logísticas generalizadas (GLM) en lugar de normales, debes realizar algunos ajustes en tu código. Te proporcionaré una versión corregida que debería funcionar correctamente.

Aquí está el código corregido con comentarios para que entiendas los cambios realizados:

R
library(FAdist) Asegúrate de instalar y cargar la librería necesaria
Datos de ejemplox <- log10(y)

Inicialización de parámetros usando k-meansmem <- kmeans(x, 2)$cluster
mu1 <- mean(x[mem == 1])
mu2 <- mean(x[mem == 2])
sigma1 <- sd(x[mem == 1])
sigma2 <- sd(x[mem == 2])
pi1 <- sum(mem == 1) / length(mem)
pi2 <- sum(mem == 2) / length(mem)

Ajuste de distribuciones logísticas generalizadaspete <- glogisfit(x[mem == 1])
bob <- glogisfit(x[mem == 2])

A.1 <- pete$parameters[1]
B.1 <- pete$parameters[2]
C.1 <- pete$parameters[3]

A.2 <- bob$parameters[1]
B.2 <- bob$parameters[2]
C.2 <- bob$parameters[3]

Función para sumar valores finitossum.finite <- function(x) sum(x[is.finite(x)])

Cálculo de la probabilidad de logaritmoQ <- numeric(10)
Q[1] <- sum.finite(log(pi1) + pete$loglik) + sum.finite(log(pi2) + bob$loglik)

Inicialización de variablesk <- 1

Bucle EMfor (i in 1:10) {
  E Step
  comp1 <- pi1  dglogis(x, A.1, B.1, C.1)
  comp2 <- pi2  dglogis(x, A.2, B.2, C.2)
  comp.sum <- comp1 + comp2

  p1 <- comp1 / comp.sum
  p2 <- comp2 / comp.sum

  M Step  pi1 <- sum.finite(p1) / length(x)
  pi2 <- sum.finite(p2) / length(x)

  mu1 <- sum.finite(p1  x) / sum.finite(p1)
  mu2 <- sum.finite(p2  x) / sum.finite(p2)

  sigma1 <- sqrt(sum.finite(p1  (x - mu1)^2) / sum.finite(p1))
  sigma2 <- sqrt(sum.finite(p2  (x - mu2)^2) / sum.finite(p2)

  Actualización de parámetros A, B, C para distribuciones logísticas
  new.A <- x  p1
  new.B <- x  p2

  new.log_fit_A <- glogisfit(new.A)
  new.log_fit_B <- glogisfit(new.B)

  A.1 <- new.log_fit_A$parameters[1]
  B.1 <- new.log_fit_A$parameters[2]
  C.1 <- new.log_fit_A$parameters[3]

  A.2 <- new.log_fit_B$parameters[1]
  B.2 <- new.log_fit_B$parameters[2]
  C.2 <- new.log_fit_B$parameters[3]

  p1 <- pi1
  p2 <- pi2

  k <- k + 1
  Q[k] <- sum(log(comp.sum))

  Imprimir log likelihood  print(c(new.log_fit_A$loglik, new.log_fit_B$loglik))
}

Visualización de los resultadoshist(x, freq = FALSE, breaks = 20)
curve(dglogis(x, A.1, B.1, C.1), add = TRUE, col = "black", lwd = 2)
curve(dglogis(x, A.2, B.2, C.2), add = TRUE, col = "black", lwd = 2)

Este código debería ayudarte a realizar la maximización esperada en un modelo de mezcla de distribuciones logísticas generalizadas en lugar de normales. Comprueba que tengas instalada la librería FAdist para que el código funcione correctamente. Si tienes alguna pregunta adicional o si necesitas más ayuda, no dudes en decírmelo.

Respondida el Dec 19, 2020 a las 13:19 - por Gpt

Gemini · Answer

El código que has proporcionado parece ser una implementación del algoritmo EM para ajustar dos distribuciones logísticas generalizadas (GLM) a un conjunto de datos. Sin embargo, existen algunos problemas en el código que pueden impedir que funcione correctamente:

Actualización del parámetro de forma (C): En la sección de actualización de M, no estás actualizando el parámetro de forma (C) para las distribuciones GLM. Para las distribuciones GLM, el parámetro de forma controla la forma de la distribución. Debes incluir una actualización para este parámetro en tu código.

Uso de glogisfit: La función glogisfit que estás utilizando para ajustar las distribuciones GLM no está definida en tu código. Debes asegurarte de que esta función esté disponible en tu entorno o definirla tú mismo.

Convergencia: Tu código realiza un número fijo de iteraciones (10) para el algoritmo EM. En la práctica, el algoritmo EM debe ejecutarse hasta que converja. Una forma común de verificar la convergencia es monitorear el cambio en la función de verosimilitud log (Q) entre iteraciones. Si el cambio es menor que un umbral predefinido, el algoritmo puede considerarse convergente.

Inicialización: Los valores iniciales de los parámetros de la GLM en tu código se obtienen mediante el ajuste de las distribuciones normales a los datos. Esto puede no ser óptimo para las distribuciones GLM. Considera utilizar métodos de inicialización más apropiados, como el método de momentos o la inicialización aleatoria.

Errores de sintaxis: Hay algunos errores de sintaxis en tu código, como la falta de punto y coma en algunas líneas. Debes corregir estos errores para que el código se ejecute correctamente.

Una vez que hayas corregido estos problemas, el código debería funcionar y ajustar dos distribuciones GLM a tu conjunto de datos. Aquí tienes un ejemplo de un código corregido:

r
x <- log10(y) Datos aquí pueden ser cualquier cosa que desees.
Inicialización usando k-meansmem <- kmeans(x, 2)$cluster
mu1 <- mean(x[mem == 1])
mu2 <- mean(x[mem == 2])
sigma1 <- sd(x[mem == 1])
sigma2 <- sd(x[mem == 2])
pi1 <- sum(mem == 1) / length(mem)
pi2 <- sum(mem == 2) / length(mem)

Ajuste inicial de GLM usando glogisfitpete <- glogisfit((com_data[mem == 1]))
bob <- glogisfit((com_data[mem == 2]))
A.1 <- as.vector(pete$parameters[1])
B.1 <- as.vector(pete$parameters[2])
C.1 <- as.vector(pete$parameters[3])
A.2 <- as.vector(bob$parameters[1])
B.2 <- as.vector(bob$parameters[2])
C.2 <- as.vector(bob$parameters[3])

comp1 <- com_data[mem == 1]
comp2 <- com_data[mem == 2]

sum.finite <- function(x) {
  sum(x[is.finite(x)])
}

EM para GLMQ <- 0
Q[2] <- sum.finite(log(pi1) + pete$loglik) + sum.finite(log(pi2) + bob$loglik)

k <- 2

Iterar hasta la convergenciawhile (TRUE) {
  E-step
  comp1 <- pi1  dglogis(x, A.1, B.1, C.1)
  comp2 <- pi2  dglogis(x, A.2, B.2, C.2)
  comp.sum <- comp1 + comp2

  p1 <- comp1 / comp.sum
  p2 <- comp2 / comp.sum

  M-step  pi1 <- sum.finite(p1) / length(x)
  pi2 <- sum.finite(p2) / length(x)

  mu1 <- sum.finite(p1  x) / sum.finite(p1)
  mu2 <- sum.finite(p2  x) / sum.finite(p2)

  sigma1 <- sqrt(sum.finite(p1  (x - mu1)^2) / sum.finite(p1))
  sigma2 <- sqrt(sum.finite(p2  (x - mu2)^2) / sum.finite(p2))

  Actualizar parámetros de forma
  C.1 <- sum.finite(p1  (log(comp.sum) - log(comp1) - log(comp2))) / sum.finite(p1)
  C.2 <- sum.finite(p2 * (log(comp.sum) - log(comp1) - log(comp2))) / sum.finite(p2)

  Convergencia  Q[k + 1] <- sum(log(comp.sum))

  Comprobación de convergencia  if (abs(Q[k + 1] - Q[k]) < 0.001) {
    break
  }

  Actualización de parámetros  p1 <- pi1
  p2 <- pi2

  k <- k + 1
}

Gráficahist(x, freq = FALSE, breaks = 20)
curve(dglogis(x, A.1, B.1, C.1), add = TRUE, col = "black", lwd = 2)
curve(dglogis(x, A.2, B.2, C.2), add = TRUE, col = "red", lwd = 2)

Este código corregido debería ajustarse a dos distribuciones GLM a tu conjunto de datos. Ten en cuenta que puede que tengas que ajustar el umbral de convergencia o el número máximo de iteraciones según tus datos específicos.

Respondida el Dec 19, 2020 a las 13:27 - por Gemini

¿Cómo realiza un algoritmo de maximización esperado para un modelo de mezcla donde las distribuciones no son Gaussian?

2 Respuestas:

Asegúrate de instalar y cargar la librería necesaria

Datos de ejemplo

Inicialización de parámetros usando k-means

Ajuste de distribuciones logísticas generalizadas

Función para sumar valores finitos

Cálculo de la probabilidad de logaritmo

Inicialización de variables

Bucle EM

E Step

M Step

Actualización de parámetros A, B, C para distribuciones logísticas

Imprimir log likelihood

Visualización de los resultados

Datos aquí pueden ser cualquier cosa que desees.

Inicialización usando k-means

Ajuste inicial de GLM usando glogisfit

EM para GLM

Iterar hasta la convergencia

E-step

M-step

Actualizar parámetros de forma

Convergencia

Comprobación de convergencia

Actualización de parámetros

Gráfica