企业官网建设流程全解析

做临床研究和公卫数据分析的朋友，几乎都遇到过这个场景：把一个连续暴露变量直接扔进回归模型，线性项不显著；换成分位数分组，趋势检验却又很漂亮。问题出在哪？大概率就是暴露与结局之间压根不是直线关系。这时候，限制立方样条（restricted cubic spline，RCS）几乎是公认的最优解之一。这几年你在SCI论文里看到的那些S形、U形、J形剂量反应曲线，十有八九都是用RCS做的。

我最早接触RCS，是给一位临床科室的医生做数据分析，他们想搞清楚某炎症指标和住院死亡风险之间的关系。一开始用logistic回归把指标当作连续变量，线性项P=0.18，不显著；把指标按四分位数分组，最高组和最低组相比P<0.001；分组趋势检验也显著。审稿人反问了一句：你凭什么确定它是线性的？这个反问直接把我送进了RCS的坑。这篇文章就把我在Stata里用RCS的完整方法、P值解读逻辑和踩过的坑摊开写清楚，尽量让第一次用的人少走弯路。

1. RCS的原理：为什么非线性关系需要限制立方样条

1.1 线性假设的陷阱：一个常见场景

假设你在研究BMI和死亡风险之间的关系。如果把BMI作为连续变量放进Cox回归，得到一个风险比HR=0.98，P=0.06，你可能会说“BMI和死亡风险没有显著关联”。但学过一点流行病学就知道，BMI和死亡的经典关系是J形或U形的：偏瘦和肥胖风险都高，正常体重风险最低。用一根直线去拟合这个U形，直线斜率可能真的不显著，甚至接近0，这就会造成“检不出关联”的尴尬。

更麻烦的是，如果你把BMI分组，比如按四分位数，结果往往非常显著，但分组本身带来两个问题：一是分组切点怎么选，选得不一样，结论可能也不一样，这是“切点依赖”；二是分组后丢失了BMI作为连续变量的剂量信息，你只能知道四个组的相对风险，无法知道BMI从22到23再到24，每增加一个单位风险变化趋势是什么样的。分组模型的拟合效果并不差，但它像一个粗筛子，能筛出有没有趋势，却筛不出曲线的具体形态。

这时候需要的是一个“让数据自己决定形状”的模型。它不应该被强加成一条直线，但又要有足够的稳定性，不能为了拟合数据而出现剧烈抖动。限制立方样条就是为这种需求设计的。

1.2 限制立方样条到底是怎么工作的

样条的本质是分段多项式。把BMI的取值范围切成几段，每一段用一条三次多项式去拟合，段与段之间保证连接点处的函数值、一阶导数、二阶导数连续，这样曲线看起来就是一条非常平滑的整体曲线，而不是几段拼凑出来的折线。这些连接点叫“结点（knot）”。

RCS的特殊之处有两个。第一，每个分段都是三次polynomial，所以它能灵活表示出U形、J形、S形这些复杂形状；第二，也是最关键的，它在第一结点之前和最末结点之后不再保持三次方形态，而是被强制约束成线性。这就是“限制”二字的含义。末尾两边的数据量少、信息稀疏，如果不强约束成线性，三次多项式会在数据两端乱甩，曲线尾部可能出现明显夸大的上翘或下坠。限制成线性之后，尾部更稳定，外推时的误差也小很多。

用RCS建模时，需要先确定结点数量和位置，然后根据结点位置生成基函数。假设你选了k个结点，那么一共会生成k-1个基变量。这些基变量本质上是原始暴露变量的某种变换，把它们全部放进回归模型，就得到了一个非线性拟合。在Stata中，最常见的rc_spline命令就是自动完成这件事的。

举个具体例子：如果你选4个结点，RCS会生成3个基变量，记为bmi1、bmi2、bmi3。其中bmi1其实就是原始BMI本身，代表线性部分；bmi2和bmi3是经过样条变换后的非线性校正项。模型变成：

