【临床研究】一个你无法逃避的问题：多元回归分析中的变量筛选

临床模型研究，说到底是做一个模型，那么模型应该如何纳入自变量，纳入哪些自变量，这都是至关重要的问题。线性回归，逻辑回归和Cox比例风险回归模型是被广泛使用的多元回归分析方法。我们在前面的几篇文章中解释过他们的统计学意义、应用及结果释义。但是我们很少讨论自变量筛选的方法，这些方法在数据分析和撰写论文时应用较为混乱，却十分重要。本集整理并总结了前沿的自变量筛选方法，我们来一观究竟。

变量筛选方法及原则

Background

在自变量筛选遇到问题时，研究者常常求助统计学家，统计学家会建议使用软件中自动筛选，例如IBM SPSS中的Logistic回归和Cox回归，给出了7种变量筛选的方法：

1. 条件参数估计似然比检验（向前：条件）；
2. 最大偏似然估计的似然比检验（向前：LR）；
3. Wald卡方检验（向前：Wald）；
4. 条件参数估计似然比检验（向后：条件）；
5. 最大偏似然估计的似然比检验（向后：LR）；
6. Wald卡方检验（向后：Wald）；
7. Enter法（变量全部进入）。

实际上，在阅读临床文献时可以发现许多作者采用了这些筛选方法。其流程是：首先逐个对变量进行单因素回归分析，把单因素回归分析p值小于0.1的纳入最终的回归方程（此处可根据样本量的大小选择性调整p值的范围，如果样本量过大，可以控制p<0.05，如果样本量过小，可以控制p<0.2, 但在通常情况下，P值的范围在0.05-0.2之间，并无统一标准）。

对于研究者而言，如何选择更好的方法实际上是因人而异。但为了更好地筛选变量，我们仍应遵守一些基本原则：

1. 足够的统计学效能。这里有一个经验性判断统计学效能是否足够的标准：即一个单变量因素至少需要对应有20个有效样本量的存在，举例来说，比如我们做Cox回归分析，如果我们收集了10个与预后相关的变量，那么至少应该有200个患者出现了我们定义的终点事件，比如死亡（应包括200名死亡患者，而非总共200名患者，因为那些没有终点事件的样本我们一般不视为有效样本）。

1. 依靠临床研究报告以往经验筛选。当不满足足够的统计学效能这一原则时，应该采用大多数临床研究报告中使用的变量筛选方法，即逐个对变量进行单因素回归分析，然后把单因素回归分析p值小于一定范围的变量纳入最终的回归方程。正如之前提到的，这种方法虽被广泛应用，但颇有争议。
2. 结合临床知识筛选。在实践中，依靠临床研究报告以往经验分析的方法也有可能无法让研究者“心满意足”。有时我们会发现某些已经确定的与某种疾病临床预后相关的变量，在单因素分析的时候并未达到我们所设定的变量筛选标准，而被排除在多因素回归模型之外。例如，在某项前列腺癌预后因素分析的研究中，作者并未发现Gleason评分与预后显著相关，而临床知识指明的是Gleason评分与前列腺癌患者的预后显著相关，此时我们应该怎样做出取舍呢？答案显而易见，对于那些已知的确定与某疾病预后显著相关的变量，即便未达到我们设定的统计学筛选标准，我们也应该纳入回归模型，这么做的考量即是从临床专业角度筛选变量。

综上，我们想着重强调第三种基本原则，筛选变量时理应充分结合临床相关知识，统筹考虑统计学上的单因素分析结果与已知临床专业知识决定纳入回归方程的变量。当然的，单因素分析结果和临床因素，样本量和统计效能理应综合考虑。下文我们将对变量筛选的进行具体阐释。

争论与共识

五个要点要记牢

有关变量筛选的争论从未停止过，统计学家有统计学家的考量，但临床论文的作者往往并未严格遵照统计学家的建议。我们暂且搁置争议，先不论对错，目前确实也很难分出对错，因为临床研究往往有很多现实的困境：比如样本量不够，比如很难根据我们已知的知识来确定哪些因素才是导致结局发生的“罪魁祸首”。

尽管如此，变量筛选也并非毫无章法可言，我们回顾顶级医学杂志发表的文章，其中有关变量筛选的方法大体考虑以下5点：

1. 结合临床专业知识。这一点应该是变量筛选最基础的考量，医学统计分析如果脱离临床可能会变得虚无缥缈。根据目前临床专业知识，已知的确定与结局发生有关的变量应该纳入回归模型，而不去过多考虑其统计学参数。如前文所述例子：Gleason评分与前列腺癌预后显著相关，这是我们的共识，那么对于评价前列腺癌预后影响因素的分析，像Gleason评分这样的变量应该参与建模。而无需去考虑其统计学参数是否有统计学意义。

