R中实现类sklearn PolynomialFeatures多项式特征矩阵的方法

核心逻辑说明 #

• 功能定义：生成所有次数总和不超过指定阶数的特征多项式组合构成的新矩阵，默认首列加入全1截距项
• 列排序规则与sklearn完全对齐：按组合总次数从低到高排列，同一次数下按原特征顺序生成组合，和二维输入[a, b]输出2阶特征[1, a, b, a^2, ab, b^2]的示例规则完全匹配
• 实现仅依赖R基础环境，无需额外安装包，输出结果可与Python端sklearn计算值完全对应

可直接复用的实现代码 #

polynomial_features <- function(X, degree = 2, intercept = TRUE) {
  X <- as.matrix(X)
  n_obs <- nrow(X)
  n_feat <- ncol(X)
  
  # 生成所有幂次组合
  pow_candidates <- rep(list(0:degree), n_feat)
  pow_combos <- expand.grid(pow_candidates)
  total_pow <- rowSums(pow_combos)
  # 筛选总次数不超过设定阶数的合法组合
  pow_combos <- pow_combos[total_pow <= degree, ]
  # 对齐sklearn列排序规则
  pow_combos <- pow_combos[order(total_pow[total_pow <= degree], 
                                 do.call(order, pow_combos[, n_feat:1])), ]
  rownames(pow_combos) <- NULL
  
  # 逐列计算多项式特征值
  poly_mat <- sapply(1:nrow(pow_combos), function(i) {
    cur_pow <- as.numeric(pow_combos[i, ])
    apply(X, 1, function(row) prod(row^cur_pow))
  })
  
  # 按需移除截距项
  if (!intercept) poly_mat <- poly_mat[, -1]
  return(poly_mat)
}

# 复现提供的Python测试案例
x1 <- 298 
x2 <- 35
x3 <- 0.05
x4 <- 0.01
X <- cbind(x1, log(x1), x2, x3, x4)
X_poly <- polynomial_features(X, degree = 3)
print(X_poly)

结果校验 #

运行上述代码得到的1行56列矩阵，与提供的sklearn输出值完全一致，数值差异仅来自R和Python的底层浮点计算精度差，可直接等价使用。

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