r - 如何使用等高线图和瓦片图可视化 GAM 结果(使用 ggplot2)

Question

我想通过 gam 结果用 ggplot2 绘制等高线图。下面是我想要的详细解释:

#packages
library(mgcv)
library(ggplot2)
library(tidyr)  
#prepare data
df <- data.frame(x = iris$Sepal.Width,
                 y = iris$Sepal.Length,
                 z = iris$Petal.Length)
#fit gam
gam_fit  <- gam(z ~
                  s(x) +
                  s(y),
                data=df,na.action = "na.fail")

为了根据 gam_fit 预测 z 值，我从 https://drmowinckels.io/blog/2019-11-16-plotting-gamm-interactions-with-ggplot2/ 找到了一种方法

#predict z values
df_pred <- expand_grid(
  x = seq(from=min(df$x), 
              to=max(df$x), 
              length.out = 100),
  y = seq(from=min(df$y), 
              to=max(df$y), 
              length.out = 100)
)
df_pred <- predict(gam_fit, newdata = df_pred, 
                      se.fit = TRUE) %>%  
    as_tibble() %>% 
    cbind(df_pred)
gg <- ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
  geom_point(data=df,aes(x=x, y=y))+
  scale_fill_distiller(palette = "YlGnBu")+
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white")
print(gg)

这给了我下面的情节

我的目标是移除没有测量 x-y 点的平铺和轮廓。例如，绘图的右上角和左上角周围没有测量点。

我想知道 mgcViz 是否可以实现这一点，但它需要包括 x 和 y 作为交互项，如下所示(我也不确定如何在下图中添加测量点):

library(mgcViz)
gamm_fit2  <- gam(z ~
                   s(x,y),
                data=df,na.action = "na.fail") #,REML=TRUE
b <- getViz(gamm_fit2)
plot(sm(b, 1))

我认为 df_pred 可能不是实现我目标的最佳格式，但我不确定该怎么做。如果您能给我任何关于 ggplot2 的解决方案，我将不胜感激。

Answer 1

可能有一个包被设计用来处理这个任务，但是如果你找不到合适的“工具”来完成这项工作，一个选择是围绕“点”绘制一个多边形并将多边形之外的所有东西都涂成灰色，例如

library(tidyverse)
library(mgcv)

#prepare data
df <- data.frame(x = iris$Sepal.Width,
                 y = iris$Sepal.Length,
                 z = iris$Petal.Length)
#fit gam
gam_fit  <- gam(z ~
                  s(x) +
                  s(y),
                data=df,na.action = "na.fail")

df_pred <- expand_grid(
  x = seq(from=min(df$x), 
          to=max(df$x), 
          length.out = 100),
  y = seq(from=min(df$y), 
          to=max(df$y), 
          length.out = 100)
)
df_pred <- predict(gam_fit, newdata = df_pred, 
                   se.fit = TRUE) %>%  
  as_tibble() %>% 
  cbind(df_pred)

ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
  geom_point(data=df,aes(x=x, y=y))+
  scale_fill_distiller(palette = "YlGnBu")+
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white") +
  coord_cartesian(xlim = c(1.9, 4.5),
                  ylim = c(4, 8))

# Get the 'hull' around all of the dots
hulls <- df[chull(df$x, df$y), ]
# Get the 'edges' of the frame, starting at the first hull point
edges <- data.frame(x = c(4.1,4.5,4.5,1.9,1.9,4.5),
                    y = c(5.2,4,8,8,4,4),
                    z = NA)
# Combine
draw_poly <- rbind(hulls, edges)

# Draw the plot, and overlay the gray polygon
ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
 geom_point(data=df, aes(x=x, y=y)) +
  scale_fill_distiller(palette = "YlGnBu") +
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white") +
  geom_polygon(data=draw_poly, aes(x=x, y=y), fill = "grey")

# Without the points
ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
#  geom_point(data=df, aes(x=x, y=y)) +
  scale_fill_distiller(palette = "YlGnBu") +
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white") +
  geom_polygon(data=draw_poly, aes(x=x, y=y), fill = "grey")

^{由 reprex package 创建于 2022-09-16 (v2.0.1)}

---

这是使用 concaveman package 的另一个示例计算凹包:

library(ggforce)
#install.packages("concaveman")
library(concaveman)

border <- concaveman(as.matrix(df[,1:2]), concavity = 2)
edges <- data.frame(V1 = c(4.5,4.5,1.9,1.9,4.5),
                    V2 = c(4,8,8,4,4))
draw_poly <- rbind(border, edges)

ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
  geom_point(data=df, aes(x=x, y=y)) +
  scale_fill_distiller(palette = "YlGnBu") +
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white") +
  geom_shape(data=draw_poly, aes(x=V1, y=V2), fill = "grey",
             expand = unit(-0.05, "cm"))

ggplot() +
  geom_tile(data=df_pred, aes(x=x, y=y, fill = fit)) +
 # geom_point(data=df, aes(x=x, y=y)) +
  scale_fill_distiller(palette = "YlGnBu") +
  geom_contour(data=df_pred, aes(x=x, y=y, z = fit), colour = "white") +
  geom_shape(data=draw_poly, aes(x=V1, y=V2), fill = "grey",
             expand = unit(-0.05, "cm"))

^{由 reprex package 创建于 2022-09-16 (v2.0.1)}