Блог

Как вы находите набор значений для предикторов модели (смесь линейных и нелинейных), которые дают наибольшее значение для ответа.

Пример модели:

 library(lme4); library(splines)

summary(lmer(formula = Solar.R ~ 1   bs(Ozone)   Wind   Temp   (1 | Month), data = airquality, REML = F))
 

Здесь меня интересует, какие условия (предикторы) приводят к наибольшему солнечному излучению (результату).

Этот вопрос кажется простым, но мне не удалось найти хороший ответ с помощью Google.

Если бы модель была простой, я мог бы использовать производные, чтобы найти максимум или минимум. Кто-то предположил, что если функция модели может быть извлечена, stats::optim() функция может быть использована. В качестве последнего средства я мог бы смоделировать все разумные варианты входных значений, подключить их к predict() функции и искать максимальное значение.

Последний упомянутый подход кажется не очень эффективным, и я полагаю, что это достаточно распространенная задача (например, поиск оптимальных клиентов для рекламы), для решения которой кто-то создал некоторые инструменты. Любая помощь приветствуется.

Здесь есть некоторые концептуальные проблемы.

• для простых терминов ( Wind и Temp ) ответ является линейной (и, следовательно, как монотонной, так и неограниченной) функцией предикторов. Таким образом, если эти термины имеют положительные оценки параметров, увеличение их значений до бесконечности ( Inf ) даст вам бесконечный ответ ( Solar.R ); значения должны быть как можно меньше (отрицательное бесконечное), если коэффициенты отрицательны. С практической точки зрения, тогда вы хотите установить для этих предикторов минимальное или максимальное разумное значение, если оценки параметров соответственно отрицательные или положительные.
• что касается bs термина, я не уверен, какие свойства B-сплайна находятся за пределами граничных узлов, но я почти уверен, что кривые уходят в положительную или отрицательную бесконечность, так что у вас та же проблема. Однако для случая bs также возможно, что существует один или несколько внутренних максимумов. В этом случае я бы, вероятно, попытался извлечь базовые термины и оценить сплайн по диапазону данных …

В качестве альтернативы, ваше упоминание optim заставляет меня думать, что это возможно:

 data(airquality)
library(lme4)
library(splines)
m1 <- lmer(formula = Solar.R ~ 1   bs(Ozone)   Wind   Temp   (1 | Month),
           data = airquality, REML = FALSE)
predval <- function(x) {
    newdata <- data.frame(Ozone=x[1],Wind=x[2],Temp=x[3])
    ## return population-averaged prediction (no Month effect)
    return(predict(m1, newdata=newdata, re.form=~0))
}
aq <- na.omit(airquality)
sval <- with(aq,c(mean(Ozone),mean(Wind),mean(Temp)))
predval(sval)
opt1 <- optim(fn=predval,
      par=sval,
      lower=with(aq,c(min(Ozone),min(Wind),min(Temp))),
      upper=with(aq,c(max(Ozone),max(Wind),max(Temp))),
      method="L-BFGS-B",  ## for constrained opt.
      control=list(fnscale=-1))  ## for maximization
## opt1
## $par
## [1] 70.33851 20.70000 97.00000
## 
## $value
## [1] 282.9784
 

Как и ожидалось, это промежуточное значение в диапазоне значений озона (1-168) и min/ max для ветра (2,3-20,7) и температуры (57-97).

Это решение методом перебора можно было бы сделать намного более эффективным, автоматически выбирая минимальные / максимальные значения для простых терминов и оптимизируя только сложные (полиномиальные / сплайновые / и т.д.) Термины.