4

Olá, o meu projeto precisa do algoritmo EM abaixo, onde está todos os códigos. O erro está no loop while, que é onde estão as etapas esperança e maximização. A mensagem de erro é:

"Error in while (abs(Elogv[r] - Elogv[r - 1]) >= 1e-06) { : missing value where TRUE/FALSE needed"**.

Como faço pra resolver esse erro, se o loop while não contém comandos TRUE e FALSE e, se já verifiquei detalhadamente que não há nenhum erro nos comandos e não há valor NAs?

Código:

n=100
u<-runif(n)
QUANTIL <- function(u){
  Q <- rep(NA, length(u))
  for (i in 1:length(u)) {
    if(u[i] <  0.2634253829){
      Q[i] <- 1*tan(pi*(0.9490353482*u[i]-0.5))+0
    }
    if(u[i]>=0.2634253829 && u[i] < 0.7365746171){
      Q[i] <-  1*qnorm(1.4428629504*u[i]-0.2214315)+0
    }
    if(u[i]>0.7365746171){
      Q[i] <- 1*tan(pi*(0.9490353482*u[i]-0.4490353))+0
    } 
  }
  return(Q)
}
x<-QUANTIL(u)
y<-c(sort(x))
i<-seq(1,n)
v<-c(i/(n+1))

t<-QUANTIL(v)
mi<-median(y)
s<-c(y[26:73])
sigma<-sqrt(sum((s-mi)^2)/(n-1))
p=0.4731492342

alpha<-(2*t^3)/(1+t^2)^2
beta<-(1-t^2)/(1+t^2)^2
eta<-(t^4-t^2)/(1+t^2)^2
lambda<-2*t/(1+t^2)^2
gama<-(-t^2)
delta<-2*t

k<-((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*((y-mi)^2)))/(((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*(y-mi)^2))+((((1-p)*1.0537015317/sigma*pi))*(1/(1+((y-mi)/sigma)^2))))
r<-2
Elogv<-sum(k*((-1/2)*((y-mi)/sigma)^2))-sum(k*log(sigma*sqrt(2*pi)))-sum((1-k)*log(sigma*pi))-sum((1-k)*log(1+((y-mi)/sigma)^2))+sum(k*log(p))+(n-sum(k))*log(1-p)+log(0.6930665173)*sum(k)+log(1.0537015317)*sum(1-k)
Elogv[1]<-0

while (abs(Elogv[r]-Elogv[r-1])>=0.000001) {

  w<-(2*beta-2*k*beta+k)
  q<-k*delta+2*lambda*(1-k)
  sigma<-(sum(y*w)*sum(q)-sum(w)*sum(y*q))/(-2*sum(alpha*(1-k))*sum(q)+sum(w)*sum(k*gama-1)+2*sum(w)*sum(eta*(1-k)))                                                  
  mi<-(sum(y*w)+2*sigma*sum(alpha*(1-k)))/sum(w)
  k<-((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*((y-mi)^2)))/(((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*(y-mi)^2))+((((1-p)*1.0537015317/sigma*pi))*(1/(1+((y-mi)/sigma)^2))))
  Elogv[r]<-sum(k*((-1/2)*((y-mi)/sigma)^2))-sum(k*log(sigma*sqrt(2*pi)))-sum((1-k)*log(sigma*pi))-sum((1-k)*log(1+((y-mi)/sigma)^2))+sum(k*log(p))+(n-sum(k))*log(1-p)+log(0.6930665173)*sum(k)+log(1.0537015317)*sum(1-k)
  r<-r+1
  

}
4
  • 1
    Eu desconheço um programa que interprete algorítimos. O nome do programa é Algorítimo EM? Se for onde eu faço o download e onde encontro a documentação? Commented 5/09/2019 às 3:04
  • Esse algoritmo EM não tem documentação, é teoria estatística, deve ter um pacote pronto mas nao está usado aí, está feito manualmente. O ponto central do problema é o loop while que está na parte final do algoritmo, que causa a mensagem de erro. Eu coloquei o código todo pra se alguém quiser rodar o loop while, pois necessitaria dos códigos desde o início. Commented 5/09/2019 às 10:06
  • 2
    Eu não sou especialista em "R", não sei nada sobre ,mas sei ler documentação e caçar bugs. Na documentação diz o erro acontece poque sua comparação retorna o valor NA e uma alternativa é utilizar o método isTRUE que retorna o valor lógico da expressão e manipula o resultado caso retorne NA.Ficaria: while (isTRUE(abs(Elogv[r]-Elogv[r-1])>=0.000001, na=FALSE)) Commented 5/09/2019 às 14:50
  • 2
    Deixo como comentário e não como uma resposta pois eu não sei nada sobre r e seria leviano da minha parte afirmar uma resposta. Então deixo um comentário com a minha opinião e espero que ajude. Commented 5/09/2019 às 14:53

1 Resposta 1

-2

Como faço pra resolver esse erro, se o loop while não contém comandos TRUE e FALSE e, se já verifiquei detalhadamente que não há nenhum erro nos comandos e não há valor NAs?

