Ошибка: неожиданное '}' в "}" в моем EM для модели смеси - PullRequest
0 голосов
/ 24 июня 2019

Я написал ЭМ для модели смеси, которая состоит из двух олов. Но это показывает:

"Ошибка: неожиданная '}' в"} "".

Я проверил свой код и не обнаружил проблем с "}".

Буду признателен, если будут обнаружены другие ошибки в моем коде.

Бумага, о которой я говорил, здесь .

И мой код выглядит следующим образом:

   #prepare data
slope1 <- -.3;slope2 <- .3;slope3 <- 1.8; slope4 <- 0.5;intercept1 <- 1.5    
age <- sample(seq(-2,2,len=201), 40) 
grade <- sample(seq(-2,2,len=201), 40) 
not_smsa <- sample(seq(-2,2,len=201), 40) 
wage <- intercept1 + slope1*age +slope2*grade + slope3*not_smsa + rnorm(length(age),0,.15) 

X <- cbind(1, age , grade , not_smsa )
Y <- wage
N<- length(Y) 
mydata <- cbind.data.frame(X,Y)


library(mclust)
em <- Mclust(Y)
summary(em, parameters = T)

ans <- lm(wage ~ age + grade + not_smsa , 
          data = mydata) 
init_params <- c(coef(ans) , coef(ans))
params <-init_params
pi1 <- 0.373321
pi2 <- 0.626679 
loglik<- rep(NA, 1000)

mu1 <- X %*% params[1:4]
sigma1 <- sqrt((t((Y-X %*% params[1:4]))%*%(Y-X %*% params[1:4]))/N)

mu2 <- X %*% params[5:8]
sigma2 <- sqrt((t((Y-X %*% params[5:8]))%*%(Y-X %*% params[5:8]))/N)

loglik[1] <- 0
loglik[2] <- sum(log(pi1*dnorm(Y,mu1,sigma1)+pi2*dnorm(Y,mu2,sigma2)))

k<-2

wls <- function(X,Y,W) {
  solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% Y
}

while (abs(loglik[k]-loglik[k-1]) >= 0.00001)(na.rm = T){

  # E step
  w1<-pi1*dnorm(Y,mu1,sigma1)/(pi1*dnorm(Y,mu1,sd=sigma1)+pi2*dnorm(Y,mu2,sigma2))
  w2<-pi2*dnorm(Y,mu2,sigma2)/(pi1*dnorm(Y,mu1,sd=sigma1)+pi2*dnorm(Y,mu2,sigma2))
  w1[is.na(w1)] <- 0.5
  w2[is.na(w2)] <- 0.5

  # M step
  pi1<-sum(w1)/length(Y)
  pi2<-sum(w2)/length(Y)

  w <- cbind(w1,w2)
  params[1:4] <- wls(X,Y,diag(w[,1]))
  params[5:8] <- wls(X,Y,diag(w[,2]))

  r1<-Y-X %*% params[1:4]
  r2<-Y-X %*% params[5:8]

  sigma1h <- sum(w1*r1)/sum(w1)
  sigma1h1 <- exp(sigma1h)
  sigma2h <- sum(w2*r2)/sum(w2)
  sigma2h1 <- exp(sigma2h)  
  #  loglik[k]
  loglik[k+1] <- sum(log(pi1*dnorm(Y,mu1,sigma1h1)+pi2*dnorm(Y,mu2,sigma2h1)),na.rm=T)
  k<- k+1
}

Спасибо за ваше терпение, и я с нетерпением жду ваших отзывов.

...