在ggplot/R中如何对一组等高线图进行分面处理？

Question

编辑：我同意Roland的看法，我不需要在Shiny特定的内容上花太多的文字。删除了这些内容，并思考后增加了数据框在应该呈现的样子。

---

编辑_2：虽然shiny的内容与问题无关，但是我根据Roland的解决方案创建了一个例子，如果路人感兴趣可以查看我所做的内容。请耐心等待图形加载；可能会有点慢。

---

我正在尝试使用R和shiny绘制一组预测建模数据。我有四个变量，我想显示它们之间的交互作用，以轮廓图的形式展示。对于每个变量，我都要求用户定义一个范围以及一个保持值。每个变量有两种情况：

• 用作轴变量之一：该范围确定了我为该变量在模型中预测新响应值的值。
• 未直接出现在图表中：保持值用于将非特征变量设置为一个常数值，以便我的另外两个变量预测结果仅给出每个x和y组合的唯一/单个表面z值。

我遇到了处理数据的问题，以便友好地生成一组轮廓图。我理想情况下希望屏幕显示四个变量之间的6个交互作用(₄C₂)。

我基本上需要两组数据：

• 一个是用于训练模型的输入数据集的原始形式(这样我就可以使用predict(model, newData)来获得用于z值的输出列)
• 原始数据集的子集/重新排列形式，用于绘制/分组

对于友好的分组版本，这是我所需要的内容(在我的想法中；也许有更好的方法)：

| x        | y        | z | col | row |
|----------+----------+---+-----+-----|
| var1_min | var2_min | z |   1 |   1 |
| var1_min | ...      | z |   1 |   1 |
| var1_min | var2_max | z |   1 |   1 |
| ...      | ...      | z |   1 |   1 |
| var1_max | var2_min | z |   1 |   1 |
| var1_max | ...      | z |   1 |   1 |
| var1_max | var2_max | z |   1 |   1 |
|----------+----------+---+-----+-----|
| var1_min | var3_min | z |   1 |   2 |
| var1_min | ...      | z |   1 |   2 |
| var1_min | var3_max | z |   1 |   2 |
| ...      | ...      | z |   1 |   2 |
| var1_max | var3_min | z |   1 |   2 |
| var1_max | ...      | z |   1 |   2 |
| var1_max | var3_max | z |   1 |   2 |
|----------+----------+---+-----+-----|
| ...      | ...      | z |     |     |
|----------+----------+---+-----+-----|
| var3_min | var4_min | z |   3 |   2 |
| var3_min | ...      | z |   3 |   2 |
| var3_min | var4_max | z |   3 |   2 |
| ...      | ...      | z |   3 |   2 |
| var3_max | var4_min | z |   3 |   2 |
| var3_max | ...      | z |   3 |   2 |
| var3_max | var4_max | z |   3 |   2 |
|----------+----------+---+-----+-----|

以这种方式，我拥有了我的x和y值，相应预测响应的列，以及用于创建facet_grid的东西（一个2x3或3x2的facet）。

对于预测数据框架，其形式必须与我的初始预测数据相匹配，并且几乎类似于上述形式的转换/宽形式：

| var1      | var2      | var3      | var4      |
|-----------+-----------+-----------+-----------|
| var1_min  | var2_min  | var3_hold | var4_hold |
| var1_min  | ...       | var3_hold | var4_hold |
| var1_min  | var2_max  | var3_hold | var4_hold |
| ...       | ...       | var3_hold | var4_hold |
| var1_max  | var2_min  | var3_hold | var4_hold |
| var1_max  | ...       | var3_hold | var4_hold |
| var1_max  | var2_max  | var3_hold | var4_hold |
| ...       | ...       | ...       | ...       |
| var1_hold | var2_hold | var3_max  | var4_min  |
| var1_hold | var2_hold | var3_max  | ...       |
| var1_hold | var2_hold | var3_max  | var4_max  |

我将这输入模型中，以获取预测响应作为等高线图中的 z 使用。

