为什么 optim() 总是给我相同的结果?
Why optim() always gives me the same result?
我正在尝试使用引导程序来计算参数 v 的 MLE 的置信区间。但是在每次迭代中,optim()
都会给我相同的结果。为什么会这样?
这是我的代码:
# Create log-likelihood function
dofloglik = function(v, x) { n <- length(x)
loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
- (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
return(-loglik) }
# Sample
y = c(0.97, 0.41, -0.13, -0.18, 1.69,
0.17, -0.46, -0.65, 2.66, -0.28,
2.30, -0.15, 0.24, 0.43, -0.54,
-1.55, 2.13, -0.48, 1.24, 1.05)
# Bootstrapping
n = length(y)
sims = 10
mle = numeric(sims)
for (i in 1:sims) {
bs = sample(y, n, replace=TRUE)
mle[i] = optim(1, dofloglik, x = bs, method = "CG")$par
print(bs)
}
mle = sort(mle)
mle
输出为:
[1] 1.05 1.69 -0.13 2.66 2.13 -0.28 -0.54 -0.65 1.69 0.43 -0.18 0.17 1.24 -0.13 -1.55 0.41 -0.18 0.24 -0.54 2.13
[1] -0.13 -0.18 2.13 -1.55 -0.65 -0.28 1.24 2.13 -0.48 2.30 -0.28 -0.54 -0.13 0.41 1.24 1.24 -0.54 -0.18 2.66 -0.54
[1] 2.30 0.43 -1.55 1.05 2.13 -0.18 1.24 0.24 -0.13 2.30 -0.48 2.66 -0.13 -1.55 -0.54 -0.13 0.43 -0.13 0.24 2.30
[1] 1.24 1.24 -0.54 2.30 -0.46 -0.65 2.13 2.66 -0.46 -0.13 -0.18 2.66 0.24 -0.48 -1.55 -1.55 1.05 0.24 0.17 -0.46
[1] 0.97 2.66 2.13 0.24 -1.55 0.43 2.13 0.43 -0.18 -1.55 -0.46 1.69 -0.48 2.66 -0.18 2.30 2.66 -1.55 -0.54 -0.54
[1] 0.41 0.43 1.24 1.69 -0.13 0.41 1.69 -0.54 -0.28 2.13 -0.46 2.30 0.17 0.97 1.24 2.30 2.66 0.43 0.43 1.24
[1] 1.24 1.24 1.05 2.66 0.17 2.13 0.17 -0.13 2.30 -0.46 -1.55 -0.28 1.24 2.30 -0.48 0.24 -0.54 0.41 -0.65 -1.55
[1] 2.66 -0.54 1.05 -0.15 0.17 -1.55 0.41 1.69 -1.55 -1.55 -0.28 -0.46 -0.48 -0.13 -0.46 0.43 1.24 0.24 -0.46 -0.28
[1] 1.05 0.24 2.13 0.97 1.69 1.05 2.13 -0.15 -0.48 -1.55 1.05 -0.15 0.43 -0.13 -0.28 0.17 2.66 -0.15 1.24 -0.28
[1] -0.48 -0.18 0.24 2.30 -0.46 -0.54 0.43 -0.54 2.66 -0.48 2.66 0.24 2.13 0.97 1.05 -0.18 2.30 -0.13 -0.46 1.24
> mle = sort(mle)
> mle
[1] 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466
> dofloglik = function(v, x) { n <- length(x)
+ loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
+ - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
+ return(-loglik) }
>
> dofloglik(1, 1:4)
[1] 2.28946
> dofloglik(1, 2:5)
[1] 2.28946
> dofloglik(1, 3:6)
[1] 2.28946
这可能不是您想要的。
我怀疑不是您想要的部分是您将某些内容分成两行。除了每一行都是一个有效的命令,所以 R 自己执行第一行。
+ loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
+ - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
所以 loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
是 运行 并且是有效的,所以它被存储到 loglik 中。第二行 - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
有效。它执行但不存储任何内容,因为我们不在顶层,它甚至不打印。
您要么需要确保第一行的末尾不会产生完整的语句(通过在第一行的末尾添加减号),要么采取其他措施以不同的方式将其拆分.
我正在尝试使用引导程序来计算参数 v 的 MLE 的置信区间。但是在每次迭代中,optim()
都会给我相同的结果。为什么会这样?
这是我的代码:
# Create log-likelihood function
dofloglik = function(v, x) { n <- length(x)
loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
- (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
return(-loglik) }
# Sample
y = c(0.97, 0.41, -0.13, -0.18, 1.69,
0.17, -0.46, -0.65, 2.66, -0.28,
2.30, -0.15, 0.24, 0.43, -0.54,
-1.55, 2.13, -0.48, 1.24, 1.05)
# Bootstrapping
n = length(y)
sims = 10
mle = numeric(sims)
for (i in 1:sims) {
bs = sample(y, n, replace=TRUE)
mle[i] = optim(1, dofloglik, x = bs, method = "CG")$par
print(bs)
}
mle = sort(mle)
mle
输出为:
[1] 1.05 1.69 -0.13 2.66 2.13 -0.28 -0.54 -0.65 1.69 0.43 -0.18 0.17 1.24 -0.13 -1.55 0.41 -0.18 0.24 -0.54 2.13
[1] -0.13 -0.18 2.13 -1.55 -0.65 -0.28 1.24 2.13 -0.48 2.30 -0.28 -0.54 -0.13 0.41 1.24 1.24 -0.54 -0.18 2.66 -0.54
[1] 2.30 0.43 -1.55 1.05 2.13 -0.18 1.24 0.24 -0.13 2.30 -0.48 2.66 -0.13 -1.55 -0.54 -0.13 0.43 -0.13 0.24 2.30
[1] 1.24 1.24 -0.54 2.30 -0.46 -0.65 2.13 2.66 -0.46 -0.13 -0.18 2.66 0.24 -0.48 -1.55 -1.55 1.05 0.24 0.17 -0.46
[1] 0.97 2.66 2.13 0.24 -1.55 0.43 2.13 0.43 -0.18 -1.55 -0.46 1.69 -0.48 2.66 -0.18 2.30 2.66 -1.55 -0.54 -0.54
[1] 0.41 0.43 1.24 1.69 -0.13 0.41 1.69 -0.54 -0.28 2.13 -0.46 2.30 0.17 0.97 1.24 2.30 2.66 0.43 0.43 1.24
[1] 1.24 1.24 1.05 2.66 0.17 2.13 0.17 -0.13 2.30 -0.46 -1.55 -0.28 1.24 2.30 -0.48 0.24 -0.54 0.41 -0.65 -1.55
[1] 2.66 -0.54 1.05 -0.15 0.17 -1.55 0.41 1.69 -1.55 -1.55 -0.28 -0.46 -0.48 -0.13 -0.46 0.43 1.24 0.24 -0.46 -0.28
[1] 1.05 0.24 2.13 0.97 1.69 1.05 2.13 -0.15 -0.48 -1.55 1.05 -0.15 0.43 -0.13 -0.28 0.17 2.66 -0.15 1.24 -0.28
[1] -0.48 -0.18 0.24 2.30 -0.46 -0.54 0.43 -0.54 2.66 -0.48 2.66 0.24 2.13 0.97 1.05 -0.18 2.30 -0.13 -0.46 1.24
> mle = sort(mle)
> mle
[1] 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466 47.16466
> dofloglik = function(v, x) { n <- length(x)
+ loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
+ - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
+ return(-loglik) }
>
> dofloglik(1, 1:4)
[1] 2.28946
> dofloglik(1, 2:5)
[1] 2.28946
> dofloglik(1, 3:6)
[1] 2.28946
这可能不是您想要的。
我怀疑不是您想要的部分是您将某些内容分成两行。除了每一行都是一个有效的命令,所以 R 自己执行第一行。
+ loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
+ - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
所以 loglik <- n*lgamma((v+1)/2) - n*lgamma(v/2)
是 运行 并且是有效的,所以它被存储到 loglik 中。第二行 - (n/2)*log(v*pi) - ((v+1)/2)*sum(log(1+x^2/v))
有效。它执行但不存储任何内容,因为我们不在顶层,它甚至不打印。
您要么需要确保第一行的末尾不会产生完整的语句(通过在第一行的末尾添加减号),要么采取其他措施以不同的方式将其拆分.