1

我编写了一个函数来计算广义帕累托分布的 MLE 估计值。当我将它与任何数据一起使用时,虽然我遇到了这样的错误

1: In log(beta * ksi) : NaNs produced
2: In nlm(loglik, theta, stepmax = 5000, iterlim = 1000) :
  NA/Inf replaced by maximum positive value

我想知道是否有人可以发现我的代码有任何错误?

MLGPD<-function(data){
    xi0 <- 1
    beta0 <- 360
    theta <- c(xi0, beta0)
    excess <- data
    assign("tmp", excess)
    loglik <- function(theta){
        ksi <- theta[1]
        beta <- theta[2]
        y <- ((tmp - 0.1)/beta)
        f <- ((1/ksi)+1)*sum(log(1+y)) + length(tmp) * log(beta*ksi)
        f
    }
    fit <- nlm(loglik, theta, stepmax = 5000, iterlim= 1000)
    return(fit)
    par.ests <- fit$x
    return(par.ests)
}

#Checking our MLE algorithm works:

rgpd<-function(n,ksi, beta){
    10000+beta*(((1-runif(n, min=0, max=1))^-ksi)-1)
}

rgpd1 <- rgpd(100, 1, 2.5)

MLGPD(rgpd1)

谢谢!

4

0 回答 0