Как написать метод формулы S3 для объекта статистической модели в R

Question

У меня есть функция, которая выполняет M-тест Box для равенства ковариации матриц в многомерной линейной модели. Я хотел бы превратить его в общую функцию S3 с помощью метода формулы, который является наиболее естественным интерфейсом.

Полный текущий код: https://gist.github.com/friendly/749b5a69a067e02b87dd. Я мог бы вставить все это здесь, но, возможно, достаточно ссылки .

Я не понимаю много магии, используемой в функциях, которые получают доступ к компонентам объекта модели. Я использовал в качестве шаблона код, который я нашел в leveneTest в пакете car, который решает аналогичную проблему для одномерных моделей.

Вот быстрый тест с использованием метода по умолчанию boxM.default:

data(iris) 
res <- boxM(iris[, 1:4], iris[, "Species"]) 
res 

, которая дает желаемый результат:

>  data(iris) 
>  res <- boxM(iris[, 1:4], iris[, "Species"]) 
>  res 

     Box's M-test for Homogeneity of Covariance Matrices 

data: iris[, 1:4] 
Chi-Sq (approx.) = 140.94, df = 20, p-value < 2.2e-16 
> 

Когда я пытаюсь вызвать метод формулу boxM.formula непосредственно, она также работает, давая тот же результат, что и выше.

boxM(cbind(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width) ~ Species, data=iris) 

Однако это испытание метода boxM.lm не удается:

> iris.mod <- lm(cbind(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width) ~ Species, data=iris) 
> boxM(iris.mod) 
Error in cbind(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width) : 
    object 'Sepal.Length' not found 
> traceback() 
8: cbind(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width) 
7: eval(expr, envir, enclos) 
6: eval(predvars, data, env) 
5: model.frame.default(form, data) 
4: model.frame(form, data) at boxM.R#59 
3: boxM.formula(formula(y), data = model.frame(y), ...) at boxM.R#76 
2: boxM.lm(iris.mod) at boxM.R#2 
1: boxM(iris.mod) 
> 

Я думаю, что я понимаю, почему он не --- что-то делать с окружающей средой, для нахождения переменных в model.frame(), но не как исправить это.

Может кто-нибудь помочь?

Answer 1

Вы создали свою функцию boxM, которая может принимать объект lm в качестве входного сигнала. Реализация пытается извлечь формулу и model.frame из lm и повторно использовать их с boxM.formula.

Похоже, причина в том, что это не сработало из-за того, что model.frame(iris.mod) не возвращает оригинал data.frame, а 2-столбцовый data.frame, где 1-й столбец содержит матрицу левых переменных, а второй вектор правой части. Вы можете проверить это

class(model.frame(iris.mod)) 
dim(model.frame(iris.mod)) 
names(model.frame(iris.mod)) 
model.frame(iris.mod)[,1] 
model.frame(iris.mod)[,2] 

Поскольку model.frame(iris.mod) уже разобраны данные в вычислимой формате, вы можете применить boxM.default вместо boxM.formula когда lm объект является входом. Например, это работает:

boxM.default(Y = model.frame(iris.mod)[,1], 
      group = model.frame(iris.mod)[,2]) 

# Box's M-test for Homogeneity of Covariance Matrices 

#data: model.frame(iris.mod)[, 1] 
#Chi-Sq (approx.) = 140.94, df = 20, p-value < 2.2e-16 

Answer 2

Мой коллега, который решил это, сказал: «Вас укусят нестандартная оценка».

Вот решение, которое работает, и в более общем плане соответствует методам S3 для объектов модели. Он находит data в среде формулы модели .

boxM.lm <- function(y, ...) { 
    data <- getCall(y)$data 
    y <- if (!is.null(data)) { 
    data <- eval(data, envir = environment(formula(y))) 
    update(y, formula(y), data = data) 
    } 
    else update(y, formula(y)) 

    boxM.formula(formula(y), data=eval(data, envir = environment(formula(y))), ...) 
}