在 R 中编写适当的法线 log-likelihood

Question

我对以下型号有疑问，

我想对μ和tau进行推断，u是已知向量，x是数据向量。 log-likelihood 是

我在 R 中编写 log-likelihood 时遇到问题。

x <- c(3.3569,1.9247,3.6156,1.8446,2.2196,6.8194,2.0820,4.1293,0.3609,2.6197)
mu <- seq(0,10,length=1000)

normal.lik1<-function(theta,x){ 
  u <- c(1,3,0.5,0.2,2,1.7,0.4,1.2,1.1,0.7)  
  mu<-theta[1] 
  tau<-theta[2] 
  n<-length(x) 

logl <-  sapply(c(mu,tau),function(mu,tau){logl<- -0.5*n*log(2*pi) -0.5*n*log(tau^2+u^2)- (1/(2*tau^2+u^2))*sum((x-mu)^2) } )

  return(logl) 
  } 

#test if it works for mu=1, tau=2
head(normal.lik1(c(1,2),x))
#Does not work..

我希望能够插入 mu 的矢量并将其绘制在 mu 上以获得固定的 tau 值，比如 2。我还想使用 optim 函数找出 tau 和 mu 的 MLE。我试过了：

theta.hat<-optim(c(1,1),loglike2,control=list(fnscale=-1),x=x,,method="BFGS")$par

但它不起作用..对我如何写可能性有什么建议吗？

Answer 1

首先，正如您问题的评论中提到的，没有必要使用sapply()。您可以简单地使用 sum()——就像在 logLikelihood 的公式中一样。

我在 normal.lik1() 中更改了这部分并将分配给 logl 的表达式乘以负 1，以便函数计算负对数似然。您想在 theta 上搜索 最小值 ，因为函数 returns 为正值。

x < c(3.3569,1.9247,3.6156,1.8446,2.2196,6.8194,2.0820,4.1293,0.3609,2.6197)
u <- c(1,3,0.5,0.2,2,1.7,0.4,1.2,1.1,0.7) 

normal.lik1 <- function(theta,x,u){ 
  mu <- theta[1] 
  tau <- theta[2] 
  n <- length(x) 
  logl <- - n/2 * log(2*pi) - 1/2 * sum(log(tau^2+u^2)) - 1/2 * sum((x-mu)^2/(tau^2+u^2))
  return(-logl) 
}

这可以使用 nlm() 来完成，例如

nlm(normal.lik1, c(0,1), hessian=TRUE, x=x,u=u)$estimate

其中 c(0,1) 是算法的起始值。

要绘制 mu 和一些固定 tau 值范围的对数似然值，您可以调整函数，使 mu 和 tau 是单独的数字参数.

normal.lik2 <- function(mu,tau,x,u){ 
  n <- length(x) 
  logl <- - n/2 * log(2*pi) - 1/2 * sum(log(tau^2+u^2)) - 1/2 * sum((x-mu)^2/(tau^2+u^2))
  return(logl) 
}

然后为 mu 定义一些范围，计算对数似然并使用 plot()。

range.mu <- seq(-10,20,0.1)

loglik <- sapply(range.mu, function(m) normal.lik2(mu=m,tau=2,x=x,u=u))

plot(range.mu, loglik, type = "l")

我敢肯定有更优雅的方法可以做到这一点，但这确实有效。

在 R 中编写适当的法线 log-likelihood

Writing a proper normal log-likelihood in R

statistics

r

normal-distribution

mle

log-likelihood