1. 根据单因素分析结果筛选变量。单因素分析p值“显著”的变量放入多元回归方程。所谓p值“显著”一般设为p<0.1，也可设为p<0.05或者p<0.2，需根据样本量大小做出调整，样本量够大可以把P值调小，样本量不足则需要更保守一点，把P值设大。这一类变量筛选方法在既往发表的临床研究论文中很常见，即便是高分论文中也很常见。尽管对于这种方法，绝大多数统计学家提出质疑。而现实情况是：如果弃用这种方法，目前是否有更加准确的、更科学的替代方法呢？答案显然是否定的，统计学家也并未找到更有说服力的新方法。我们习惯于轻易的否定一个既有的方案，但在没有找到新的解决方案之前，否定既有的方案是无意义的。
2. 根据混杂因素“Z”对试验因素或暴露因素“X”的影响大小筛选变量。具体说来，先观察，调整“Z”与不调整“Z”，“X”对因变量“Y”的作用是否有变化。先运行仅纳入“X”的基本模型，记录回归系数β1，再在该模型中加入“Z”，看β1变化多大，通常认为β1变化超过10%则需要调整该变量，否则不需要。这种方法与第2种根据单因素分析结果筛选变量的差别在于：这里把混杂因素对试验因素的影响量化。这种方法也并非是完美的，“Z”和“X”对“Y”的影响也同样可能受到其他混杂因素的影响。如果沿着这个思路继续思考下去，我们会陷入一个思维的“怪圈”。剩下的更复杂的方法学难题留给更有智慧的人，我们暂且把这个方法认为是一个可供参考的变量筛选方法，尤其是对于那些研究目的非常明确：明确探究“X”对因变量“Y”的作用，而这种作用很有可能是客观存在的，我们要做的无非是去调整这些混杂因素。
3. 决定最终纳入模型的变量个数。这点至关重要。如果样本量足够大，统计效能足够，我们可以借助统计软件提供的变量筛选方法自动筛选变量，并可以筛选出适合的独立影响结果的变量。但“理想很饱满，现实很骨感”。往往我们考虑变量很多的时候，样本量却很小。在统计效能和变量筛选之间我们得做出折中，折中才能带来更好的结果。
4. 其他方法。以上列举了四种变量筛选方法，此外还有很多其他变量筛选方法，比如根据模型参数：决定系数R^2，AIC，似然比对数、C-Statistics等等。筛选变量的方法越多，恰恰证实了目前变量筛选没有公认的最好的方法。我们写这篇文章的目的也不在于讨论孰优孰劣，只有基于客观条件的最合适研究人员的选择才是好的，因此，根据实际情况选择最合适的筛选变量方法才是本文的目标。

不同类型变量的纳入方法

Different values, different ways

为了能够帮助大家掌握不同类型变量的筛选方法，我们针对不同类型的变量进行了以下的总结，我们逐一看过来。

01

连续变量

对于连续变量，有一个很好的处理方法可供参考。如果变量与结果之间的关系是线性的，则可以在回归公式中包含连续变量。如果不是，可以将其转换为二分法变量或序数分类变量，然后将它们放入回归公式中。我们已经用这种方法把原来的连续变量变成了分类变量。我们进行这种转换是因为变量可能与结果不是线性的。可能会出现一些其他关系，而不是线性关系。

连续变量转换总结

回归模型中包含连续变量时，应尽可能将原始变量包含在回归模型中，并考虑实际需要。变量可以根据一些规则进行转换。为了更好的专业解释，可以进行二类分组，等分分组，等距分组和临床临界值分组。通过最优截断点分析，将连续变量转化为分类变量，并将其作为哑元变量引入回归模型。在回归模型中，连续变量可以以不同的方式表示。我们将于下文中举具体的例子。无论以何种方式呈现，总的原则是，这种变化更有利于专业人士的解读和理解。

1. 普通转换

对于服从正态分布的连续变量，这不是问题。然而，当我们面对不符合正态分布的数据时，我们可以根据某个函数进行变换，然后对这些数据进行归一化处理。并与回归模型进行拟合。原始数据可以根据其自身特点用不同的函数进行归一化，如平方根法、LNX法、Log10X法和(1/X)法等。如果对原始数据进行了归一化处理，则应对正态变换后的变量进行解释，而不是回归模型中的原始变量，也可以根据变换中使用的函数来估算原始自变量对原始因变量的影响。

