Ошибка в effectfun и parres для модели взаимодействия области Гиббса

#r #spatstat

#r #spatstat

Вопрос:

Я использую spatstat v.1.59-0 для построения некоторых моделей точечных процессов, но у меня возникают проблемы с некоторыми инструментами проверки, в частности effectfun и с остаточным графиком parres

Компонентами модели являются, rings_pp которая состоит из точечных местоположений для 52 археологических объектов и elev которая представляет собой ЦМР, преобразованную в spatstat пиксельное изображение. Я пытаюсь оценить соответствие модели точечного процесса Гиббса параметру взаимодействия области ( AreaInter ) и elev

Я подгоняю модель к данным, используя следующий код:
rr1 <- data.frame(r=seq(100, 2000, by=50))
p1 <- profilepl(rr1, AreaInter, rings_pp~elev, aic=T)
ppm5 <- as.ppm(p1)

Кажется, что все работает нормально (например, другие диагностики показывают, что модель разумно соответствует данным), за исключением того, что когда я пытаюсь использовать effectfun и parres для оценки эффекта ковариации высот, они не будут работать.

parres выдает ошибку Error in m[, d$relevant, drop = FALSE] : (subscript) logical subscript too long При обратном отслеживании ошибки я получаю следующее:

 `> trace(parres(ppm5, "elev"))
Error in m[, d$relevant, drop = FALSE] : 
  (subscript) logical subscript too long
> traceback()
4: effectFun.can(x)
3: effectFun(xxx)
2: parres(ppm5, "elev")
1: trace(parres(ppm5, "elev"))`
  

effectfun работает, когда se.fit=FALSE , но не будет работать, когда se.fit=TRUE и выдает следующую ошибку: Error in quadform(mm, vc) : ncol(x) == nrow(v) is not TRUE Когда я отслеживаю ошибку, я получаю следующее:

 `> trace(plot(effectfun(ppm5, "elev", se.fit=T)))
Error in quadform(mm, vc) : ncol(x) == nrow(v) is not TRUE
> traceback()
8: stop(simpleError(msg, call = sys.call(-1)))
7: stopifnot(ncol(x) == nrow(v))
6: quadform(mm, vc)
5: predict.ppm(orig.model, locations = fakeloc, covariates = fakecov, 
   se = se.fit)
4: predict(orig.model, locations = fakeloc, covariates = fakecov, 
   se = se.fit)
3: effectfun(ppm5, "elev", se.fit = T)
2: plot(effectfun(ppm5, "elev", se.fit = T))
1: trace(plot(effectfun(ppm5, "elev", se.fit = T)))`
  

Когда я подгоняю неоднородную модель, ppm1 <- ppm(rings_pp, ~elev) оба effectfun и parres работают нормально и показывают хорошую подгонку (хотя использование остаточной функции K Kres предполагает, что модель не учитывает кластеризацию). Похоже, это связано с тем, как я применил модель Гиббса AreaInter .

Любой совет был бы высоко оценен.

Комментарии:

1. Спасибо за хороший вопрос со многими деталями. Можете ли вы либо предоставить свои данные, либо попытаться воссоздать проблему с помощью имитируемого набора данных или с использованием одного из встроенных наборов данных? Гораздо проще отлаживать проблему, которая может быть воспроизведена локально.

2. Привет, доктор Рубак. Спасибо за ваш ответ. Мне не удалось воссоздать проблему с одним из примеров наборов данных (bei), и я не думаю, что возможно загрузить мои данные в stackoverflow. Ковариация DEM довольно большая (380 Мб). Тем не менее, я был бы рад отправить вам данные и R-скрипт по электронной почте, если вы не возражаете.

3. Если это вообще возможно, разместите общедоступную ссылку, например, на Dropbox, чтобы любой мог загрузить и протестировать пример. У меня не будет времени запустить код и отладить его на следующей неделе или больше, но я, вероятно, когда-нибудь этим займусь. Я не могу получить 380 МБ по электронной почте. Я думаю, что мой лимит составляет 10 МБ.

Ответ №1:

Еще раз спасибо за четкий вопрос.

Это ошибка, которую можно продемонстрировать на простом примере

 fit <- ppm(cells ~ x, AreaInter(0.07))
plot(effectfun(fit, se.fit=TRUE))
plot(parres(fit))
  

Я скоро проведу расследование и исправлю это.

Postscript: теперь это исправлено в версии для разработки spatstat (version 1.59-0.020 ), доступной в репозитории GitHub

.

Комментарии:

1. Привет, доктор Бэддели. Хорошо, большое спасибо! Ваша книга по R с Рубаком и Тернером просто фантастическая. Есть планы на 2-е издание?