logit(P) = b0 + b1*bmi1 + b2*bmi2 + b3*bmi3 + 其他协变量

如果b2和b3的系数都是0，那么模型自动退化为普通线性模型。所以检验非线性，本质上就是检验这两个非线性基变量的系数是否同时为0。

1.3 “限制”两个字，为什么不能省略

很多人第一次听到RCS会在脑子里想：既然三次多项式拟合能力强，是不是结点越多越好，限制不限制也无所谓？我的经验是：如果没有尾部线性约束，模型在数据边缘的表现往往非常难看。

举个例子，假设你的BMI数据主要集中在20到35之间，但有一小部分人BMI高达42。没有限制的普通三次样条为了拟合那少数高BMI个体的数据，可能在35到42这段区间上做出一个非常陡峭的上升或下降，等到了42以上又开始反向回头。这种振荡在少数异常值的驱动下会非常夸张。加了限制约束之后，最末结点以外的部分保持线性，相当于在尾部收住了缰绳，曲线不会乱跑。

当然，限制也有代价。如果暴露和结局之间的真实关系在极端区域是非线性的，RCS会低估这种尾部非线性，因为尾部被强制拉直了。但实际数据分析中，极端尾部往往是数据稀疏区，把效应强行拟合出来也没有可信度，不如承认“我不知道尾部是什么样的，只能按近似线性来推断”。所以RCS的尾部线性约束在大多数医学研究里是合理且推荐的。

1.4 与几种替代方案的对比

做非线性建模不只有RCS一种方案。多项式回归、分数多项式、普通样条、分组哑变量都是选项。我用一张表把它们的优缺点理清楚，这也是我给不同场景做选型时的依据。

从我处理过的实际项目来看，RCS在医学与公共卫生论文里的接受度最高。审稿人看到“restricted cubic spline”比看到“quadratic term”更放心，因为RCS绕开了全局多项式那种“为了弯曲所有地方都被拉弯”的尴尬，同时比分组分析保留了更多信息。

2. Stata全流程实操：从生成样条变量到检验非线性

2.1 准备环境与安装命令

RCS在Stata里的实现有好几条路。官方命令mkspline可以做，但更常用的是Stata Journal发布的rc_spline命令。这个命令的最大好处是自动把结点放在指定百分位位置，并且直接生成一组基变量，后续建模不需要再手工构造。

安装方法很简单：

ssc install rc_spline

如果网络环境下SSC源不可用，也可以直接从Stata Journal软件包安装，本质上是一个文件。装完之后可以用which rc_spline验证一下。

我的建议是：如果纯粹做RCS，直接用rc_spline即可；如果你还想同时控制其他样条变换或者需要自定义结点位置，再考虑用mkspline的cubic选项。

2.2 数据准备与变量检查

建模之前，先确认数据结构是否满足RCS的基本要求。RCS本身对数据没有特殊要求，它只是给连续自变量做了一组变换，因变量可以是二分类、生存数据和连续变量。但有几个点需要提前检查。

第一，暴露变量必须是数值型的连续变量，不要先做标准化或中心化，中心化应该放到结果显示阶段再做。第二，检查缺失值。RCS基变量的计算是基于结点百分位的，如果暴露变量有较多缺失，最好先评估缺失机制，别盲目插补。第三，确认协变量类型。分类变量建议用i.前缀，Stata会自动生成哑变量，这样后续margins和检验都方便。

以一个示例数据集来说明，假设变量有bmi（连续）、age（连续）、sex（1男0女）、death（0存活1死亡）、time（随访月数）：

use "cohort.dta", clear summarize bmi age tab death sex

看连续性变量的分布范围和分位数，这决定了后续结点的选择是否合理。如果BMI的P50在24附近，而P5在18、P95在35，那么RCS结点大概率会落在18到35之间，曲线在这些位置之间的形态估计最可靠。

2.3 生成RCS基变量的两种方式

使用rc_spline生成基变量，最核心的参数是结点数nk()。以4结点为例：

