Пакет R spatstat: как получить стандартные ошибки для негладких членов в модели процесса Пуассона (функция: ppm), когда use.gam = TRUE?

В пакете R spatstat (я использую текущую версию 1.31-0) есть опция use.gam. Когда вы устанавливаете это значение true, вы можете включать сглаженные члены в линейный предиктор, так же, как вы это делаете с пакетом R mgcv. Например,

g <- ppm(nztrees, ~1+s(x,y), use.gam=TRUE)

Теперь, если мне нужен доверительный интервал для перехвата, вы обычно можете использовать summary или vcov, что работает, когда вы не используете gam, но не работает, когда вы используете гам.

vcov(g)

что дает сообщение об ошибке

Error in model.frame.default(formula = fmla, data = 
    list(.mpl.W = c(7.09716796875,  :invalid type (list) for variable 's(x, y)'

Я знаю, что это стандартное приближение ошибки здесь не оправдано при использовании gam, но это фиксируется в сообщении предупреждение:

In addition: Warning message: model was fitted by gam();
            asymptotic variance calculation ignores this

Меня это не беспокоит - я готов оправдать использование этих стандартных ошибок для той цели, для которой я их использую - мне просто нужны числа и я бы не хотел, чтобы это делалось «самим писать».

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

Есть ли у кого-нибудь какие-нибудь советы или идеи по этому поводу? Это ошибка? К моему удивлению, мне не удалось найти обсуждение этого в Интернете. Любая помощь или подсказки приветствуются.

Изменить: Хотя я окончательно ответил на свой вопрос здесь, я пока не приму свой ответ. Таким образом, если кто-то заинтересован и готов приложить усилия, чтобы найти "обходной путь" для этого, не дожидаясь следующего выпуска spatstat, я могу наградить его / нее. В противном случае я просто приму свой ответ в конце периода вознаграждения.


r statistics spatstat
person Macro    schedule 25.02.2013    source источник
comment
Я получаю такое же сообщение об ошибке, если просто пытаюсь вывести g на экран. Похоже, проблема возникает при вызове функции model.frame. При простой отладке кажется, что ошибка возникает в строковых данных ‹- .Internal (model.frame (formula, rownames, variables, varnames, extras, extranames, subset, na.action))   -  person Jouni Helske    schedule 28.02.2013
comment
Привет, @Hemmo, я это вижу. Но модель по-прежнему оценивается (например, coef(g) работает), и вы можете построить прогнозируемые значения и т. Д. (Хотя, когда вы пытаетесь получить стандартные ошибки для прогнозов, вы снова возвращаетесь к этой ошибке). Какие-нибудь советы?   -  person Macro    schedule 28.02.2013
comment
Бегло взглянув на код ppm и mpl.engine, я бы сказал, что ppm и его подфункции не используют подход model.frame. Он сохраняет формулу и данные в выходных данных (g $ internal), но анализ формулы по умолчанию /model.frame.default не может обрабатывать список s (x, y), поскольку во фрейме данных такого нет. Я предполагаю, что это ошибка, и вам следует спросить об этом у автора пакета. Вы также можете проверить это с более старой версией пакета и посмотреть, появится ли у вас такая же ошибка.   -  person Jouni Helske    schedule 28.02.2013

Ответы (2)


arrow_upward
4
arrow_downward

Я связался с одним из авторов пакета, Адрианом Баддели, по этому поводу. Он оперативно ответил и сообщил мне, что это действительно ошибка программного обеспечения, и, по-видимому, я первый, кто с ней столкнулся. К счастью, ему потребовалось совсем немного времени, чтобы отследить проблему и исправить ее. Исправление будет включено в следующий выпуск spatstat, 1.31-1.

Изменить: выпущена обновленная версия spatstat, в которой больше нет этой ошибки:

g <- ppm(nztrees, ~1+s(x,y), use.gam=TRUE)
sqrt( vcov(g)[1,1] ) 
[1] 0.1150982
Warning message:
model was fitted by gam(); asymptotic variance calculation ignores this

Другие примечания к выпуску см. На веб-сайте spatstat. Спасибо всем, кто читал и участвовал в этой ветке!

person Macro    schedule 01.03.2013

arrow_upward
1
arrow_downward

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

The default formula, ~1, indicates the model is stationary and no trend is to be fitted.

Но вместо этого вы можете указать модель так:

g <- ppm(nztrees, ~x+y, use.gam=TRUE)   
#Then to extract the coefficientss:
>coef(g)
(Intercept)             x             y 
-5.0346019490  0.0013582470 -0.0006416421 
#And calculate their se:
vc <- vcov(g)
se <- sqrt(diag(vc))
> se
(Intercept)           x           y 
0.264854030 0.002244702 0.003609366

Есть ли в этом смысл / ожидаемый результат? Я знаю, что авторы пакета очень активно работают над r-sig-geo рассылка lsit, поскольку они помогали мне в прошлом. Вы также можете разместить свой вопрос в этом списке рассылки, но вам следует указать здесь свой вопрос, когда вы это сделаете.

person Simon O'Hanlon    schedule 28.02.2013
comment
Привет @Simon, спасибо за внимание. s(x,y) - это gam синтаксис для определения эффекта x,y непараметрической сглаживающей функцией (которая оценивается). См. документацию по игре. Обратите внимание, что модель оценивается как функциональный параметр, когда вы запускаете мой код (например, строите прогнозируемую поверхность - plot(predict(g))), но кажется, что связь с gam является неполной, поскольку SE для негладких условия недоступны. Модель, которую вы подходите, является обычной лог-линейной в x и y, и мне удалось получить эти стандартные ошибки и т. Д. - person Macro; 28.02.2013
comment
@Macro ах, мне следовало уделить больше внимания. Я изучу его дальше, чтобы узнать, могу ли я помочь, но я настоятельно рекомендую публиковать сообщения на r-sig-geo, и я настоятельно рекомендую прочитать сначала руководство по публикации. Раньше я был назначен для плохо сформированного вопроса в списке рассылки! - person Simon O'Hanlon; 28.02.2013
comment
Нелокальная документация по игре для всех, кто следит за этим - person Simon O'Hanlon; 28.02.2013