Great question! These two functions both implement natural (restricted) cubic splines, but their underlying basis functions, boundary handling, and extrapolation behavior have key differences that explain the dramatic gaps you're seeing outside your training data range. Let's break this down step by step.

Core Differences #

1. Basis Function Type #

• splines::ns() uses B-spline bases: These are local, piecewise functions that only have non-zero values between adjacent knots. They're optimized for local smoothing and are the standard spline implementation in base R's toolkit.
• rms::rcs() uses truncated power bases for restricted cubic splines: This approach builds on polynomial terms (like $x, x^2, x^3$) plus truncated terms, with an explicit constraint that the spline is linear beyond the outermost knots (the "restricted" part of its name).

2. Extrapolation Behavior (The Critical Gap) #

Both are natural splines (forced linear outside boundary knots), but how they define those boundary knots is totally different:

• ns(): By default, it uses the full range of your training data (range(x)) as boundary knots. You can override this with the Boundary.knots parameter, but out of the box, extrapolation starts from the min/max of your observed x values.
• rcs(): When paired with datadist() (required for most rms functionality), it uses the x range from your datadist object. If you manually specify knots, it treats the first and last knot as boundary points—so extrapolation starts at those knot values, not your data's actual min/max.

This is exactly why your extrapolation to x=-10 and x=10 looks so different: ns() extends linearly from your data's true extremes (~-3.5 to 3.2), while rcs() extends linearly from your specified knots (-2 and 2). Their linear slopes are completely different as a result.

3. Parameter Interpretability & Tooling #

• ns() parameters are hard to interpret directly, since B-spline bases are local and don't map to intuitive polynomial terms.
• rcs() parameters are tied to truncated power terms, which are slightly easier to reason about. Plus, the rms package gives you extra tools for spline visualization (like Predict()), model validation, and calibration that integrate seamlessly with rcs().

Verifying With Your Simulation #

Your code perfectly demonstrates these differences—let's walk through what's happening:

1. Prepare Non-Linear Data #

library(tidyverse)
library(splines)
library(rms)
set.seed(100)
xx <- rnorm(1000)
yy <- 10 + 5*xx - 0.5*xx^2 - 2*xx^3 + rnorm(1000, 0, 4)
df <- data.frame(x=xx, y=yy)

2. Fit Models (Note Boundary Knot Differences) #

# ns() uses df$x's full range (~-3.5 to 3.2) as boundary knots, plus your specified knots
ns_mod <- lm(y ~ ns(x, knots=c(-2, 0, 2)), data=df)

# rms requires datadist; rcs uses df$x's range as boundaries, with your knots as internal knots
ddist <- datadist(df)
options("datadist" = "ddist")
trunc_power_mod <- ols(y ~ rcs(x, knots=c(-2, 0, 2)), data=df)

3. Training Set Fit (Nearly Identical!) #

Within your training data range, both models perform almost exactly the same—their MSE values will be nearly identical:

mean(ns_mod$residuals^2)  # ~15.7
mean(trunc_power_mod$residuals^2)  # ~15.7

Plotting the fitted values confirms they overlap completely where you have data.

4. Extrapolation (Huge Differences) #

When you predict outside your training data, the boundary knot difference becomes obvious:

newdata <- data.frame(x=seq(-10, 10, 0.1))
pred_ns_new <- predict(ns_mod, newdata=newdata)
pred_trunc_new <- predict(trunc_power_mod, newdata=newdata)
newdata$pred_ns_new <- pred_ns_new
newdata$pred_trunc_new <- pred_trunc_new
newdata_melted <- newdata %>% gather(key="model", value="predictions", -x)
ggplot(newdata_melted, aes(x=x, y=predictions, group=model, linetype=model)) + geom_line()

You'll see ns() extends linearly from your data's true min/max, while rcs() starts its linear extrapolation at -2 and 2—leading to wildly different predictions at x=-10 or x=10.

Can They Be Interchanged? #

It depends entirely on your use case:

1. Only predicting within training data range: Yes—if you manually align the boundary knots (e.g., set Boundary.knots=c(-2,2) for ns() to match rcs()'s knot boundaries), their fit will be nearly identical.
2. Need to extrapolate: Absolutely not. Their extrapolation logic is fundamentally different, and using them interchangeably here will give you incorrect results.
3. Toolchain compatibility: If you're using other rms features (like model validation, calibration, or complex nested models), stick with rcs(). For simpler GLMs/LMs, ns() is lighter weight and integrates seamlessly with base R workflows.

Final Takeaway #

ns() and rcs() are both natural cubic spline implementations, but their boundary handling and extrapolation rules make them non-interchangeable for most real-world use cases—especially if you need to predict outside your training data. Always double-check your boundary knot settings if you're trying to align their behavior!

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