* 生成基于BMI的RCS基变量 rc_spline bmi, nk(4) * 查看新生成了哪些变量 describe bmi1 bmi2 bmi3

运行之后，Stata会默认在BMI的第5、第35、第65、第95百分位放置结点，并生成bmi1、bmi2、bmi3三个变量。其中bmi1严格等于原始BMI。你可以用list bmi bmi1 bmi2 bmi3 in 1/10验证一下。

如果你希望自己指定结点位置，可以在rc_spline里用knots()选项，比如：

rc_spline bmi, knots(18 22 25 30 35)

但我个人很少在第一次分析时手工指定结点，除非有明确的临床临界值依据。原因很简单：手工指定结点容易引入主观偏差，审稿人也容易追问“为什么选这几个切点”。使用百分位自动放置结点，至少是多数文献认可的标准做法。

另一种生成方式是官方命令mkspline：

mkspline bmi_rcs = bmi, cubic nknots(4)

这个命令也会生成bmi_rcs1、bmi_rcs2、bmi_rcs3等基变量。两个命令生成的基函数数值不完全一样，但拟合出的曲线形状是等价的。你选哪个都行，但报告时要写清楚用的是什么命令，方便复现。

2.4 拟合回归模型并输出核心检验

生成基变量之后，接下来的用法和普通回归没有本质区别。先看Logistic回归的情形：

* Logistic回归：结局为二分类death logistic death bmi1 bmi2 bmi3 age i.sex estimates store m_rcs

如果想做Cox回归，需要先设置生存数据：

stset time, failure(death==1) stcox bmi1 bmi2 bmi3 age i.sex

模型报告里会出现bmi1、bmi2、bmi3三个系数。单独看任何一个系数都没太大意义，因为RCS的效应是由三个基变量联合表达的。必须做联合检验。

整体关联P值：

test bmi1 bmi2 bmi3

这个检验的原假设是三个系数均为0，等价于“BMI与死亡风险没有任何关联”。这个P值是对整个暴露效应的总体检验，可以理解为“BMI整体上是否和结局有关”。

非线性P值：

test bmi2 bmi3

这个检验的原假设是“非线性基变量的系数为0”，如果P<0.05，说明仅仅用线性项不够，曲线存在明显的非线性变化。这是RCS分析里最关键的P值。

我见过不少人在这一步犯错，把test bmi1 bmi2 bmi3当作非线性检验来报告，实际上报告的是整体关联。这两个P值含义完全不同，后面会单独讲清楚。

3. P值解读的完整逻辑：整体关联P、非线性P与临床意义

3.1 两个P来自哪里

RCS分析中会出现至少两个P值：一个是整体关联检验（global association test），一个是非线性检验（non-linearity test）。它们的区别可以用一个生活类比来理解：你去医院做体检，测得血压值偏高，医生需要回答两个层面的问题——第一，你的血压异常是否和心血管风险有关？第二，这种关系是线性趋势还是非线性曲线？整体P回答第一个问题，非线性P回答第二个问题。

从统计原理上看，test bmi1 bmi2 bmi3是在检验所有RCS基变量的联合显著性，相当于在检验一个自由度较高的模型是否比空模型好。test bmi2 bmi3则是比较完整RCS模型和只保留线性项的简化模型，它检验的是“增加非线性项是否显著提升了拟合”。

实际操作中，如果用似然比检验来做非线性检验，可以这样跑：

logistic death bmi1 bmi2 bmi3 age i.sex estimates store full logistic death bmi age i.sex lrtest full

最后一行出来的似然比检验P值，和test bmi2 bmi3的Wald检验P值略有差异，但结论通常一致。文章里报告哪种都可以，建议全文统一，不要一会儿Wald一会儿LRT。

3.2 一个完整的案例解读

假设我在某队列数据里分析BMI与全因死亡的关系。样本量1.2万人，随访8年，死亡事件1800例。使用4结点RCS拟合Cox回归，得到以下结果：

