[R-br] Estimativas de parâmetros de uma Weibull utilizando duas metodologias: Newton-Raphson "vs" Nelder Mead (estimativas muito diferentes) o que pode estar acontecendo??

Romero Luiz M. Sales Filho romero.sfilho em gmail.com
Quarta Maio 6 19:20:27 BRT 2015


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.
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20150506/f9b8921a/attachment.html>


Mais detalhes sobre a lista de discussão R-br