Нелинейная ограниченная оптимизация в 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")

PS: Я нашел свою ошибку.Вектор S должен был быть меньше 1, а функция должна была быть x * a, а не (x / a).Еще раз спасибо за помощь мне!

Answer 1

Это просто код, предложенный в принятом ответе с правильными результатами выживания (ограниченный интервалом [0-1] и исправленными результатами:

S<-c(81,78,72,65)/100
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)
[1] 0.5975635

png(); 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") ; dev.off()

Answer 2

В вашем случае S и x являются фиксированными, и вы хотите найти a и b, которые минимизируют сумму (i = от 1 до 4) exp (a * x [i] ** b) -S [i]) ** 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. Поскольку я понятия не имел, что может быть хорошего, я просто положил а = 0 и б = 0 Результат имеет компонент оценки, который дает лучшие параметры, найденные алгоритмом

ab <- minim$estimate

Тогда вы можете извлечь a и b из ab

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

Вы можете создать свою функцию прогнозирования

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

Прогнозируемая выживаемость через 20 лет составляет приблизительно 63%