R中同一画布绘制自定义lmer模型分组双3D曲面图的实现方法

核心逻辑为手动生成两个分组的预测曲面数据，再叠加绘制即可，不需要依赖绘图包的内置拟合功能，操作逻辑如下：

• 第一步：训练你自定义的lmer混合效应模型
• 第二步：生成覆盖两个预测变量取值范围的规则网格，分别代入分组变量的两个水平，通过训练好的模型的predict函数计算每个网格点的响应值
• 第三步：用plot3D包的persp3D函数叠加两个曲面，可自定义颜色、透明度区分分组，也可按需叠加原始散点数据。

完整可运行示例代码如下，可直接基于你的现有数据调整：

library(lme4)
library(car)
library(dplyr)
library(plot3D)

# 1. 生成示例数据（和你提供的逻辑一致）
X <- seq(76, 135) + rnorm(sd = 2, n = 60) 
Y <- seq(65, 365, length.out = 60) + rnorm(sd = 4, n = 60)
Test.grid <-  expand.grid(X = X, Y = Y)
Test.grid$A <- 1
Test.grid$Z <- 2*X + 0.5*Y 
df1 <- sample_n(Test.grid, 60)
df2 <- df1 %>% mutate(A = 2, Y = Y + 50, Z = Z + 10)
Test <- rbind(df1, df2)

# 2. 训练自定义lmer模型（随机效应部分可按需补充你的原有参数）
model <- lmer(Z ~ poly(X,2)*poly(Y,2)*as.factor(A) + (1|A), data = Test)

# 3. 生成预测网格
pred_grid <- expand.grid(
  X = seq(min(Test$X), max(Test$X), length.out = 50),
  Y = seq(min(Test$Y), max(Test$Y), length.out = 50),
  A = c(1,2) # 两个分组水平
)
pred_grid$pred_Z <- predict(model, newdata = pred_grid, allow.new.levels = T)

# 4. 拆分两个分组的预测数据，转为曲面绘制要求的矩阵格式
pred_A1 <- filter(pred_grid, A == 1)
z_mat1 <- matrix(pred_A1$pred_Z, nrow = 50, ncol = 50)
pred_A2 <- filter(pred_grid, A == 2)
z_mat2 <- matrix(pred_A2$pred_Z, nrow = 50, ncol = 50)
x_mat <- matrix(pred_A1$X, nrow = 50, ncol = 50)
y_mat <- matrix(pred_A1$Y, nrow = 50, ncol = 50)

# 5. 叠加绘制两个曲面
# 先绘制第一个分组的曲面
persp3D(x_mat, y_mat, z_mat1, col = adjustcolor("red", 0.6), border = NA, xlab = "pred1", ylab = "pred2", zlab = "response", ticktype = "detailed")
# 叠加第二个分组的曲面，add=T表示在原有画布上叠加
persp3D(x_mat, y_mat, z_mat2, col = adjustcolor("black", 0.6), border = NA, add = T)
# 可选：叠加原始散点区分分组
points3D(Test$X, Test$Y, Test$Z, col = ifelse(Test$A ==1, "darkred", "gray20"), pch = 16, add = T)

如果需要类似car包scatter3d的交互式操作效果，可以替换为rgl包的persp3d函数，同样使用add = T参数叠加曲面即可，逻辑和上述代码一致。如果需要调整曲面平滑度，修改预测网格的length.out参数即可，数值越大曲面越平滑。

内容的提问来源于stack exchange，提问作者T.Tree