由于我需要将变量排列成一个公共轴比例尺，跨越面板行或下面板列（可以是任意一个，不需要两个都有），所以情况也变得棘手起来。我会像这样排列组合：

| x    | y    | row | column |
|------+------+-----+--------|
| var1 | var2 |   1 |      1 |
| var1 | var3 |   2 |      1 |
| var2 | var3 |   1 |      2 |
| var2 | var4 |   2 |      2 |
| var4 | var3 |   1 |      3 |
| var4 | var1 |   2 |      3 |

现在我可以拥有三列和两行的外观，其中第一列具有共享的var1轴，第二列具有var2，第三列具有var4。

我正在考虑手动使用expand.grid来创建六个变量的唯一组合。完成后，我意识到每一行都会使用其保持值设置了两个变量，所以也许我可以创建这些六个组合的列表，然后将非保持值变量提取到两个新列中以用于绘图数据框？

有什么建议吗？

下面是一个恶劣的例子，我尝试着运用三个变量，着重研究var1和c(var2, var3)之间的交互作用：

# the min/max arguments to `seq()` are like the user-defined range
# take the second argument to `c()` is to be user-defined hold value

library(ggplot2)

var1 <- seq(0, 25, length.out = 10) # hold value = 11.1
var2 <- seq(5, 45, length.out = 10) # hold value = 17
var3 <- seq(55, 90, length.out = 10) # hold value = 72

# create combinations between var1 and var2, with var3 held
test_data <- expand.grid(var1 = var1, var2 = var2, var3 = 72)

# same, but for var1 vs. var3, with var2 held
test_data <- rbind(test_data,
    expand.grid(var1 = var1, var2 = 17, var3 = var3))

# create response; analog to using predict() in real life
test_data$resp <- (test_data$var1 + test_data$var2) / test_data$var3

# facet variable placeholder and filling in
test_data$facet <- rep("", nrow(test_data))
test_data[test_data$var2 == 17, "facet"] <- "var1 vs. var3"
test_data[test_data$var3 == 72, "facet"] <- "var1 vs. var2"

# now I melted
test_data2 <- melt(test_data, id.vars = c("var1", "resp", "facet"))

不幸的是，这给我留下了一堆情况，其中value被填充了来自var2和var3的所有保留值，因此我必须将它们删除：

test_data2 <- test_data2[test_data2$value != 72 & test_data2$value != 17, ]

现在，我能做到这个：

ggplot(test_data2, aes(x = var1, y = value, z = resp)) +
    stat_contour() + facet_grid(~ facet)

已经得到了我想要的数字范围。 现在我猜我需要一种优雅的方式来进行组合并保存值，而不会产生丑陋的结果。

这是一个更新版本，现在我知道如何在同一坐标轴上绘制行/列（因为我有两列和一行，所以需要在两个面板中将y轴设置为相同的变量var1）：

ggplot(test_data2, aes(x = value, y = var1, z = resp)) + 
    stat_contour() + facet_grid(~ facet, scales = "free_x")

- Hendy

---

很遗憾，您将此与闪亮的代码混合在一起。这对于您真正的问题并不重要。我不清楚您的输入应该是什么样子的，比如四个变量。关于您最后一段的内容：看起来您想设置facet_grid的scales参数。 - Roland

@Roland 同意，对此很抱歉。已删除与 Shiny 相关的内容。我添加了一个概念性的变量交互示例，以尝试说明数据框需要看起来像什么。关于比例尺，可以实现，但不能在上面的实际示例中使用；scales = "free_y" 对于该情况没有任何作用。我需要做的是在同一 y 轴上绘制（使用 var1，公共变量），然后使用 scales = "free_x"。我添加了那个图。 - Hendy

从 "edit 2" 中出现了无效的链接到闪亮的例子。 - r2evans

@r2evans 感谢您的提醒。我想当时我是在一个闪亮的 beta 服务器上，现在已经迁移到了 shinyapps.io。已更新！ - Hendy

Answer 1

我真的希望我理解你的意思是正确的。

