r - R中的非线性约束优化

Question

我正在尝试预测特定疾病的存活率。我唯一掌握的信息是诊断后 1 年、3 年、5 年和 10 年的生存率。

例如:

S，即存活率

S<-c(81,78,72,65)

x, 诊断后的时间

x<-c(1,3,5,10)

我正在尝试测试一些函数，这些函数可以让我估计 20 年后的存活率。

我的一个函数被定义为

f(x)= exp(ax^b)，a 和 b 未知但需要为正数。我非常友善地使用了 fmarm 提供的代码，但使用另一个函数进行了测试。

f(x) = (1 + (x/a)^b)^-1

但是，我得到的值非常奇怪，都低于 1，而且我似乎不明白为什么。我错过了什么吗？

S<-c(81,78,72,65)
x<-c(1,3,5,10)

f<-function(ab)
{
  a <- ab[1]
  b <- ab[2]
  return(sum((((1+(x/a)**b)**-1)-S)**2))
}

minim <- nlm(f,p=c(1,1))

ab <- minim$estimate

a_opt <- ab[1]
b_opt <- ab[2]

prediction_exp <- function(x){
  return((1+(x/a_opt)**b_opt)**-1)
}
prediction_exp(20)

plot(prediction_exp(1:20), type="l", col="blue", xlab="Nb d'années après diagnostic", ylab="survie nette en %")
lines(x,S,col="black")

P.S:我发现了我的错误。 S 向量需要小于 1，函数应该是 x*a 而不是 (x/a)。再次感谢 fmarm 帮助我!

Answer 1

在你的例子中，S 和 x 是固定的，你想找到最小化 sum(i=1 to 4) exp(a*x[i]**b)-S[i])** 的 a 和 b 2

你可以创建一个函数

f <- function(ab){
  a <- ab[1]
  b <- ab[2]
  return(sum((exp(a*x**b)-S)**2))
}  

ab 是一个长度为 2 的向量，第一个是 a，第二个是 b

要最小化此功能，您可以使用 nlm

minim <- nlm(f,p=c(0,0))

你必须给 p : ab 的起始参数。因为我不知道什么是好的，所以我只是把 a=0 和 b=0结果有一个估计组件，可以为您提供算法找到的最佳参数

ab <- minim$estimate

然后你可以从ab中提取a和b

a_opt <- ab[1]
b_opt <- ab[2]

您可以创建您的预测函数

prediction_exp <- function(x){
  return(exp(a_opt*x**b_opt))
}
prediction_exp(20)

预计20年后的生存率约为63%