大神们，我在用R做极大似然估计时用optim函数，遇到下面情况。

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

该怎么解决呢？用加上什么参数限制吗？

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏1 回帖

关键词：极大似然估计 Optim 似然估计 极大似然 Tim

时间是一个疯子 发表于 2016-3-8 18:39
该怎么解决呢？用加上什么参数限制吗？

是否有无效数据？nan表示求根时出现了无效的数或者是负数。

lowbrow 发表于 2016-3-8 18:49
是否有无效数据？nan表示求根时出现了无效的数或者是负数。

我把代码发给你看下。
lmdata <- function(n, theta){
  
  sigma <- sqrt(theta[length(theta)])
  bet <- theta[-length(theta)]
  
  x <- runif(n, -2, 2)  #  covariates
  e <- rnorm(n, 0, sigma)  # random error
  
  ## including intercept
  y <- as.vector(cbind(1, x) %*% bet) + e
  
  dat <- cbind(y, x)
  return(dat)
}   这个函数是为了产生随机数

下面是产生似然函数。
likfun <- function(dat, theta){
  y <- dat[, 1]
  x <- dat[, -1]
  sigma2 <- theta[length(theta)]
  bet <- theta[-length(theta)]
  xb <- as.vector(cbind(1, x) %*% bet)
  lik <- sum(log(1/sqrt(2 * pi * sigma2) *exp(-(y-xb)^2/(2*sigma2))))
  return(lik)
}


mleFun <- function(dat, theta){
  
  theta.ini <- c(0, 0, 1)
  fit <- optim(
    par = theta.ini,
    fn = function(theta){-likfun(dat, theta)},
    method = "BFGS"
  )
  
  theta.est <- fit$par
  
}
n <- 200
theta0 <- c(1, 2, 1^2)  假定三个待估参数分别为斜率，截距以及误差项的方差。
dat=lmdata(n,theta0)
mleFun(dat,theta0)

lowbrow 发表于 2016-3-8 18:49
是否有无效数据？nan表示求根时出现了无效的数或者是负数。

然后就出现了那个结果，搞不懂，，，好像也没什么错误。

可以在Rstudio中逐句执行。

lowbrow 发表于 2016-3-8 19:15
可以在Rstudio中逐句执行。

在Rstudio中也是如此，，，好无语。。。

log(1/sqrt(2 * pi * sigma2) *exp(-(y-xb)^2/(2*sigma2)
这句产生的向量可能有问题，产生的向量组中有负值。

lowbrow 发表于 2016-3-8 19:29
log(1/sqrt(2 * pi * sigma2) *exp(-(y-xb)^2/(2*sigma2)
这句产生的向量可能有问题，产生的向量组中有负值 ...

我要估计的是误差项的方差，设定的也是正数。。可能R识别不出来，这个要怎么解决呢？

把mleFun改成

mleFun <- function(dat, theta){
  
  theta.ini <- c(0, 0, 1)
  fit <- constrOptim(
    theta=theta.ini,
    f = function(theta){-likfun(dat, theta)}, ui=c(0,0,1), ci=c(0,0,0),
    method = "Nelder-Mead"
  )
  
  theta.est <- fit$par
  return(theta.est)
}

> n <- 200
> theta0 <- c(1, 2, 1^2)
> dat=lmdata(n,theta0)
> mleFun(dat,theta0)
[1] 0.9750839 2.0498282 0.8936549

{:3_41:}