ordisurf vs mgcv: игровая модель - PullRequest
Новая
       8

ordisurf vs mgcv: игровая модель

2 голосов
/ 21 июня 2019

Я хотел обновить модель GAM за функцией ordisurf в R vegan (добавить случайный коэффициент). Для этого я сначала построил модель, которая должна быть идентична той, которая стоит за ordisurf, как объяснено здесь: https://www.fromthebottomoftheheap.net/2011/06/10/what-is-ordisurf-doing/

Тем не менее, модель ordisurf и GAM, которую я сделал в mgcv, дают действительно разные результаты. Пример кода с данными дюны ниже. Кто-нибудь может объяснить? Как должна выглядеть модель GAM, чтобы работать аналогично обычному серфингу? (Важно знать, прежде чем начинать улучшать его ...:)) 

ord <-metaMDS(dune, k=3) 

surf <- ordisurf(ord, dune.env$A1, choices = c(2,3))

scrs <- data.frame(scores(ord, display = "sites", choices = c(1,2,3)))

dat <- with(dune, cbind(scrs, dune.env$A1))

mod_23 <- gam(dune.env$A1 ~ s(NMDS2, NMDS3, k = 10), data = dat)

plot.gam(mod_12, se=FALSE, cex = 1, pch = 1, col="blue")

1 Ответ

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

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

Первое из них означает, что мы подгоняем модель к mgcv::gam(), используя method = 'REML', чтобы воспользоваться преимуществом выбора гладкости REML. Выбор плавности GCV (в настоящее время) используется по умолчанию в mgcv::gam(), но это не рекомендуется, и Саймон Вуд, разработчик mgcv , намекает, что какая-то будущая версия mgcv может изменить по умолчанию от method = 'GCV.Cp' до method = 'REML'. Причиной этого является то, что выбор плавности GCV наблюдался недостаточно гладким во многих приложениях. Это происходит, когда профиль функции GCV становится очень плоским относительно глобального оптимального показателя GCV. Незначительные различия в данных могут привести к тому, что будут выбраны очень разные параметры сглаживания, некоторые из которых приводят к значительному недостаточному сглаживанию. Было показано, что выбор плавности REML и ML менее подвержен этой проблеме; профиль оценки REML имеет тенденцию быть более изогнутым и иметь ярко выраженный минимум.

Второй пункт и то, почему мы сейчас используем select = TRUE, является результатом наблюдения, что штраф за гладкость в GAM, который определяет степень волнистости расчетной поверхности, влияет только на волнистые базисные функции. Базовое расширение содержит пару базовых функций, которые, с точки зрения штрафа, являются совершенно гладкими; есть две линейные, двумерные базисные функции, представляющие две линейные плоскости, которые по определению имеют нулевую кривизну (нулевую вторую производную) и, следовательно, вообще не вносят вклад в штраф за вихревость (по умолчанию в любом случае, который заключается в использовании штраф на кривизну расчетного сплайна). Конечным результатом является то, что штраф за гладкость может оштрафовать полностью до предполагаемой плоской поверхности, но не дальше. На практике это означает, что пользователи могут интерпретировать результат модели ordisurf() как линейную поверхность, даже если эта линейная поверхность не была статистически значимой (и многие пользователи просто смотрят на график, а не на основную GAM, которая сообщит их, если их самолет был значительным или нет). select = TRUE добавляет второе сглаживание к сглаживанию, которое влияет только на идеально гладкие базисные функции. Это приводит к тому, что GAM может оштрафовать / сжать модель целиком на поверхности, которую мы оцениваем. Другими словами, он может сжать поверхность до плоской горизонтальной поверхности с нулевым эффектом === нулевой взаимосвязи между конфигурацией ординации и переменной отклика.

Вместе я чувствовал, что эти параметры лучше всего защищены от ложных срабатываний пользователей.

Если вы измените свой gam() вызов на: 

mod_23 <- gam(dune.env$A1 ~ s(NMDS2, NMDS3, k = 10), data = dat,
              method = 'REML', select = TRUE)

тогда вы должны получить тот же вывод, что и ordisurf(): 

r$> mod_23                                                                       

Family: gaussian 
Link function: identity 

Formula:
dune.env$A1 ~ s(NMDS2, NMDS3, k = 10)

Estimated degrees of freedom:
0.285  total = 1.28 

REML score: 43.25057     

r$> surf                                                                         

Family: gaussian 
Link function: identity 

Formula:
y ~ s(x1, x2, k = 10, bs = "tp", fx = FALSE)

Estimated degrees of freedom:
0.285  total = 1.28 

REML score: 43.25057
...