Вычисление R ^ 2 для нелинейного метода наименьших квадратов

Предположим, что у меня есть значения x значения y и ожидаемые значения y f (из некоторой нелинейной кривой наилучшего соответствия).

Как я могу вычислить R ^ 2 в R? Обратите внимание, что эта функция не является линейной моделью, но нелинейные наименьшие квадраты (nls) подходят, поэтому не подходят lm.

Ответ 1

Вы просто используете функцию lm для соответствия линейной модели:

x = runif(100)
y = runif(100)
spam = summary(lm(x~y))
> spam$r.squared
[1] 0.0008532386

Обратите внимание, что квадрат r не определен для нелинейных моделей или, по крайней мере, очень сложной цитаты из R-help:

Существует веская причина, что модель nls, соответствующая R, не обеспечивает r-квадрат. R-squared не имеет смысла для общей модели nls.

Один из способов мышления r-квадрата - это сравнение остаточной суммы квадратов для подгонной модели с остаточной суммой квадратов для тривиальной модели, состоящей только из константы. Вы не можете гарантировать, что это сравнение вложенных моделей при работе с моделью nls. Если модели не вложены, это сравнение не имеет особого смысла.

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

Если вы хотите получить рецензируемые доказательства, см. Эту статью, например; это не значит, что вы не можете вычислить значение R ^ 2, просто это может означать не то же самое/иметь те же желательные свойства, что и в случае с линейной моделью.

Ответ 2

Похоже, что f - ваши предсказанные значения. Таким образом, расстояние от них до фактических значений, разделенных n * дисперсией y

так что-то вроде

1-sum((yf)^2)/(length(y)*var(y))

должен дать вам квазиразрешенное значение, если ваша модель достаточно близка к линейной модели, а n довольно большая.

Ответ 3

Еще один квази-R-квадрат для нелинейных моделей - это квадрат корреляции между фактическими значениями y и предсказанными значениями y. Для линейных моделей это регулярное R-квадрат.

Ответ 4

В качестве альтернативы этой проблеме я использовал несколько раз следующую процедуру:

вычислить соответствие данных с помощью функции nls
используя полученную модель, делают прогнозы
Трассировка (график...) данных против значений, предсказанных моделью (если модель хорошая, точки должны быть вблизи биссектрисы).
Вычислите R2 линейной регрессии.

С наилучшими пожеланиями. Патрик.

Ответ 5

Как прямой ответ на поставленный вопрос (а не утверждать, что R2/псевдо R2 не являются полезными) в nagelkerke функцию в rcompanion пакете будет сообщать различные значения псевдо R2 для нелинейных наименьших квадратов (NLS) модели, предложенный Mcfadden, Кокс и Снелл и Нагелькерке, например

require(nls)
data(BrendonSmall)
quadplat = function(x, a, b, clx) {
          ifelse(x  < clx, a + b * x   + (-0.5*b/clx) * x   * x,
                           a + b * clx + (-0.5*b/clx) * clx * clx)}
model = nls(Sodium ~ quadplat(Calories, a, b, clx),
            data = BrendonSmall,
            start = list(a   = 519,
                         b   = 0.359,
                         clx = 2304))
nullfunct = function(x, m){m}
null.model = nls(Sodium ~ nullfunct(Calories, m),
             data = BrendonSmall,
             start = list(m   = 1346))
nagelkerke(model, null=null.model)

Пакет soilphysics также сообщает Efron pseudo R2 и скорректировал значение псевдо R2 для моделей nls как 1 - RSS/TSS:

pred <- predict(model)
n <- length(pred)
res <- resid(model)
w <- weights(model)
if (is.null(w)) w <- rep(1, n)
rss <- sum(w * res ^ 2)
resp <- pred + res
center <- weighted.mean(resp, w)
r.df <- summary(model)$df[2]
int.df <- 1
tss <- sum(w * (resp - center)^2)
r.sq <- 1 - rss/tss
adj.r.sq <- 1 - (1 - r.sq) * (n - int.df) / r.df
out <- list(pseudo.R.squared = r.sq,
            adj.R.squared = adj.r.sq)

который также является pseudo R2 рассчитанным функцией accuracy в пакете rcompanion. В основном, это R2 измеряет, насколько лучше ваше соответствие сравнивается с тем, если вы просто нарисуете плоскую горизонтальную линию через них. Это может иметь смысл для моделей nls если ваша нулевая модель является той, которая позволяет использовать только модель перехвата. Также для конкретных других нелинейных моделей это может иметь смысл. Например, для модели мошенничества, в которой используется строгий рост сплайнов (bs = "mpi" в слоте сплайна), приспособленная модель для наихудшего возможного сценария (например, где ваши данные были строго снижены) была бы плоской линией, и, следовательно, R2 равным нулю. Затем скорректированный R2 также оштрафовал модели с более высокими значениями установленных параметров. Использование скорректированного значения R2 уже рассмотрит многие критические замечания в документе, приведенном выше, http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2892436/ (кроме того, если клянутся, используя информационные критерии для выполнения выбор модели - вопрос, какой из них использовать - AIC, BIC, EBIC, AICc, QIC и т.д.).

Просто используя

r.sq <- max(cor(y,yfitted),0)^2
adj.r.sq <- 1 - (1 - r.sq) * (n - int.df) / r.df

Я думаю, что это также имеет смысл, если у вас есть нормальные гауссовские ошибки, т.е. Корреляция между наблюдаемым и установленным у (обрезанная на нуле, так что отрицательная связь будет означать нулевую прогностическую мощность) в квадрате, а затем скорректирована для nr установленных параметров в скорректированная версия. Если y и yfitted идут в одном направлении, это будет R2 и adjusted R2 значение adjusted R2 как сообщается для регулярной линейной модели. Для меня это, по крайней мере, имело бы прекрасный смысл, поэтому я не согласен с откровенным отрицанием полезности pseudo R2 значений pseudo R2 для моделей nls поскольку, как представляется, подразумевается ответ выше.

Для ненормальных структур ошибок (например, если вы использовали GAM с ненормальными ошибками) McFadden pseudo R2 определяется аналогично

1-residual deviance/null deviance

См. Здесь и здесь для некоторого полезного обсуждения.