例如，作者在其发表于JACC，2016的文章中做了正态性检验。原始表达式如下：用Kolmogorov-Smirnov检验评估连续变量的正态性。正态性检验方法包括使用数据分布参数（偏度值和峰度值）和使用数据分布图（直方图、P-P图、Q-Q图）。或者采用一些非参数检验方法（Shapiro-Wilk检验、Kolmogorov-Smirnov检验）来帮助我们评估数据的正态性。在他们的研究中，肌钙蛋白I、NT-proBNP或corin等变量符合非正态分布。因此，作者用中位数（四分位数-三分位数）来描述这些招募对象的基线特征。例如，用肌钙蛋白I的中位数对corin进行多元线性回归分析。原表达如下：采用多元线性回归分析确定影响corin水平的因素。肌钙蛋白I、NT-proBNP和corin的水平通过Log10转换进行了标准化。用log10函数对肌钙蛋白I、NT-proBNP、corin等变量进行了归一化，并将其纳入多元线性回归。然后进行Cox回归分析。虽然对Cox回归没有特别的要求，但使用log10函数将肌钙蛋白I、NTproBNP和Corin归一化。所有这三个变量都被纳入多元线性回归模型，以保持与原始模型的一致性。

2. 对于固定增量的每一次变化进行变换

如果连续变量以其原始形式直接引入模型，则回归参数被解释为因变量因每个单位变化而产生的变化的影响。然而，有时这种变化的影响可能是微弱的。因此，我们可以将连续的自变量按固定的区间，以等距分组的方式转化为一个分类变量，然后将它们引入到模型中进行分析。这一分组有助于更好地理解和应用于患者。例如，我们包括年龄在31到80岁之间的患者。我们可以根据10岁的年龄间隔将其分为10-40、41-50、51-60、61-70、71-80组。然后五个已经设置好的哑元变量将被包含到模型中进行分析。但是，如果变量的范围很大，按照前面提到的方法进行分组会导致分组和哑元变量太多，这在分析过程中是相当冗余的，临床上也很难解释。相反，有些数据的范围很小，不能再分组，也不能转换成分类变量。那么，面对这两种情况，我们该怎么办？

在这里，我们可以参考发表在JACC，2016的一篇文章。我们发现，在模型中，作者使用了很多“per”，如每5%变化、每0.1U每100ml/min等，这是连续变量在每次变化中以固定增量进行的变换，这种变换一直存在于“per+区间+单位”中。我们将在本文中演示两个示例。吸氧效率的平均斜率为1655U，5-95%的患者将从846 U变化到2800U，这确实是一个很大的范围。如果将原始数据放入公式中，每1U的变化将导致风险率非常微小的变化，这在临床上是没有意义的。而如果将其转换为分类变量，将出现许多组。因此，作者将每100U的变化纳入模型中，发现当氧摄取效率斜率以每100U为单位增加时，死亡风险率将降低9%（HR=0.91，95%CI:0.89–0.93）。另一个例子是可变峰值呼气交换比。中位数为1.08 U，5-95%的患者将从0.91-1.27 U变化。这是一个很小的范围。如果将原始数据放入公式中，每1U的变化将导致风险率发生很大变化。在临床实践中，1U的改变是非常罕见的，这种结果会使实用性大打折扣。由于数据的范围较小，其分类变量变换也会非常困难。因此，作者将每0.1U的变化纳入模型中，发现当峰值呼气交换比每增加0.1U时，死亡风险率将降低6%（HR=0.94，95%CI:0.86–1.04），虽然没有统计学意义。

那么，我们怎样才能实现这种转变呢？如果我们想把因子从每1个单位改为100个单位，它会是原来的100倍，我们只需要将原始变量除以100，然后将其包含到模型中即可。类似地，如果我们想将因子从1个单位改为0.1个单位，它将是原来的1/10。只需将原始变量乘以10并将其纳入回归模型中。

3. 标准差的变换

在临床研究中，我们得到了另一种转换方法：每增加一个SD时的自变量变化。让我们看看2016年在JACC上发表的一篇文章。模型中年龄和收缩压按SD增加值计算。年龄每增加一个标准差，动脉粥样硬化性心脏病（ASCVD）的风险增加70%（HR=1.70，95%CI:1.32-2.19）。收缩压（SBP）每升高一个标准差，ASCVD风险增加25%（HR=1.25，95%CI:1.05-1.49）。本文将连续变量以每SD递增的形式引入到模型中。假设变量服从正态分布，则均值±1SD区间内的面积为68.27%，平均值为±1.96，SD区间内的面积为95%。如果平均值为±2.58，则SD区间内的面积为99%。我们可以说，如果数据范围在4SD以内，大约95%的样本将被覆盖。因此，新的变量，特别是那些临床解释尚不清楚的罕见变量，我们可以把每个SD放入模型中。这可以引导患者根据自己的实际测量结果，判断自己在人群分布水平的几个标准差范围内，然后评估相应的风险会发生多大的变化。

