Pessoal, boa noite!

Tenho um conjunto de dados que provavelmente segue uma distribuição Weibull. 

Os dados são os seguintes:

y<-c(295000,

869000 ,

869900 ,

1573335 ,

151400 ,

152000 ,

183700 ,

218000 ,

30200 ,

45100 ,

46900 ,

47300

)

Efetuei a estimação dos parâmetros de duas formas, uma utilizando as expressões analíticas onde tive que usar o método de Newton-Raphson... E outra utilizando a função optim em que selecionei o método de Nelder Mead... Os dois algoritmos seguem abaixo:

###### Newton Paphson Method ##########

n=length(y)

x=vector()

##### Função que quero obter a raiz (estimador de lambda)

f <- quote((sum((y^delta)*log(y))/sum(y^delta))-(1/delta)-(1/n)*sum(log(y)))

fx0 <- function(delta){ eval(f) }

##### Primeira derivada da função (estimador de lambda)

f1 <-quote((1/delta^2)-((sum((log(y))*(y^delta)))^2/(sum(y^delta))^2)+(sum((log(y))^2*(y^delta))/sum(y^delta)))

fx1 <- function(delta){ eval(f1) }

### Obtendo a estimativa para delta

i <- 1       # valor inicial para o passo

dif <- 1   # valor inical para a diferença entre duas avaliações

tol <- 10^-6 # diferência mínima entre duas avaliações (tolerância)

passo <- 1000   # número máximo de passos

x <- 0.3      # valor inicial para a raiz            

while(i<passo & dif>tol){

  x[i+1] <- x[i]-(fx0(x[i])/fx1(x[i]))

  dif <- abs(x[i+1]-x[i])

d<-x[i] 

i<-i+1

}

### Obtendo a estimativa para alpha

alpha=((sum(y^d))/n)^(1/d)

### d é o parâmetro de forma e alpha é o parâmetro de escala

d

alpha

########## Teste de aderência 

ypop<-rweibull(1000,d,alpha)

ks.test(y,ypop)

#####################################################

###### Nelder-Mead através da função Optim ##########

n=length(y)

vero <- function(par,x){

k = par[1] #### Shape parameter

lambda = par[2] #### Scale Parameter

   saida<-(n*log(k/lambda) + (k-1)*sum(log(x/lambda)) - sum((x/lambda)^k))

return(-saida)

}

saiday<-optim(par=c(1.8,4),fn=vero,

      method= "Nelder-Mead",x=y

      )

k<-saiday[1]$par[1]

lambda<-saiday[1]$par[2]

### k é o parâmetro de forma e lambda e o parâmetro de escala

k

lambda

########## Teste de aderência 

ypop<-rweibull(1000,k,lambda)

ks.test(y,ypop)

###############################################

O que acontece é que as estimativas através dos dois métodos estão dando muito diferentes e através de um KS.test concluo que as estimativas através de Newton-Raphson são muito boas, pois o p-valor fica realmente bem alto; já as estimativas utilizando Nelder-Mead ficam péssimas e o ajuste, claro, fica péssimo com p-valor muito pequeno....

Gostaria de saber se existe alguma explicação para que isso esteja ocorrendo... O método de Nelder-Mead não é totalmente indicado?? O que estaria provocando essa disparidade nas estimativas... Sei que são dois métodos de otimização diferentes, mas acredito que os dois eram pra fornecer vaores bem próximos, concordam???

Observe que se eu gerar uma weibull, por exemplo:

Y<-rweibull(12,5,6)

As estimativas serão bem próximas e o ks.test me dará um p-valor bem condizente com a realidade (isso para os dois métodos de otimização)...

O que pode estar acontecendo?? Vcs tem alguma idéia??? O problema está nos dados que estou analizando?? que tipo de problema seria esse??

Estou confuso pois quando gero dados no R os dois métodos apresentam boas estimtivas e p-valor do ks.test é bem condizente... Entretanto se uso os dados que tenho, os parâmetros só são bem estimados pelo método de Newton-Raphson o outro método fica péssimo... 

Agradeço qualquer dica!

Abraço,

Romero.

_______________________________________________
R-br mailing list
R-br@listas.c3sl.ufpr.br
https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça código mínimo reproduzível.