我创建了自己的示例，实际上符合一个模型。

#some data
set.seed(42)
x1 <- rnorm(20)
x2 <- runif(20)
x3 <- rpois(20,10)
x4 <- rexp(20)
y <- 10 + 2*x1 + 3*x2^2 + 4*x3 +5*x4 + rnorm(20, sd=0.1)

dat <- data.frame(x1, x2, x3, x4, y)

#fit the model
fit <- lm(y~x1+I(x2^2)+x3+x4, data=dat)
summary(fit)

#ranges and fixed values
fix_x <- c(0.3, 0.4, 15, 1)
min_x <- c(-3, 0, 5, 0)
max_x <- c(3, 1, 20, 7)

#all combinations
combis <- combn(seq_len(ncol(dat)-1),2)
#number of x-values 
#(warning! don't make too large since expand.grid is used)
n <- 100

#create new data and predict for each combination
newdat <- lapply(seq_len(ncol(combis)),
                 function(i) {
                   gr <- expand.grid(seq(from=min_x[combis[1,i]],
                                         to=max_x[combis[1,i]],
                                         length.out=n),
                                     seq(from=min_x[combis[2,i]],
                                         to=max_x[combis[2,i]],
                                         length.out=n))

                   newdat <- as.data.frame(matrix(nrow=nrow(gr), ncol=ncol(dat)-1))
                   newdat[,combis[,i]] <- gr
                   newdat[,-combis[,i]] <- matrix(rep(fix_x[-combis[,i]],each=nrow(gr)), nrow=nrow(gr))

                   newdat <- as.data.frame(newdat)
                   names(newdat) <- head(names(dat),-1)

                   newdat$y <- predict(fit, newdata=newdat)

                   newdat$comb <- paste(combis[,i],collapse=" vs. ")
                   #rename so rbind works as needed
                   names(newdat)[combis[,i]] <- c("xa","xb")
                   names(newdat)[-combis[,i]] <- c(paste0("fix",letters[seq_len(ncol(dat)-3)]), "y", "comb")
                   newdat
                 })

newdat <- do.call(rbind,newdat)

#plot
library(ggplot2)
ggplot(newdat, aes(x=xa, y=xb, z=y)) + 
  stat_contour() + 
  facet_wrap(~comb, scales="free", ncol=2) +
  xlab("") +
  ylab("")

- Roland

---

你理解得非常准确。比我预期的要花费更多的精力，但我喜欢巧妙地使用了一种“索引表”的方法（当你发布时我正在尝试这个想法），并且使用facet_wrap而不是facet_grid来解决我需要指定类似于var4 vs. var1这样的东西的问题（你的解决方案使顺序无关紧要）。谢谢！ - Hendy

我没想到会花这么多精力，它吸引了我。也许有一天我自己可以用它（虽然我想知道是否真的没有一个包实现类似的东西）。你打算将你的闪亮应用程序公开吗？ - Roland

我从你的例子中学到了一个很棒的东西：我没想到可以只使用相同的列名来rbind数据框，我以为它们还必须按照相同的顺序。太棒了！不幸的是，我不能公开这个。这是为了交互式地可视化DOE数据，以帮助确定最佳条件。我可以做的是可能找到一些类似的数据，并创建一个类似于我正在处理非敏感数据的模拟。我认为这是一个很酷的想法，交互式部分应该非常棒。非常感谢你的帮助。 - Hendy

我刚刚收到了账户确认邮件，这是在提交请求以获得Shiny beta服务器访问权限后的回复。为了回报您的慷慨，我将使用您的模型和我最终使用的代码，并将其放在公共服务器上。我会在接下来的一两天内回来发布链接，这样您就可以看到它的运行情况 :) - Hendy

shiny 应用已经启动，使用您的示例模型！再次感谢。这太棒了。我创建了另一个问题，关于平铺和分面，因为我想创建一个带有叠加白色等高线的瓷砖背景，这样可能更容易看到 z 值。不确定。欢迎对可视化进行建议。 - Hendy