Найти максимум функции в R - PullRequest
2 голосов
/ 18 июня 2019

У меня есть следующая функция.

enter image description here

Пусть F (.) - накопительная функция распределения gamma с shape = 1и rate =1.Знаменателем является функция выживания S(X) = 1 - F(X).g(x) - это средняя функция остаточного ресурса.

Я написал следующую функцию в r.

x = 5
denominator = 1 -pgamma(x, 1, 1)
numerator = function(t) (1 - pgamma(t, 1, 1))

intnum  = integrate(numerator , x, Inf)

frac = intnum$value/denominator
frac

Как найти максимум функции g(x) для всех возможных значенийX >= 0?Могу ли я сделать это в r?Большое спасибо за вашу помощь.

1 Ответ

1 голос
/ 18 июня 2019

Перед началом я определил функцию, которую вы сделали

surviveFunction<-function(x){
  denominator = 1 -pgamma(x, 1, 1)
  numerator = function(t) (1 - pgamma(t, 1, 1))

  # I used sapply to get even vector x
  intnum  = sapply(x,function(x){integrate(numerator , x, Inf)$value})

  frac = intnum/denominator
  return(frac)
}

Тогда давайте подгоним нашу функцию к функции, называемой «кривая», она нарисует график с непрерывными данными. И результат показан ниже

df = curve(surviveFunction, from=0, to=45)
plot(df, type='l')

enter image description here И отрегулируйте xlim, чтобы найти максимальное значение

df = curve(surviveFunction, from=0, to=45,xlim = c(30,40))
plot(df, type='l')

enter image description here

И теперь мы можем догадаться, что глобальный максимум находится в районе 35

.

Я предлагаю два варианта, чтобы найти глобальный максимум.

Сначала используйте данные df, чтобы найти максимум

> max(df$y,na.rm = TRUE)
 1.054248 #maximum value

> df$x[which(df$y==(max(df$y,na.rm = TRUE)))]
 35.55 #maximum value of x 

Второе использование optimize.

> optimize(surviveFunction, interval=c(34, 36), maximum=TRUE)

$maximum
[1] 35.48536

$objective
[1] 1.085282

Но функция optimize находит не глобальное максимальное значение, которое я думаю.

Если вы видите ниже

optimize(surviveFunction, interval=c(0, 36), maximum=TRUE)

$maximum
[1] 11.11381

$objective
[1] 0.9999887

Приведенный выше результат не является глобальным максимумом, я думаю, это локальный максимум.

Итак, я предлагаю вам использовать первое решение.

...