O While é uma estrutura de repetição que necessita de que uma condição lógica seja atendida para ser executado. Para ser mais claro, ele basicamente só vai repetir enquanto a condição apresentada, no seu caso, (abs(Elogv[r]-Elogv[r-1])>=0.000001), for verdadeira.

Para mais informação sobre While: (https://cran.r-project.org/doc/manuals/r-devel/R-lang.html#while)

Para sermos mais precisos, podemos dizer que o seu problema está na lógica de programação.

A variável Elogv foi iniciada recebendo valores de um somatório que na sequência foi substituído pelo valor 0 indicado na posição [1], Elogv[1]<-0, transformando também essa variável num vetor neste momento.

De acordo com a sua condição While, (abs(Elogv[r]-Elogv[r-1])>=0.000001), enquanto o valor absoluto de Elogv na posição [r] menos Elogv na posição [r-1] for maior ou igual a 0.000001 será feito o que estiver dentro das chaves{}.

Portanto, ela nunca será atendida, pois, sendo r igual a 2, logo Elogv[2] é igual a Na, porque em nenhum momento foi atribuído algum valor na posição [2], e Elogv[2-1], isto é, Elogv[1] é igual a 0, por isso nunca serão maior que 0.

Mas cuidado, embora pareça que o resultado deveria ser falso ao invés de emitir erro, você precisa se atentar que como não é possível subtrair um valor por outro que está como Na (Not available/Não disponível) na linguagem R, o resultado é Na. Por isso que acusou o erro dizendo que não foram encontrados valores onde TRUE (verdadeiro) ou FALSE (falso) são necessários. A condição não resultou em verdadeiro ou falso, apenas em Na.

Para resolver o problema você precisará entender a lógica necessária para o desenvolvimento do seu algoritmo, porque se prosseguirmos com a análise do script é possível notar que não bastará acabar com o alerta de erro acertando a atribuição inicial da variável, observe que dentro do While ocorrerá mais um Na devido a variável r ser acrescentada sem que Elogv receba atribuição de valor, logo Elogv[3] será igual a Na.

Para mais informação sobre Na: (https://cran.r-project.org/doc/manuals/r-devel/R-lang.html#NA-handling)

Além dessas observações, outro ponto importante são os cálculos, você precisa analisar se eles não resultarão num dado contínuo que fará com que o seu While nunca tenha fim.

Segue abaixo uma alternativa para acabar com o erro, lembrando que você precisa entender a lógica desse programa para poder fazer os devidos ajustes corretamente. Não é porque não emitirá erro que o código estará cumprindo 100% com o seu objetivo.

#inicia o vetor
Elogv <- vector()
#atribui valor na posição 1 do vetor
Elogv[1]<-0
#atribui valor na posição 2 do vetor
Elogv[2]<-sum(k*((-1/2)*((y-mi)/sigma)^2))-sum(k*log(sigma*sqrt(2*pi)))-sum((1-k)*log(sigma*pi))-sum((1-k)*log(1+((y-mi)/sigma)^2))+sum(k*log(p))+(n-sum(k))*log(1-p)+log(0.6930665173)*sum(k)+log(1.0537015317)*sum(1-k)

while (abs(Elogv[r]-Elogv[r-1])>=0.000001) {
  w<-(2*beta-2*k*beta+k)
  q<-k*delta+2*lambda*(1-k)
  sigma<-(sum(y*w)*sum(q)-sum(w)*sum(y*q))/(-2*sum(alpha*(1-k))*sum(q)+sum(w)*sum(k*gama-1)+2*sum(w)*sum(eta*(1-k)))                                                  
  mi<-(sum(y*w)+2*sigma*sum(alpha*(1-k)))/sum(w)
  k<-((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*((y-mi)^2)))/(((p*0.6930665173/sigma*sqrt(2*pi))*exp((-1/2*sigma^2)*(y-mi)^2))+((((1-p)*1.0537015317/sigma*pi))*(1/(1+((y-mi)/sigma)^2))))
  #foi invertido a ordem do r e Elogv
  r<-r+1
  Elogv[r]<-sum(k*((-1/2)*((y-mi)/sigma)^2))-sum(k*log(sigma*sqrt(2*pi)))-sum((1-k)*log(sigma*pi))-sum((1-k)*log(1+((y-mi)/sigma)^2))+sum(k*log(p))+(n-sum(k))*log(1-p)+log(0.6930665173)*sum(k)+log(1.0537015317)*sum(1-k)
}
2
  • 1
    Sua resposta pode ser melhorada adicionando informações fundamentais. Basta editar para adicionar mais detalhes, como citações ou documentação para que outras pessoas possam confirmar que sua resposta está correta. Você pode encontrar mais informações sobre como escrever boas respostas na Central de Ajuda.
    – Comunidade Bot
    Commented 3/05/2022 às 10:44
  • 1
    Como está escrito atualmente, sua resposta não é clara. Por favor, pode editar para adicionar mais detalhes, que ajudarão outras pessoas a entender como isso resolve a pergunta feita. Você pode encontrar mais informações sobre como escrever boas respostas na Central de Ajuda.
    – Comunidade Bot
    Commented 3/05/2022 às 23:48

Você deve fazer log-in para responder a esta pergunta.

Esta não é a resposta que você está procurando? Pesquise outras perguntas com a tag .