在医学研究中，我们经常构建回归模型来分析自变量和因变量之间的关系。事实上，大多数的回归模型有一个重要的假设就是自变量和因变量呈线性关联，这个条件实际很难满足。常见的解决方法是将连续变量分类，但类别数目和节点位置的选择往往带有主观性，并且分类往往会损失信息。因此，一个更好的解决方法是拟合自变量与因变量之间的非线性关系，限制性立方(Restricted cubic spline，RCS)就是分析非线性关系的最常见的方法之一。
近年来在Lancet、BMJ等杂志经常见到利用限制性立方样条来拟合非线性关系。
什么是立方样条?
回归样条(regression spline)本质上是一个分段多项式, 但它一般要求每个分段点上连续并且二阶可导，这样可以保证曲线的平滑性。而限制性立方样条是在回归样条的基础上附加要求：样条函数在自变量数据范围两端的两个区间内为线性函数。
在利用限制性立方样条绘制曲线关系时，通常需要设置样条函数节点的个数(k)和位置(ti)。绝大多数情况下, 节点的位置对限制性立方样条的拟合影响不大, 而节点的个数则决定曲线的形状, 或者说平滑程度。当节点的个数为2时, 得到的拟合曲线就是一条直线，大多数研究者推荐的节点为3-5个。
在《Regression Modeling Strategies》这本书中，Harrell建议节点数为4时，模型的拟合较好，同时可以兼顾曲线的平滑程度和避免过拟合造成的精度降低。而当样本量较大时，例如因变量为未删失的连续变量并且大于100时，5个节点是更好的选择。小样本(如n<30)可以选择3个节点。以下是Harrell推荐的节点数和相应的节点位置，大家可以参考。
案例说明(模拟数据)
目前SAS、STATA、R等软件都可以进行限制性立方样条分析。基于画图的方便，我们以R语言为例进行说明。首先参照rms包，生成一个模拟数据集，包括性别(sex)，年龄(age)以及生存时间(time)和结局变量(death)。
若想分析年龄和生存率之间关系，传统的方法可以在Cox回归中将年龄作为连续变量处理，也可以对年龄进行分组，这样的做法都无法更直观的呈现年龄与死亡风险之间的关联。以下我们用限制性立方样条来分析年龄与死亡风险之间的关系：
可以看到age整体是有意义的(包括线性或者非线性关联)，然后看P-Nonlinear =0.0168<0.05，这里我们可以说年龄与死亡风险之间存在非线性关联。
如果自变量与关注的结局变量存在非线性关系，如何在文章中对结果更详细的描述呢，建议大家可以参照上文中提到的Lancet的文章。
R语言中的RCS表示生存分析之限制性立方样条；NCS是优化小区邻区关系的工具。线性回归和二元逻辑回归的非线性关系研究本文主要介绍生存分析Cox中，采用限制性立方样（Restricted cubic spline，RCS）拟合连续性自变量和事件风险之间的关系，拐点改为截断点，截断点是默认中位数。NCS可以测量是否有邻小区缺失，也可以对已有的邻小区进行检查。如果一个小区的邻区缺失或不全，必然会引起手机信号强度不足或者掉话。模型拟合对于人口模型可以采用Logistic增长函数形式，它考虑了初期的指数增长以及总资源的限制。其函数形式如下。
首先载入car包以便读取数据，然后使用nls函数进行建模，其中theta1、theta2、theta3表示三个待估计参数，start设置了参数初始值，设定trace为真以显示迭代过程。nls函数默认采用Gauss-Newton方法寻找极值，迭代过程中第一列为RSS值，后面三列是各参数估计值。然后用summary返回回归结果。
library(car)
pop.mod1 <- nls(population ~ theta1/(1+exp(-(theta2+theta3*year))),start=list(theta1 = 400, theta2 = -49, theta3 = 0.025), data=USPop, trace=T)
summary(pop.mod)
在上面的回归过程中我们直接指定参数初始值，另一种方法是采用搜索策略，首先确定参数取值范围，然后利用nls2包的暴力方法来得到最优参数。但这种方法相当费时。
还有一种更为简便的方法就是采用内置自启动模型(self-starting Models)，此时我们只需要指定函数形式，而不需要指定参数初始值。本例的logistic函数所对应的selfstarting函数名为SSlogis
pop.mod2 <- nls(population ~ SSlogis(year,phi1,phi2,phi3),data=USPop)
二、判断拟合效果
非线性回归模型建立后需要判断拟合效果，因为有时候参数最优化过程会捕捉到局部极值点而非全局极值点。最直观的方法是在原始数据点上绘制拟合曲线。
library(ggplot2)
p <- ggplot(USPop,aes(year, population))
