如何在R语言中实现Lagrange插值公式函数？考虑使用双循环

实现思路 #

拉格朗日插值的核心是构造基函数 ( l_i(x) = \prod_{\substack{j=1 \ j \neq i}}^n \frac{x - x_j}{x_i - x_j} )，最终插值结果为 ( L(x) = \sum_{i=1}^n y_i \cdot l_i(x) )。这里用两层循环实现：

• 第一层循环：遍历所有需要插值的新点，逐个计算其插值结果
• 第二层循环：遍历每个已知数据点，计算对应基函数的值并加权累加

完整代码 #

lagrange_interp <- function(x, y, x_new) {
  # 输入合法性校验
  if (length(x) != length(y)) {
    stop("x和y的长度必须一致")
  }
  if (any(duplicated(x))) {
    stop("x中不能包含重复的自变量值")
  }
  
  n <- length(x)
  m <- length(x_new)
  y_interp <- numeric(m)  # 初始化插值结果向量
  
  # 第一层循环：处理每个待插值点
  for (k in seq_along(x_new)) {
    current_x <- x_new[k]
    total <- 0
    
    # 第二层循环：计算每个拉格朗日基函数并加权求和
    for (i in seq_along(x)) {
      # 计算基函数分子：(current_x - x_j)，j≠i
      numerator <- current_x - x[-i]
      # 计算基函数分母：(x_i - x_j)，j≠i
      denominator <- x[i] - x[-i]
      # 计算基函数值
      l_i <- prod(numerator / denominator)
      
      total <- total + y[i] * l_i
    }
    
    y_interp[k] <- total
  }
  
  return(y_interp)
}

测试示例 #

# 构造测试数据：y = x²的几个点
x <- c(1, 2, 3, 4)
y <- x^2

# 待插值的点
x_new <- c(2.5, 3.2)

# 执行插值
result <- lagrange_interp(x, y, x_new)

# 输出结果
cat("插值结果：", result, "\n")
# 预期输出：插值结果：6.25 10.24

代码说明 #

• 输入校验：先检查x和y长度是否一致，以及x是否有重复值（避免分母为0的情况）
• 两层循环逻辑：外层循环遍历每个待插值点，内层循环计算每个已知点对应的基函数，通过向量运算（x[-i]）替代了第三层循环，既满足两层循环的要求，又保证了效率
• 结果返回：将所有插值结果存储在向量中返回，适配多个待插值点的场景

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