[R-br] EMV

Walmes Zeviani walmes em ufpr.br
Segunda Junho 23 12:24:49 BRT 2014


Ana Paula,

A optim() otimiza qualquer função. Nem precisa ser uma função de
verossimilhança. Ela pode ser usada para definir delineamentos ótimos,
máximo de funções como curvas de lactação, enfim, serve à propósitos
gerais. Sendo assim, depende do usuário escrever a função objetivo e
passar. Abaixo seguem exemplos da gama com a optim() e da weibull com
optim(), nlm(), survival::survreg() e bbmle::mle2().
​
##-----------------------------------------------------------------------------
## Gama.

n <- 100;
shape <- 5; scale <- 3
curve(dgamma(x, shape, scale), 0, 10)
y <- rgamma(n, shape=shape, scale=scale)

ll <- function(theta, y){
    sum(dgamma(y, shape=theta[1], scale=theta[2], log=TRUE))
}

start <- c(5,3)
op <- optim(start, fn=ll, y=y, hessian=TRUE,
            method="BFGS", control=list(fnscale=-1))
op$par
op$value
op$hessian

##-----------------------------------------------------------------------------
## Weibull.

n <- 100;
shape <- 5; scale <- 3
curve(dweibull(x, shape, scale), 0, 10)
y <- rweibull(n, shape=shape, scale=scale)

ll <- function(theta, y){
    sum(dweibull(y, shape=theta[1], scale=theta[2], log=TRUE))
}

start <- c(5,3)
op <- optim(start, fn=ll, y=y, hessian=TRUE,
            method="BFGS", control=list(fnscale=-1))
op$par
op$value
op$hessian

##-----------------------------------------------------------------------------
## Usando a nlm (só minimiza).

ll <- function(theta, y){
    -sum(dweibull(y, shape=theta[1], scale=theta[2], log=TRUE))
}

start <- c(5,3)
op <- nlm(p=start, f=ll, y=y, hessian=TRUE)
op

##-----------------------------------------------------------------------------
## Usando a survreg (para a gama seria a glm(..., family=gamma)).

require(survival)

m0 <- survreg(Surv(time=y)~1, dist="weibull")
summary(m0)

exp(coef(m0)) ## = scale.
1/m0$scale    ## = shape.

##-----------------------------------------------------------------------------
## Usando a bblme.

require(bbmle)

## Forma de especificar I.
ml0 <- mle2(y~dweibull(shape=sh, scale=sc),
             start=list(sh=5, sc=3), data=list(y=y))
class(ml0)

## Vantagem é a disponibilidade de métodos.
summary(ml0)
vcov(ml0)
plot(profile(ml0))
confint(ml0)
confint(ml0, method="quad")

## Forma de específicar II.
ll <- function(shape, scale, y){
    -sum(dweibull(y, shape=shape, scale=scale, log=TRUE))
}

start <- list("shape"=5, "scale"=3)
ml1 <- mle2(minuslogl=ll, start=start, method="BFGS", data=list(y=y))

summary(ml1)
vcov(ml1)
plot(profile(ml1))
confint(ml1)
confint(ml1, method="quad")

##-----------------------------------------------------------------------------


À
​ disposição
Walmes.​
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20140623/cb2e13ec/attachment.html>


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