结论应该是：BMI与全因死亡存在显著关联（整体P<0.001），且这种关联不是简单的线性关系（非线性P=0.006），曲线形态需要按非线性方式解读。然后结合图形看具体形状，是U形还是J形，再报告关键节点的HR。

如果反过来，整体P<0.001，但非线性P=0.42，应该怎么报告？这时应该写“在本次数据中，BMI与死亡风险呈线性正相关，未观察到显著的非线性趋势”。曲线应该基本接近直线，不需要强行描述成“倒J形”之类。很多人在线性P不显著时试图在图形里找出弯曲感，这是过度解读。

还有一种情况：整体P=0.20，非线性P=0.01。这看起来有点矛盾，但实际可能出现。因为非线性基变量的检验探测的是曲线形状，而整体P检验的是总效应的存在性。一个典型的S形曲线可能左右两侧效应方向相反，平均值抵消，导致整体关联不显著，但非线性项显著。这种情况下，能不能写“没有关联”？我的建议是谨慎。你应该看曲线图，如果曲线确实在某个区域穿过风险基线，就需要表述为“并未发现总体上的单调正相关，但曲线显示非单调变化，需在其他人群中验证”。不要一看到整体P>0.05就直接写“无关联”。

3.3 什么时候可以认定非线性关系成立

我的标准是三条同时满足：第一，非线性P<0.05；第二，曲线图呈现肉眼可辨的弯曲，不是只有统计上显著但视觉上近乎直线；第三，不同结点数设置下曲线形态不变或基本一致。

第三点特别重要。RCS的结点选择会影响曲线形状，如果你只跑4个结点就报告结果，审稿人可能会要求你用3、4、5个结点各跑一遍做敏感性分析。如果4结点和5结点跑出来的曲线都呈现同一个U形，结论就非常稳；如果4结点是U形，5结点变成S形，那说明曲线形状不稳定，可能是数据噪声驱动的，要谨慎下结论。

3.4 别只报P：效应量、置信区间与临床意义

P值只能回答“有没有统计学证据”，不能回答“效应有多大”“临床上重不重要”。临床试验和流行病学评审现在非常反感“P<0.05就宣布胜利”的写法。RCS结果报告必须包含具体效应量和置信区间。

比如报告BMI=30相对于参考值BMI=22的HR和95%CI。这需要你在图形或表格中给出明确的对比点。常见做法是选取一个有临床意义的参考值，比如BMI=22，然后计算其他BMI取值下相对该参考值的HR或OR。

效应量的解读还需要结合置信区间宽度。如果尾部的置信区间宽到横跨整条曲线，例如BMI=35对应的HR为2.1，但95%CI是0.7到6.0，这说明尾部估计精度很差，不能据此说“肥胖显著增加风险”。

3.5 亚组与交互中的RCS

这个话题在群里经常被问，关键词是“亚组分析”和“交互”。如果你要做性别分层下BMI与死亡风险的RCS，正确的做法不是分别筛出男性和女性各跑一遍RCS。那样虽然操作简单，但没法直接给出“男女之间曲线是否不同”的统计检验。

更好的做法是拟合一个带交互项的完整模型。以logistic回归为例，如果分组变量是sex（0女1男），可以这样：

logistic death c.bmi1##i.sex c.bmi2##i.sex c.bmi3##i.sex age test 1.sex#c.bmi1 1.sex#c.bmi2 1.sex#c.bmi3

上面test检验的是“BMI与结局的RCS曲线形态在男女之间是否有显著差异”，也就是交互P值。如果交互P显著，再分别画男女两组的曲线；如果交互P不显著，不建议分男女各画一组不同曲线去强行解读。亚组分析最忌讳的就是“只在某个亚组里显著，另一个亚组不显著”就宣称存在亚组差异，因为差异要用交互项检验来验证，不能只看组内P值是否跨过0.05。

4. 发表级图表制作与结果汇报

4.1 绘制剂量反应曲线的两种做法

