Разделение на основе модели с «двухуровневым взаимодействием» (сегментированные модели)

Нет, но да ;-)

Нет, lmtree() не может быть легко использован в качестве установщика для mob().

1. Размер внутреннего дерева (lmtree()) не фиксирован, то есть вы можете получить дерево без какого-либо раздела или со многими подгруппами, и это может сбить с толку внешнее дерево (mob()).

2. Даже если бы кто-то обошел проблему размерности или исправил ее, всегда заставляя один разрыв, потребовалось бы больше работы, чтобы настроить правильный вектор коэффициентов, матрицу оценочных функций и т. д. Это также непросто, потому что скорость сходимости (и, следовательно, вывод) отличается, если заданы точки останова (например, для двоичного фактора) или должны быть оценены (например, для ваших числовых переменных b).

3. Как вы настроили функцию fit(), внутренний lmtree() не знает, где найти b. Все, что у него есть, это числовой вектор y и числовая матрица x, но не исходные данные.

Но да, я думаю, что все эти проблемы можно решить, изменив представление с подгонки «двухслойного» дерева к подгонке «сегментированной» модели внутри дерева. У меня сложилось впечатление, что вы хотите соответствовать модели y ~ x (или a ~ b в вашем примере), где в x используется кусочно-линейная функция с дополнительной точкой останова. Если в x предполагается, что кусочно-линейная функция является непрерывной, то можно легко использовать пакет segmented. Если нет, то можно использовать strucchange. Предполагая, что вы хотите первое (поскольку вы смоделировали свои данные таким образом), я включаю рабочий пример segmented ниже (а также немного изменил ваш вопрос, чтобы отразить это).

Немного изменив имена и код, ваши данные d имеют сегментированную кусочно-линейную зависимость y ~ x с коэффициентами, зависящими от групповой переменной g.

set.seed(321)
d <- data.frame(
  x = runif (200),
  g = factor(sample(1:2, 200, replace = TRUE))
)
d$y <- jitter(ifelse(d$g == "1",
  pmax(0, 2 * d$x - 1),
  pmax(0, 4 * d$x - 1.2)
), amount = 0.1)

Затем в каждый узел дерева я могу поместить модель segmented(lm(y ~ x)), которая поставляется с подходящими экстракторами для coef(), logLik(), estfun() и т. д. Таким образом, функция мафиози проста:

segfit <- function(y, x, start = NULL, weights = NULL, offset = NULL, ...)
{
  x <- as.numeric(x[, 2])
  segmented::segmented(lm(y ~ x))
}

(Примечание: я не пробовал, будет ли segmented() также поддерживать объекты lm() с weights и offset.)

Таким образом, мы можем получить полное дерево, которое просто разбивается на g в этом базовом примере:

library("partykit")
segtree <- mob(y ~ x | g, data = d, fit = segfit)
plot(segtree, terminal_panel = node_bivplot, tnex = 2)

Практическое введение в segmented доступно в: Muggeo VMR (2008). «сегментированный: пакет R для соответствия моделям регрессии с ломаной связью». Новости R, 8 (1), 20-25. https://CRAN.R-project.org/doc/Rnews/

Для базовой методологической основы см .: Muggeo VMR (2003). «Оценка регрессионных моделей с неизвестными точками разрыва». Статистика в медицине, 22 (19), 3055-3071. DOI: 10.1002 / sim.1545