Сообщение об ошибке в mgcv GAM при использовании family = gaulss ()

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

Следовал рекомендациям по HGAM из недавней статьи (https://peerj.com/articles/6876/) вместе с указаниями из примечаний в описании пакета lmvar (https://cran.r-project.org/web/packages/lmvar/vignettes/Intro.html).


library(mgcv); library(datasets)

# Using the CO2 dataset for illustration

data <- datasets::CO2

# Changing the 'Plant' factor from ordered to unordered

data$Plant <- factor(data$Plant, ordered = FALSE)

# This model fits with no errors

CO2_modG_1 <- gam(log(uptake) ~ s(log(conc), k = 5, bs = 'tp') + s(Plant, k = 12, bs = 're'), data = data, method = 'REML', family = 'gaussian')

# This model fails with error

CO2_modG_2 <- gam(log(uptake) ~ s(log(conc), k = 5, bs = 'tp') + s(Plant, k = 12, bs = 're'), data = data, method = 'REML', family = gaulss())

Это возвращает сообщение об ошибке:


Error in while (sqrt(mean(dlb/(dlb + lami * dS * rm)) * mean(dlb)/mean(dlb +  : 
  missing value where TRUE/FALSE needed

r gam mgcv
person momcanally    schedule 19.08.2019    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Вы должны передать список из двух объектов формулы в gam() при подгонке гауссовой модели в масштабе местоположения с семейством gaulss().

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

CO2_modG_2 <- gam(list(log(uptake) ~ s(log(conc), k = 5, bs = 'tp') + 
                         s(Plant, k = 12, bs = 're'), 
                       ~ 1), 
                  data = data, method = 'REML', family = gaulss())

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

CO2_modG_3 <- gam(list(log(uptake) ~ s(log(conc), k = 5, bs = 'tp') + 
                         s(Plant, k = 12, bs = 're'), 
                       ~ s(Plant, k = 12, bs = 're')), 
                  data = data, method = 'REML', family = gaulss())

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

list(response ~ x1 + s(x2), # location linear predictor, left & right hand sided
              ~ x1 + s(x3)  # scale linear predictor, right-hand side only
    )

следовательно, согласно первому примеру, который я показал выше, если вам нужен постоянный член только в линейном предикторе для масштаба, вам нужно использовать нотацию R для перехватывающего или постоянного члена; ~ 1

list(response ~ x1 + s(x2), # location linear predictor, left & right hand sided
              ~ 1           # scale linear predictor, constant
    )
person Gavin Simpson    schedule 19.08.2019