R语言factoextra包绘制PCA双图时遇环境对象不可子集化错误

问题背景 #

使用factoextra包绘制PCA双图时反复触发上述错误，数据框结构如下：

head(pcadata)
  log_response_ratio water_depth_m maxbottomtemp BOsalinity
1         -0.4514584            10      17.09915     35.066
2         -0.1408754            10      17.09915     35.066
3         -0.6519815            10      17.09915     35.066
4         -0.5682017            15      17.09915     35.066
5         -0.2082759            15      17.09915     35.066
6         -0.3052670            15      17.09915     35.066

复现代码 #

pcadata <- read.csv("pca.csv")
pcadata$water_depth_m <- as.numeric(pcadata$water_depth_m)
log_response_ratio <- complete.cases(pcadata$log_response_ratio)
pcadata <- pcadata[log_response_ratio, ]

# remove lnRR for the PCA
pcanolrr <- pcadata %>% dplyr::select(-log_response_ratio)

# create scaled df
pcadatafinal <- scale(pcanolrr, center = TRUE, scale = TRUE)

col_pca <- ncol(pcadatafinal)

pca_out <- PCA(pcadatafinal, scale.unit = FALSE, ncp = col_pca, graph = FALSE) # since already scaled, not need to do it again

# scree plot
fviz_screeplot(pca_out, addlabels = TRUE, ylim = c(0, 50))

# biplot
range(pcadata$log_response_ratio, na.rm=TRUE)
length(pcadata$log_response_ratio)

round(seq(-2, 2.25, length=1397), 4)

(var_breaks <- round(seq(-2, 2.25, length=1397), 4))

var_limits <- c(-2, 2.25)

fviz_pca_biplot(pca_out,
                geom.ind = "point",
                pointsize = 1.5,
                pointshape=21,
                col.ind = "#666666",
                stroke.ind = 1,
                fill.ind = pcadata$log_response_ratio,
                alpha.ind = 1,
                label ="var",
                labelsize = 2,
                col.var = "#666666",
                repel = TRUE,
                invisible="quali"
) +
  coord_fixed(ratio = 1) +
  labs(title="", x="Dim1", y="Dim2") +
  scale_x_continuous(breaks = scales::pretty_breaks(6), limits = xlim) +
  scale_y_continuous(breaks = scales::pretty_breaks(6), limits = ylim) +
  guides(fill = guide_colorbar(title.position = "top", title.hjust = 0.5)) +
  scale_fill_gradient2(low="#2166AC", mid="white", high="#B2182B", midpoint = median(pcadata$log_response_ratio)) +
  theme_bw() +
  theme(panel.grid.major = element_line(linewidth = 0.3, linetype = "dashed"),
        panel.grid.minor = element_blank(),
        axis.title = element_text(size = 10),
        axis.text = element_text(size = 10),
        axis.title.x = element_text(margin = margin(t=10)),
        axis.title.y = element_text(margin = margin(r=10)),
        legend.title = element_text(size=10),
        legend.text = element_text(size=10),
        legend.position = c(0.18, 0.1),
        legend.key.height = unit(2, units = "mm"),
        legend.key.width = unit(8, units = "mm"),
        legend.direction = "horizontal",
        aspect.ratio = 1
  )

错误信息 #

Error in range[1] : object of type 'environment' is not subsettable

回溯信息：

17.
oob(...)
16.
self$oob(x, limits)
15.
map(..., self = self)
14.
scales[[i]][[method]](data[[var]][scale_index[[i]]])
13.
FUN(X[[i]], ...)
12.
lapply(seq_along(scales), function(i) {
scales[[i]][[method]](data[[var]][scale_index[[i]]])
})
11.
FUN(X[[i]], ...)
10.
lapply(vars, function(var) {
pieces <- lapply(seq_along(scales), function(i) {
scales[[i]][[method]](data[[var]][scale_index[[i]]])
}) ...
9.
scale_apply(layer_data, x_vars, "map", SCALE_X, self$panel_scales_x)
8.
FUN(X[[i]], ...)
7.
lapply(data, function(layer_data) {
match_id <- match(layer_data$PANEL, layout$PANEL)
x_vars <- intersect(self$panel_scales_x[[1]]$aesthetics,
names(layer_data)) ...
6.
map_position(..., self = self)
5.
layout$map_position(data)
4.
ggplot_build.ggplot(x)
3.
ggplot_build(x)
2.
print.ggplot(x)
1.
(function (x, ...)
UseMethod("print"))(x)

问题原因与解决方案 #

核心原因 #

代码中scale_x_continuous和scale_y_continuous的limits参数调用了未定义的xlim和ylim，而R中xlim是ggplot2的内置函数（属于环境对象），直接尝试将其作为参数传递会触发"无法子集环境对象"的错误。

解决步骤 #

方法1：定义PCA维度的坐标范围

从PCA结果中提取样本点在维度1和维度2的坐标范围，可适当扩展范围提升图形美观度：

# 在绘制双图前添加以下代码
xlim <- range(pca_out$ind$coord[, 1]) * 1.1  # 扩展10%范围
ylim <- range(pca_out$ind$coord[, 2]) * 1.1

方法2：删除limits参数

如果不需要手动限制坐标范围，直接删除scale_x_continuous和scale_y_continuous中的limits = xlim、limits = ylim参数，factoextra会自动适配合适的坐标范围。

修改后的完整双图代码片段

# 先定义xlim和ylim
xlim <- range(pca_out$ind$coord[, 1]) * 1.1
ylim <- range(pca_out$ind$coord[, 2]) * 1.1

# 绘制双图
fviz_pca_biplot(pca_out,
                geom.ind = "point",
                pointsize = 1.5,
                pointshape=21,
                col.ind = "#666666",
                stroke.ind = 1,
                fill.ind = pcadata$log_response_ratio,
                alpha.ind = 1,
                label ="var",
                labelsize = 2,
                col.var = "#666666",
                repel = TRUE,
                invisible="quali"
) +
  coord_fixed(ratio = 1) +
  labs(title="", x="Dim1", y="Dim2") +
  scale_x_continuous(breaks = scales::pretty_breaks(6), limits = xlim) +
  scale_y_continuous(breaks = scales::pretty_breaks(6), limits = ylim) +
  guides(fill = guide_colorbar(title.position = "top", title.hjust = 0.5)) +
  scale_fill_gradient2(low="#2166AC", mid="white", high="#B2182B", midpoint = median(pcadata$log_response_ratio)) +
  theme_bw() +
  theme(panel.grid.major = element_line(linewidth = 0.3, linetype = "dashed"),
        panel.grid.minor = element_blank(),
        axis.title = element_text(size = 10),
        axis.text = element_text(size = 10),
        axis.title.x = element_text(margin = margin(t=10)),
        axis.title.y = element_text(margin = margin(r=10)),
        legend.title = element_text(size=10),
        legend.text = element_text(size=10),
        legend.position = c(0.18, 0.1),
        legend.key.height = unit(2, units = "mm"),
        legend.key.width = unit(8, units = "mm"),
        legend.direction = "horizontal",
        aspect.ratio = 1
  )

内容的提问来源于stack exchange，提问作者Liesl Blackman