RCS分析的图形是灵魂。一张好看的RCS曲线图，能让人一眼看出暴露与结局的关系形态。Stata里最简单的做法是用margins加marginsplot。

以logistic回归为例：

logistic death bmi1 bmi2 bmi3 age i.sex margins, at(bmi=(15(1)40)) at(age=60 sex=1) predict(pr) marginsplot

这个图给出的是在不同BMI取值下，一个60岁男性个体的事件预测概率。它能在形状上反映非线性关系，但纵轴是绝对概率，不是OR或HR。如果想看相对风险，更标准的是计算相对于参考值的OR或HR。

对于Cox模型，可以先用margins计算线性预测值predict(xb)，再手动计算相对参考值的HR。最简单的示例如下：

stcox bmi1 bmi2 bmi3 age i.sex margins, at(age=60 sex=1 bmi=(15(1)40)) predict(xb) post matrix b = e(b) matrix V = e(V) * 假设选定BMI=22为参考，它对应at()列表中的第8个位置 forvalues i = 1/26 { scalar hr`i' = exp(b[1,`i'] - b[1,8]) scalar diffvar`i' = V[`i',`i'] + V[8,8] - 2*V[`i',8] scalar ll`i' = exp(b[1,`i'] - b[1,8] - 1.96*sqrt(diffvar`i')) scalar ul`i' = exp(b[1,`i'] - b[1,8] + 1.96*sqrt(diffvar`i')) }

这段代码的核心思路是：把BMI网格上每个点的线性预测值和参考点的线性预测值做差，取指数得到HR，方差用两点的方差和协方差合并计算。代码是逻辑清晰的手工计算，不依赖额外第三方命令。缺点是代码略繁琐，但胜在完全可控，你能理解每一步到底在算什么。

如果实在不想手算，Stata社区也有专门的RCS绘图命令，但核心原理都是同样的“参考点做差取指数”。

4.2 参考值选择与曲线标注

参考值的选择直接决定图表的可读性。多数人会选暴露的临床正常值、中位数或最低风险点。比如BMI选22或23，血压选120，血红蛋白选120g/L。选择时要在论文方法部分写清楚：“以BMI=22作为参照”。

图形上一般用垂直辅助线标出参考值位置，y轴对应HR=1的横线也要画上。Stata里可以用xline(22) yline(1)加上去。如果参考点不是曲线最低点，比如你选的是中位数，但曲线最低点在另一个BMI位置，那么有些读者会误以为参考点就是风险最低点。这种情况需要在图注里说明：参考点用于标准化HR显示，不代表最低风险点。

4.3 置信区间、直方图辅助线等细节

发表级RCS图通常有两个关键元素：95%CI的带状区域和暴露变量的分布直方图或地毯图。置信区间的意义不用多说，它能告诉你哪些区域估计可靠。直方图或地毯图的作用是展示变量在哪个区间有足够的样本支持，避免读者把尾部曲线当作真实效应。

Stata里可以在twoway中把直方图和RCS曲线叠加。经验做法是：先画直方图，透明度和颜色调淡一点；再画带状置信区间；最后画HR曲线和参考线。图层顺序很重要，曲线必须最显眼，直方图只是背景信息。

4.4 论文表格如何整理RCS结果

除了图形，表格也需要呈现RCS的核心结果。推荐采用这种结构：

这样既展示主分析结果，又展示敏感性分析结果。审稿人能看到不同结点设置下结论是否一致。

5. 常见问题与避坑经验

5.1 结点应该选多少？如何做敏感性分析

结点数不是越多越好。结点多，模型灵活度高，但容易跟着噪声走；结点少，容易错过关键弯曲。文献里使用最多的方案是4个结点或5个结点。4个结点可以拟合出U形、J形等常见形态，5个结点能捕捉更多细节。如果样本量很大，比如超过1万，可以考虑5个结点；如果样本量只有几百，4个结点甚至3个结点更稳。

