for(i in 1:k){ x[i] <- rgamma(n,alp,lamb) alphat[i] <- lambhat[i] <- (n*alp)/x[i] } 极大似然估计参数 alp 的R代码

假设样本 $X_1, X_2, ..., X_n$ 来自参数为 $\alpha$ 和 $\lambda$ 的 Gamma 分布，那么其似然函数为：

$L(\alpha, \lambda) = \prod_{i=1}^{n} f(X_i|\alpha, \lambda) = \prod_{i=1}^{n} \frac{\lambda^{\alpha}}{\Gamma(\alpha)} X_i^{\alpha - 1} e^{-\lambda X_i}$

其中 $\Gamma(\cdot)$ 表示 Gamma 函数。对数似然函数为：

$l(\alpha, \lambda) = \log L(\alpha, \lambda) = n \log \lambda \alpha - n \log \Gamma(\alpha) + (\alpha - 1) \sum_{i=1}^{n} \log X_i - \lambda \sum_{i=1}^{n} X_i$

极大似然估计就是要找到最大化 $l(\alpha, \lambda)$ 的参数 $\alpha$ 和 $\lambda$。我们可以利用数值优化算法来求解。在 R 中，可以使用 optimize 函数来实现。具体的代码如下：

# 定义对数似然函数
log_likelihood <- function(par, x) {
  alpha <- par[1]
  lambda <- par[2]
  n <- length(x)
  log_l <- n*log(lambda*alpha) - n*log(gamma(alpha)) + (alpha-1)*sum(log(x)) - lambda*sum(x)
  return(-log_l) # optimize 函数的要求是最小化函数，因此需要取相反数
}

# 生成样本数据
set.seed(123)
n <- 200
alpha <- 2
lambda <- 1
x <- rgamma(n, alpha, lambda)

# 利用 optimize 函数求解极大似然估计
res <- optimize(log_likelihood, c(0.1, 10), x=x)
alpha_hat <- res$minimum
lambda_hat <- res$objective / alpha_hat / n

其中 optimize 函数的第一个参数是要最小化的函数，第二个参数是参数的初始值，x 是样本数据。最终的估计结果存储在 alpha_hat 和 lambda_hat 中。