seria bom mudar a escala dos dados para evitar problemas numericos
### estima para y/1000000
o <- optim(c(1.2, 1), function(p)
-sum(dweibull(y/1e6, p[1], p[2], log=TRUE)),
method='BFGS', hessian=TRUE)
### estima para 1000 reamostras (bootstrap)
oo <- t(replicate(1000, optim(o$par, function(p)
-sum(dweibull(sample(y, replace=TRUE)/1e6,
p[1], p[2], log=TRUE)), method='BFGS')$par))
(tb <- rbind(sd.assintotico=sqrt(diag(solve(o$hess))),
sd.bootstrap=apply(oo,2,sd)))
### visualiza estimativas bootstrap e adiciona SD assintotico
par(mfrow=c(2,2), mar=c(3,3,1,1), mgp=c(2,1,0))
hist(y/1e6, main='dados')
plot(function(x) dweibull(x, o$par[1], o$par[2]), 0, 2, add=T, col=2)
hist(oo[,1], main='shape')
abline(v=o$par[1]+c(-2,0,2)*tb[1,1], col=2, lty=c(2,1,2))
hist(oo[,2], main='scale')
abline(v=o$par[2]+c(-2,0,2)*tb[1,2], col=2, lty=c(2,1,2))
Elias
_______________________________________________
R-br mailing list
[email protected]
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.