我的日常流程是：主分析用4个结点，敏感性分析跑3个和5个结点，对比曲线形态是否一致。如果3个、4个、5个结点下曲线都呈现同一个趋势，结论基本可以放心。如果只有某个结点数下出现明显曲线，其他结点数接近直线，别报喜不报忧，直接按线性关系报告可能更诚实。

5.2 自动结点位置与实际百分位的关系

Stata的rc_spline自动放置结点时，使用的是第5、35、65、95百分位（4结点），或第5、27.5、50、72.5、95百分位（5结点）。这意味着样本的分布会直接影响结点所在的具体数值。

如果你的暴露变量分布很不均匀，比如大多数样本集中在某个窄范围内，自动结点可能挤在一起，导致曲线在某些区间几乎没有信息。这时候我建议先画一下暴露变量的直方图，如果发现极端尾部样本很少，就要警惕尾部RCS估计的可靠性。论文方法部分可以如实写“结点位于暴露变量的特定百分位”，这样读者就能判断你的结点放置是否合理。

5.3 尾部风险：外推和过拟合

RCS最容易被质疑的地方就是尾部。数据稀少导致尾部置信区间异常宽，但很多人仍然会把尾部曲线当作重点来解读。比如BMI=40以上的样本只有几十个人，曲线却显示HR飙到5.0，这基本不能信。

处理方式有三种：第一，把图形范围限制在数据支持比较充分的区间，比如BMI的P1到P99，不要从14画到50；第二，报告时明确说明尾部估计的置信区间较宽，需谨慎解读；第三，如果尾部的极端值确实有临床意义，考虑用更稳健的方法，比如增加尾部样本量或使用稳健方差估计。过度外推是审稿人最喜欢的攻击点，一定要主动回避。

5.4 审稿人经常提的几个问题与应对思路

我在实际投稿和帮人改稿过程中，经常见到审稿人对RCS提出以下问题：

为什么不用分组分析？回答思路：分组会损失连续性剂量反应信息，切点选择主观；RCS在保留连续信息的同时允许数据驱动判断曲线形态。

为什么选择4个结点？回答思路：参考已有文献建议，4个结点足以拟合常见非线性形态，并且补充敏感性分析验证结果在不同结点数下一致。

非线性P不显著，怎么解释曲线看起来有弯曲？回答思路：曲线形态的视觉判断不能替代统计检验。非线性P不显著时，曲线上的起伏可能是抽样误差造成的，应按线性关系报告，不要硬解释。

RCS的置信区间在尾部很宽，结论是否可靠？回答思路：这是RCS的固有特性，因为尾部数据稀疏，估计精度低。我们已经在图中用直方图展示数据分布，结论主要基于中部数据范围，并在讨论中限定了外推边界。

5.5 我自己的几条实操经验

最后说几条实际跑多了才总结出来的经验。

第一，建模前先把暴露变量分布和事件数摸清楚。事件数太少时，RCS自由度相对较大，容易过拟合。如果事件数不足，优先用3个结点或者直接把变量按线性处理，不要硬上复杂模型。

第二，RCS基变量生成后，不要手动修改它们的值。基变量是由结点位置决定的，改动任何一点都会改变整个样条结构，后续检验就全乱了。

第三，报告时一定要有图形。光文字描述“非线性P=0.006”是空洞的，审稿人和读者都需要看到曲线，理解弯曲的方向和位置。

第四，用rc_spline生成基变量后，如果要作图，别忘了把协变量固定在一个有意义的参照组。不同协变量水平下，绝对风险或预测概率会不同，但相对风险曲线形态通常差别不大。报告时写清楚这些协变量取值，能大大提升可复现性。

我对RCS的总体感受是：它是一把快刀，能帮你切出漂亮的剂量反应曲线，但刀好不好用，取决于你懂不懂曲线背后的统计逻辑。拿到一个显著的非线性P，先别急着兴奋，做几次敏感性分析，看看置信区间宽度，想想尾部数据能不能支撑结论，最后再决定怎么往论文里写。这个方法我用了几十个数据集，始终有效。