做这种转换很简单，我们可以通过这两种方式做到这一点：

（1）在建立回归模型之前，需要对原始连续变量进行归一化处理，并将归一化后的自变量纳入回归模型。得到的回归系数是因变量对各因变量SD的影响。(注意：此处仅对自变量进行归一化)。

（2）如果原始变量没有归一化，可以直接将原始变量带入模型，得到非标准化系数，然后用自变量的标准差乘以非标准化系数来计算自变量的标准差系数，也称为标准化系数。这是因变量对自变量的每个附加SD的影响。

02

等级变量

等级变量非常常见。它是一种有序的多类别变量。通常，多个数据可以出现在同一变量中，并且这些数据彼此等级相关。如高血压分级(0=正常，1=高正常，2=1级，3=2级，4=3级)，尿蛋白水平(0=−，1=±，2=+，3=++，4=+++，5=++++)，药物疗效(无效，改善，治愈)，均为等级变量。它不同于无序多类别变量。有序多类别变量呈现单调递增或递减。当Logistic回归模型中存在有序的多类别变量时，不建议将这些变量直接作为连续变量引入，除非每单位变化会导致相同的风险比变化。然而，大多数情况下，它不会那么理想地改变。因此，我们建议将有序的多类别变量当作哑元变量来对待，这样就可以将每一级与另一级进行比较。当结果不是线性相关时，应使用最优尺度回归来探索效应拐点。

03

无序多分类变量

无序多分类变量是一种非常常见的变量类型。通常，多类别变量中有几个可能的值，而彼此之间没有层次关系。如种族(1=白色，2=黑色，3=黄色，4=其他)，给药方式(1=口服，2=皮下注射，3=静脉注射，4=其他)，都是无序的多类别变量。当无序多类变量在Logistic或Cox回归模型中时，需要设置哑元变量才能将其引入模型。下面我们将介绍哑元变量的设置方法。

哑元变量的设置方法

（1）Indicator:该方法用于指定分类变量的参考水平。这里计算的参数指的是变量的最后一级或第一级。这取决于您是选择以下参考类别中的第一个还是最后一个。

（2）Simple:该方法可以计算分类变量的各个水平与参考水平的比值。

（3）Helmert:我们将分类变量的水平与以下水平的平均值进行比较。如果某一水平的系数增加且具有统计学意义，则表明分类变量从该水平开始对风险率产生影响。它还可以用在有序分类变量中。

（4）Difference:该方法可以将分类变量与各级平均数进行比较。由于与Helmert完全相反，因此它也被称为反向赫尔默特(Reverse Helmert)。例如，级别2的平均值可以与级别1的平均值进行比较；级别3的平均值可以分别与级别1和级别2的平均值进行比较，以此类推。如果系数在一定的水平上变小，并且在统计上不显著，则分类变量对风险比的影响达到平台期。这一通常用于有序分类变量，如吸烟量等。假设研究者把它们作为独立的无序多分类变量来分析，那就没有意义了。

（5）Repeated：将分类变量的级别与它们相邻的级别进行比较（第一级别除外），其中“前一级别”作为参考级别。

小结

Take-home messages

我们已经对筛选方法和变量变换方法进行了总结。虽然事实上并没有绝对完美的方法，因此只能选择最贴切的方式。我们可以构建多个模型，并根据以前发表的，特别是影响因子较高的临床试验类文章，得出每个模型的客观结果。这是一种比较好的筛选方式。在构建预测模型的过程中，除了对所有可能的变量进行变量筛选外，还会有具体的考虑。例如，某些恶性肿瘤的TNM分期虽然对预后的价值不大，但因其易于在临床实践中应用而被广泛应用。我们如何评估模型的准确性和简洁性？虽说变量越多，模型的预测精度越高，但临床应用的难度也就会相应的增大。总之，我们在构建模型时应该选择好一个平衡点。

参考文献：Zhi-Rui Zhou, Wei-Wei Wang, Yan Li, et al. In-depth mining of clinical data: the construction of clinical prediction model with R.Annals of Translational Medicine.