Как интегрировать (AUC) модель nls и доверительный интервал Монте-Карло в R

Я пытаюсь интегрировать (разрешить область под) нелинейную функцию (от nls()) от x = 0 до бесконечности в R. Однако функция интеграции R вызывает функцию (f).

Короче, хотелось бы сделать что-то приближенное:

integrate(my.nls, lower = 0L, upper = Inf)

Но, к сожалению, my.nls на самом деле является подобранным объектом модели, а не функцией. Я рассмотрел возможность использования сглаживающего сплайна для интерполяции, а затем интегрирования полученной функции. Но я бы предпочел использовать истинную функцию nls вместо приближения. Кроме того, учитывая бесконечный характер интеграции, я должен быть очень осторожен с экстраполяцией в положительном направлении.

Если возможно, идеальным методом была бы возможность интегрировать результаты nls и других функций, например, область под смоделированным доверительным интервалом 97,5%, вычисленным на основе функции predictNLS пакета распространения.

Я новичок в R, и это только мой второй пост на SO, я думаю, поэтому, пожалуйста, простите меня, если это банальный или глупый вопрос или я совершил какой-то другой грех. Пока что неправильное использование as.function или function(){predict(my.nls()} никуда меня не привело, и я буду очень признателен за любую помощь.

Ниже приведен короткий пример, который должен служить иллюстрацией моей проблемы:

### Make up some data
x <- seq(from = 10, to = 1, length.out = 15)+(rnorm(15)+2)
y <- seq(from = 1, to = 10, length.out = 15)+(rnorm(15)+2)

### Fit an nls model, in this case, just a plain linear one.
my.nls <- nls(y~m*x+b, start = c(m=-1, b=100))

### Get confidence intervals from propagate package, might take a couple 
#seconds to run. Only serves to illustrate the type of values, the 
#function of which, I'd like to integrate (see my.preds$summary)

library(propagate) 

my.preds <- predictNLS(my.nls, newdata = data.frame("x" = x))

### Integrate (totally not right, just (hopefully) illustrating 
#the idea of what I'd like to do)

#exact.fn.auc <- integrate(my.nls, lower = 0L, upper = Inf)
#upperCI.fn.auc <- integrate(predictNLS(my.nls)$summary$Sim.97.5%, lower = 0L, upper = Inf)

PS: Я понимаю, что синтаксис в последних двух строках очень неправильный, я просто пытаюсь показать, откуда взялись бы значения, представленные функцией, если бы они были рассчитаны отдельно. Если есть какие-либо вопросы относительно того, к чему я клоню, спрашивайте, и я постараюсь перефразировать свою проблему.

PPS: Очень вероятно, что я иду к этому совершенно не в том направлении (хотя типы моделей, которым я должен соответствовать, на самом деле нелинейны (в отличие от модели, показанной выше), и я хотел бы получить область ниже mean и ее доверительные интервалы каким-либо образом), если у вас есть предложения по другим подходам, их тоже можно приветствовать. Проблема, с которой я сталкиваюсь со сплайнами, заключается в том, что мои настоящие модели становятся асимптотическими по мере приближения к y = 0, и, учитывая, что я перехожу к Inf, небольшие аберрации при экстраполяции позволяют разрешить некоторые действительно очень разные значения под кривой.


person m.s.bolton    schedule 13.02.2019    source источник


Ответы (1)


Основная проблема действительно заключается в том, что integrate нужна функция, а это не то, что вы пытались предоставить. Другая проблема, по крайней мере в этом примере, заключается в том, что интеграл расходится при увеличении до Inf.

Ограничивая внимание [0, 10], для первого случая имеем

integrate(function(p) 
  predict(my.nls, data.frame(x = p)),
  lower = 0, upper = 10)
# 102.0578 with absolute error < 1.1e-12

Во втором

integrate(function(p) 
  predictNLS(my.nls, newdata = data.frame(x = p), do.sim = FALSE)$summary$`Prop.97.5%`,
  lower = 0, upper = 10)
# 113.9549 with absolute error < 1.7e-06

где я также добавил do.sim = FALSE, чтобы не использовать Монте-Карло, поскольку это занимает немного больше времени, но вы, конечно, можете настроить параметры (например, количество итераций Монте-Карло nsim).

person Julius Vainora    schedule 13.02.2019
comment
Очень полезно и должно быть легко адаптировано к моей реальной проблеме. В своих экспериментах я был близок к тому, чтобы определить функцию из предсказания, но у меня был неправильный синтаксис. Спасибо! - person m.s.bolton; 14.02.2019