Поверхность регрессии участка

Question

Я читаю книгу Коэна, Коэна, Айкена и Уэста (2003) «Прикладной множественный регрессионный корреляционный анализ для поведенческих наук» и наткнулся на трехмерный график поверхности регрессии, показывающий взаимодействие и отсутствие взаимодействия (стр. 259).Графики выглядят так, как будто они были созданы с использованием R. Мне нравятся графики в качестве учебного пособия, и я хотел бы воспроизвести их.Графики выглядят примерно так:

Единственное дополнение к Coehn et al.графики были линиями на плоскостях в среднем, + 1sd и = 1sd для x2.Это было бы отличным дополнением, если это возможно (как правило, большинство вещей возможно с R)

Я предоставил примерный набор данных ниже с IV, 2 предикторами и центрированными предикторами.Как бы я использовал R для создания графика поверхности (плоскости) регрессии, показывающего взаимодействие и аддитивную модель как для центрированных, так и для нецентрированных данных (я предполагаю, что техника будет одинаковой, но хочу убедиться).

Всего 4 сюжета: 1. нецентрированное отсутствие взаимодействия 2. нецентрированное взаимодействие 3. центрированное отсутствие взаимодействия 4. центрированное взаимодействие

DF<-structure(list(y = c(-1.22, -1.73, -2.64, -2.44, -1.11, 2.24, 
3.42, 0.67, 0.59, -0.61, -10.77, 0.93, -8.6, -6.99, -0.12, -2.29, 
-5.16, -3.35, -3.35, -2.51, 2.21, -1.18, -5.21, -7.74, -1.34), 
    x1 = c(39.5, 41, 34, 30.5, 31.5, 30, 41.5, 24, 43, 39, 25.5, 
    38.5, 33.5, 30, 41, 31, 25, 37, 37.5, 24.5, 38, 37, 41, 37, 
    36), x2 = c(61L, 53L, 53L, 44L, 49L, 44L, 57L, 47L, 54L, 
    48L, 46L, 59L, 46L, 61L, 55L, 57L, 59L, 59L, 55L, 50L, 62L, 
    55L, 55L, 52L, 55L), centered.x1 = c(5.49702380952381, 6.99702380952381, 
    -0.0029761904761898, -3.50297619047619, -2.50297619047619, 
    -4.00297619047619, 7.49702380952381, -10.0029761904762, 8.99702380952381, 
    4.99702380952381, -8.50297619047619, 4.49702380952381, -0.50297619047619, 
    -4.00297619047619, 6.99702380952381, -3.00297619047619, -9.00297619047619, 
    2.99702380952381, 3.49702380952381, -9.50297619047619, 3.99702380952381, 
    2.99702380952381, 6.99702380952381, 2.99702380952381, 1.99702380952381
    ), centered.x2 = c(9.80357142857143, 1.80357142857143, 1.80357142857143, 
    -7.19642857142857, -2.19642857142857, -7.19642857142857, 
    5.80357142857143, -4.19642857142857, 2.80357142857143, -3.19642857142857, 
    -5.19642857142857, 7.80357142857143, -5.19642857142857, 9.80357142857143, 
    3.80357142857143, 5.80357142857143, 7.80357142857143, 7.80357142857143, 
    3.80357142857143, -1.19642857142857, 10.8035714285714, 3.80357142857143, 
    3.80357142857143, 0.803571428571431, 3.80357142857143)), .Names = c("y", 
"x1", "x2", "centered.x1", "centered.x2"), row.names = c(NA, 
25L), class = "data.frame")

Заранее спасибо.

РЕДАКТИРОВАТЬ: Следующий код строит график, но не будет работать, когда у вас есть взаимодействие (это действительно то, что меня интересует).Кроме того, я не знаю, как построить максимум (+ 1sd), минимум (-1sd) и среднее значение для x2.

x11(10,5)
 s3d <- scatterplot3d(DF[,c(2,3,1)], type="n", highlight.3d=TRUE,
      angle=70, scale.y=1, pch=16, main="scatterplot3d")

   # Now adding a regression plane to the "scatterplot3d"
    my.lm <-  with(DF, lm(y ~ x1 + x2))
 s3d$plane3d(my.lm, lty.box = "solid")

Попытка построить плоскость взаимодействия (см. Здесь):

s3d <- scatterplot3d(DF[,c(2,3,1)], type="n", highlight.3d=TRUE,
      angle=70, scale.y=1, pch=16, main="scatterplot3d")

    my.lm <-  with(DF, lm(y ~ x1 + x2 + x1:x2 ))
 s3d$plane3d(my.lm, lty.box = "solid")

выдал следующую ошибку:

Error in segments(x, z1, x + y.max * yx.f, z2 + yz.f * y.max, lty = ltya,  : 
  cannot mix zero-length and non-zero-length coordinates

Answer 1

Вот как я бы это сделал (добавив немного цвета) с пакетами 'rms' и 'lattice':

require(rms)  # also need to have Hmisc installed
require(lattice)
ddI <- datadist(DF)
options(datadist="ddI")
lininterp <- ols(y ~ x1*x2, data=DF)
 bplot(Predict(lininterp, x1=25:40, x2=45:60), 
       lfun=wireframe,  # bplot passes extra arguments to wireframe
       screen = list(z = -10, x = -50), drape=TRUE)

И модель отсутствия взаимодействия:

 bplot(Predict(lin.no.int, x1=25:40, x2=45:60), lfun=wireframe, col=2:8, drape=TRUE, 
screen = list(z = -10, x = -50),
 main="Estimated regression surface with no interaction")