【深度学习】图形模型基础(5)：线性回归模型第五部分：多变量线性回归模型

1.引言

当我们从基础的线性模型 $\text{error}$ 转向更复杂的模型 $\beta_0 + \beta_1 x_1 + \beta_2 x_2 + \ldots + \text{error}$ 时，我们面临了诸多挑战。这些挑战包括但不限于：

预测变量的选择：确定哪些变量 $x$ 应该被纳入模型中，这需要对数据进行深入分析和理解。
系数的解释：理解每个系数 $\beta_i$ 的意义，以及它们如何共同影响模型的输出 $y$ 。
系数的相互作用：分析不同系数之间的相互作用，以及它们如何共同作用于模型的预测能力。
新预测变量的构建：从现有的变量中创造新的预测变量，以捕捉数据中的离散性和非线性特征。

为了应对这些挑战，我们必须学会在引入新的预测变量时如何构建和理解模型。我们将通过一系列实例来探讨这些问题，这些实例将通过R代码、数据图表以及拟合模型的图形来展示。这些实例不仅有助于我们理解模型构建的过程，也有助于我们更好地把握模型的预测能力和局限性。通过这种方式，我们可以更深入地理解数据，更准确地预测结果，并更有效地应对模型中的复杂性。

2.在模型中增加变量

在模型中添加预测变量是一个复杂的过程，因为当多个预测变量存在时，回归系数的解释变得更加复杂。每个系数的解释部分依赖于模型中的其他变量。例如，系数 $\beta_k$ 表示在所有其他预测变量保持不变的情况下，两个个体在预测变量 $x_k$ 上相差一个单位时，结果变量 $y$ 的平均或预期差异。我们通过一系列例子来说明这一点，这些例子包括单一预测变量和多个预测变量的组合。

2.1.双变量预测模型

首先，我们考虑一个二元预测变量，即母亲是否高中毕业的指标，用1表示高中毕业，0表示没有。我们使用R语言来拟合这个模型，代码如下：

model_1 <- lm(kid_score ~ mom_hs, data = kidiq)
summary(model_1)

这个模型的公式可以表示为：

$kid_score = 78 + 12 × mom_hs + error (5.1) \text{kid\_score} = 78 + 12 \times \text{mom\_hs} + \text{error}\tag{5.1}$

其中，截距78表示母亲没有完成高中教育的儿童的平均测试分数。如果母亲完成了高中教育，儿童的平均测试分数预计会增加12分。
在这里插入图片描述
图5.1 展示了儿童的测试分数与母亲是否完成高中教育的指标之间的关系。图中叠加了回归线，这条线穿过由母亲教育水平定义的每个子群体的平均值。表示高中完成情况的指示变量已经进行了抖动处理；也就是说，每个 $x$ 值上都加上了一个随机数，以防止这些点相互重叠。

2.2.连续单变量模型

接下来，我们尝试使用一个连续预测变量，即母亲的智商测试分数。拟合的模型如下：

model_2 <- lm(kid_score ~ mom_iq, data = kidiq)
summary(model_2)

模型公式为：

$kid_score = 26 + 0.6 × mom_iq + error (5.2) \text{kid\_score} = 26 + 0.6 \times \text{mom\_iq} + \text{error}\tag{5.2}$

这意味着，如果两个儿童的母亲智商相差1分，儿童的测试分数平均相差0.6分。如果智商相差10分，儿童的测试分数平均相差6分。

2.3.双变量模型的拟合

最后，我们考虑同时包括两个预测变量的线性回归模型。在R中，我们这样拟合模型：

model_3 <- lm(kid_score ~ mom_hs + mom_iq, data = kidiq)
summary(model_3)

Median MAD_SD
(Intercept) 25.7 5.9
mom_hs 6.0 2.4
mom_iq 0.6 0.1
Auxiliary parameter(s):
Median MAD_SD
sigma 18.2 0.6

2.4.双变量线性回归模型的理解

这个模型的公式为：

$kid_score = 25.7 + 6 × mom_hs + 0.6 × mom_iq + error (5.3) \text{kid\_score} = 25.7 + 6 \times \text{mom\_hs} + 0.6 \times \text{mom\_iq} + \text{error}\tag{5.3}$

在这里插入图片描述
图5.2 展示了儿童的测试分数与母亲的智商之间的关系，并叠加了回归线。线上的每个点可以被理解为，对于具有相应智商的母亲的孩子，预测的儿童测试分数，或者作为具有该智商的母亲的子群体儿童的平均分数。

在这个模型中，儿童的测试分数对母亲智商的回归斜率对于每个母亲的教育水平子群体是相同的。截距25.7在实际中没有太大意义，因为它假设母亲的智商为0，这在现实中是不存在的。在模型(5.3)中， $mom_hs \text{mom\_hs}$ 的系数6表示，在智商相同的情况下，母亲是否完成高中教育对儿童测试分数的预期影响是6分。而 $mom_iq \text{mom\_iq}$ 的系数0.6表示，在母亲是否完成高中教育相同的情况下，母亲的智商每增加1分，儿童的测试分数预计增加0.6分。
在这里插入图片描述
图5.3 展示了儿童的测试分数与母亲的智商之间的关系。浅点代表那些母亲完成了高中教育的儿童，而深点代表那些母亲没有完成高中教育的儿童。图中叠加了基于儿童测试分数对母亲智商和母亲是否完成高中教育指标进行回归的线（深色线代表母亲没有完成高中教育的儿童，浅色线代表母亲完成了高中教育的儿童）。

3.多元线性回归系数的含义

在多元线性回归分析中，解释回归系数的过程并不总是直观的，因为我们通常假设在比较不同个体时，可以改变一个预测变量而保持其他变量不变。然而，在现实中，某些情况下这种假设并不成立。例如，如果模型中同时包含智商（IQ）和智商平方（IQ^{2）作为预测变量，那么在保持IQ}2不变的情况下单独改变IQ是没有意义的。同样，如果模型中包括母亲是否完成高中教育（mom_hs）、母亲的智商（mom_IQ）以及它们的交互项（mom_hs:mom_IQ），单独考虑任何一个变量而保持其他两个变量不变也是没有意义的。

我们对回归系数的解释通常有两种视角：预测性解释和反事实解释。

预测性解释：这种解释关注的是，当比较两组在某一预测变量上相差1个单位但在所有其他预测变量上相同的个体时，结果变量平均上的差异。在线性模型中，这个系数代表了这两个项目之间预期的y值差异。这是我们之前描述过的解释方式。
反事实解释：这种解释关注的是个体内部的变化，而不是个体之间的比较。在这里，系数表示的是在模型中保持所有其他预测变量不变，向相关预测变量增加1个单位所预期的y值变化。例如，如果我们说“将母亲的智商从100提高到101，预计会导致孩子的测试分数增加0.6分”，这就是一种反事实的解释。

尽管有时统计和回归的入门教材会警告不要使用反事实解释，但它们仍然会使用类似的表述，例如“母亲的智商每增加10分，孩子的测试分数就与增加6分有关”。然而，仅凭数据本身，回归分析只能告诉我们不同单位之间的比较，并不能告诉我们单位内部的变化。

因此，对回归系数最准确的解释应该是基于比较的，例如：“当比较两个母亲教育水平相同的孩子时，如果一个母亲的智商比另一个高出x分，那么预计这个孩子的测试分数平均会高出6x分。”或者，“如果两个项目i和j在预测变量k上相差x个单位，但在所有其他预测变量上相同，那么预测的y值差异 $y_i - y_j$ 平均是 $\beta_k \times x$ 。”

这种表述可能显得笨拙，但它有助于解释为什么人们通常更喜欢使用简单的表述，如“ $x_k$ 的变化1个单位导致或与y的变化β个单位有关”，尽管这些表述可能极具误导性。我们必须接受，回归分析虽然是一个强大的数据分析工具，但其解释有时可能相当复杂。我们将在后续的文章中进一步探讨在何种条件下可以对回归进行因果解释。

4.变量的交互作用

在模型(5.3)中，我们假定儿童测试分数对母亲智商的回归斜率在母亲是否完成高中教育的各个子群体中是相同的。然而，图5.3的数据观察表明，这些斜率实际上有显著的差异。为了解决这个问题，我们可以引入mom_hs和mom_iq之间的交互项，即这两个变量乘积的新预测变量，从而允许斜率在不同子群体中变化。在R中，我们可以通过以下代码来拟合并展示包含交互项的模型：

fit_4 <- lm(kid_score ~ mom_hs * mom_iq, data = kidiq)
summary(fit_4)

拟合后的模型可以表示为：

$kid_score = − 11 + 51 × mom_hs + 1.1 × mom_iq − 0.5 × mom_hs × mom_iq + error \text{kid\_score} = -11 + 51 \times \text{mom\_hs} + 1.1 \times \text{mom\_iq} - 0.5 \times \text{mom\_hs} \times \text{mom\_iq} + \text{error}$

这个模型的图形表示分别在图10.4a和图10.4b中展示。图10.4a为不同母亲教育水平子群体的回归线，而图10.4b则将x轴扩展到零，以展示截距，即线在y轴上穿过零点的位置。这突出了一个重要问题：当预测变量的可能值范围远离零时，截距没有直接的解释意义。
在这里插入图片描述

图5.4 (a) 展示了儿童测试分数对母亲智商的回归线，其中用不同的符号表示母亲完成了高中教育的儿童（浅圈）和母亲未完成高中教育的儿童（深点）。交互作用允许每个组别有不同的斜率，浅线和深线分别对应浅点和深点。 (b) 是相同的图表，但是水平和垂直轴都扩展到零，以显示截距。

在解释这个模型的系数时，我们需要谨慎。我们通过检视特定子群体内和跨特定子群体的平均或预测测试分数来从拟合模型中获得含义。一些系数仅对特定子群体可解释：

截距项代表母亲没有完成高中教育且智商为0的儿童的预测测试分数，但这是一个不具现实意义的假设场景。如果我们在将变量纳入回归预测之前先对输入变量进行中心化处理，截距将更具解释性。
mom_hs的系数可以被理解为，在智商为0的情况下，未完成高中教育的母亲的儿童与完成高中教育的母亲的儿童预测测试分数之间的差异。这个系数的解释需要通过代入适当的数值并比较方程来实现。由于智商为0的母亲是不可想象的，这个系数的解释并不直观。
mom_iq的系数可以被看作是未完成高中教育的母亲的儿童在智商上每增加1点时，平均测试分数的变化。这对应于图5.4中深色线的斜率。
交互项的系数表示了mom_iq斜率的差异，即完成了与未完成高中教育的母亲的儿童之间的比较，也就是图5.4中浅色线和深色线斜率之间的差异。

我们可以通过为完成和未完成高中教育的母亲的儿童分别查看回归线来理解模型：

对于未完成高中教育的母亲（ $mom_hs = 0 } \text{mom\_hs} = 0 \}$ ）：
$kid_score = − 11 + 1.1 × mom_iq \text{kid\_score} = -11 + 1.1 \times \text{mom\_iq}$

对于完成高中教育的母亲（ $mom_hs = 1 } \text{mom\_hs} = 1 \}$ ）：
$kid_score = 40 + 0.6 × mom_iq \text{kid\_score} = 40 + 0.6 \times \text{mom\_iq}$

这里，对于未完成高中教育的母亲的儿童，估计的斜率是1.1；对于完成高中教育的母亲的儿童，斜率是0.6。这些斜率是直接可解释的。然而，截距仍然面临一个问题，即仅在母亲的智商为0时才具有可解释性。

我们应该何时寻找交互作用？交互作用可能非常重要，我们通常首先在那些在没有交互作用时具有大系数的预测变量中寻找它们。例如，吸烟与癌症有很强的关联。在其他致癌物的流行病学研究中，调整吸烟作为一个未交互的预测变量和作为一个交互项是至关重要的，因为其他风险因素与癌症之间的关联强度可能取决于个体是否吸烟。我们通过家庭氡暴露的例子在图1.7中说明了这种交互作用：高水平的氡与更大的癌症可能性相关，但这种差异对于吸烟者来说比非吸烟者要大得多。包括交互作用是一种方式，允许模型以不同的方式适应数据的不同子集。

在存在交互作用的情况下解释回归系数时，如果通过将每个输入变量中心化到其均值或其他方便的参考点来预处理数据，通常可以更容易地解释模型。我们将在线性变换的背景下讨论这一点。

5.交互变量

在先前的章节中，我们探讨了利用指示变量在回归分析中进行比较的方法。现在，我们将此概念进一步应用于基于调查数据的模型拟合中，以预测身高并考虑体重和其他变量的影响。以下是来自数据集earnings的样本数据，其中包括身高、体重、性别、收入、种族、教育水平、步行量、运动量、是否吸烟、紧张程度、愤怒程度和年龄等变量。

为了从身高预测体重，我们首先建立了一个简单的线性回归模型：

fit_1 <- lm(weight ~ height, data = earnings)
summary(fit_1)

模型估计结果如下：

截距（Intercept）：(-172.9)
身高（height）的系数：(4.9)

辅助参数，即预测误差的标准差（sigma）为：(29.1)

得到的回归方程可以表示为：
[ \text{weight} = (-172.9) + 4.9 \times \text{height} ]

这个方程告诉我们，身高每增加一英寸，体重平均增加4.9磅。然而，对于身高为0英寸的人，模型预测的体重是-172.9磅，这在实际情况下是没有意义的。鉴于美国成年人的平均身高约为66英寸，我们可以重新表述模型对于平均身高人群的预测：

coefs_1 <- coef(fit_1)
predicted_1 <- coefs_1["(Intercept)"] + coefs_1["height"] * 66

或者，使用posterior_predict函数来预测特定身高的体重，并得到预测值的分布情况：

new_data <- data.frame(height = 66)
pred <- posterior_predict(fit_1, newdata = new_data)
cat("Predicted weight for a 66-inch-tall person is", round(mean(pred)),
    "pounds with a sd of", round(sd(pred)), "\n")

无论是通过直接使用模型系数进行计算，还是通过posterior_predict函数获取预测的分布，对于66英寸高的个体，我们都可以预测其体重，并且得到一个更加符合实际的预测结果。

5.1.中心化预测变量

为了增强模型预测结果的可解释性，我们对身高这一预测变量执行了中心化处理。中心化是一种常用的数据预处理方法，通过从原始数据中减去一个基准值（通常是均值）来实现。在我们的数据集earnings中，我们首先计算出所有个体的平均身高，然后从每个个体的身高中减去这个均值，得到中心化身高c_height：

earnings$c_height <- earnings$height - 66  # 假设66是身高的均值

接下来，我们使用这个中心化后的身高变量来拟合一个新的回归模型fit_2，以预测体重weight：

fit_2 <- stan_glm(weight ~ c_height, data = earnings)

这个模型的估计结果提供了以下信息：

截距（Intercept）：153.4，这表示在中心化身高为0时，模型预测的体重。
中心化身高（c_height）的系数：4.9，这表示中心化身高每增加一个单位，模型预测的体重将平均增加4.9磅。

此外，模型还提供了一个辅助参数，即预测误差的标准差sigma：

sigma：29.1，这表示预测体重值围绕平均预测值的离散程度。

通过这种方式，中心化身高不仅使我们能够更直观地理解身高对体重影响的大小，还有助于我们更准确地进行模型预测。例如，如果我们知道某个人的身高是72英寸，而平均身高是66英寸，那么他们的中心化身高就是6英寸。根据模型，我们可以预测这个人的体重会比平均体重高出4.9 * 6磅。

5.2.二元回归分析

在构建更全面的回归模型时，我们可能会考虑将二元变量纳入分析，例如性别。在下面的文本中，我们将演示如何在回归模型中包含性别这一指示变量，并解释其对预测结果的影响。

我们首先通过引入性别变量male来扩展我们的模型，其中男性用1表示，女性用0表示。以下是在R语言中实现的代码：

fit_3 <- stan_glm(weight ~ c_height + male, data = earnings)

模型估计的结果如下：

截距：149.6
c_height（以英寸为单位的中心化身高）的系数：3.9
male（性别）的系数：12.0

辅助参数包括估计的标准差sigma，其值为28.8。

性别系数12.0的值表明，在控制了身高变量后，男性相较于女性，其体重平均高出12磅。这一发现提供了有关性别对体重影响的重要信息。

为了具体计算一个身高为70英寸的女性的预测体重，我们可以使用模型系数直接进行计算：

coefs_3 <- coef(fit_3)
predicted_female <- coefs_3["(Intercept)"] + coefs_3["c_height"] * 4.0 + coefs_3["male"] * 0

或者，使用posterior_predict函数来获取基于模型的预测分布，并计算平均值：

new_female <- data.frame(c_height = 4.0, male = 0)
pred_female <- posterior_predict(fit_3, newdata = new_female)
cat("Predicted weight for a 70-inch-tall woman is", round(mean(pred_female)),
    "pounds with a sd of", round(sd(pred_female)), "\n")

无论使用哪种方法，我们得到的女性预测体重都是165磅。对于相同身高的男性，相应的预测体重计算方法类似，但将性别变量male的值设为1：

new_male <- data.frame(c_height = 4.0, male = 1)
pred_male <- posterior_predict(fit_3, newdata = new_male)
cat("Predicted weight for a 70-inch-tall man is", round(mean(pred_male)),
    "pounds with a sd of", round(sd(pred_male)), "\n")

最终，我们发现对于70英寸高的男性，预测体重为177磅，比女性的预测体重高出12磅，这与性别系数的估计值一致。

5.3.分类变量的处理方法

在统计建模中，我们经常需要处理分类变量，尤其是当这些变量具有多个水平时。以种族为例，我们的数据集earnings中的ethnicity变量就包含四个不同的水平：黑人（Black）、西班牙裔（Hispanic）、其他（Other）和白人（White）。通过在R控制台中输入table(earnings$ethnicity)，我们可以得到每个种族的观测数量：

Black Hispanic Other White
177    103     37  1473

为了在回归模型中包含种族这一分类变量，我们可以将其作为因子（factor）处理。以下是如何将种族变量纳入回归模型的代码：

fit_4 <- stan_glm(weight ~ c_height + male + factor(ethnicity), data = earnings)
print(fit_4)

模型输出结果包括了每个变量的估计系数和标准差：

中位数 MAD_SD
(Intercept) 154.1 2.2
c_height 3.8 0.3
male 12.2 2.0
factor(ethnicity)Hispanic -5.9 3.6
factor(ethnicity)Other -12.6 5.2
factor(ethnicity)White -5.0 2.3

辅助参数：

中位数 MAD_SD
sigma 28.7 0.5

在输出结果中，我们注意到虽然种族变量有四个水平，但只有三个系数被显示出来：西班牙裔、其他和白人。黑人种族没有系数，因为在R中，默认将第一个水平作为基线类别，本例中即为黑人种族。

这些系数的解释如下：

截距（154.1）是当c_height（中心化处理后的身高）为0、male（性别，1表示男性）为0，且种族为基线类别黑人时，预测的体重值。
c_height的系数（3.8）表示身高每增加一个单位，体重预期增加3.8磅。
male的系数（12.2）表示在相同身高的情况下，男性的体重比女性重12.2磅。
factor(ethnicity)Hispanic的系数（-5.9）意味着与黑人相比，西班牙裔在相同身高和性别条件下，体重平均轻5.9磅。
factor(ethnicity)Other的系数（-12.6）表示与黑人相比，其他种族在相同条件下体重平均轻12.6磅。
factor(ethnicity)White的系数（-5.0）意味着与黑人相比，白人在相同条件下体重平均轻5磅。

通过这种方式，我们可以量化不同种族间体重的差异，并将这些差异纳入回归模型中进行分析。

5.4.重新设置基线因子

在进行回归分析时，因子变量的基线水平（也就是参照组）可以任意设定。在R语言中，默认情况下，因子水平会按照字母顺序进行排序。例如，如果我们有一个表示种族的因子变量ethnicity，其水平包括"Black"（黑人）、“White”（白人）等，R会将"Black"作为排序上的第一个类别，并将其作为基线水平。

如果我们想要将"White"或其他任何类别设置为基线，我们可以通过在创建因子变量时指定水平的顺序来实现。

# 假设我们有一个名为 ethnicity 的向量，包含多个种族名称
# 我们想要将 "White" 设置为基线水平

# 首先，定义种族的顺序，将 "White" 放在第一位
ethnicity_levels <- c("White", "Black", "Hispanic", "Other")

# 然后，根据这个顺序创建因子变量
earnings$ethnicity_factor <- factor(earnings$ethnicity, levels = ethnicity_levels)

# 现在，我们可以在回归模型中使用这个新的因子变量
# 例如：
fit_ethnicity <- lm(weight ~ height + male + ethnicity_factor, data = earnings)

这种方法允许我们更直观地比较不同种族与白人之间的体重差异。

这个模型使用“White”（白人）作为基线类别，因为在设置因子变量eth时，我们首先列出了它。详细审视系数如下：

在之前的模型拟合中，截距154.1代表的是c_height（中心化身高）为0、male（性别，男性）为0，且种族为“Black”（黑人）时预测的体重。在新版本中，截距变为149.1，这是在c_height为0、male为0，且种族为“White”（白人）时的预测值，即新的基线类别。变化的5.0对应于原始回归中“White”系数的负值。
身高和男性的系数没有变化。
“Black”的系数，之前由于是基线类别而隐含地为0，现在变为5.0，这是相对于新的基线类别“White”的差异。
“Hispanic”（西班牙裔）的系数从-5.9增加到-0.9，变化了5.0，对应于基线从“Black”变为“White”。
“Other”（其他）的系数同样增加了5.0。
“White”的系数从-5.0增加到了隐含的0值。

另一种方法是直接为这四个种族群体创建指示变量：

earnings$eth_White <- ifelse(earnings$ethnicity == "White", 1, 0)
earnings$eth_Black <- ifelse(earnings$ethnicity == "Black", 1, 0)
earnings$eth_Hispanic <- ifelse(earnings$ethnicity == "Hispanic", 1, 0)
earnings$eth_Other <- ifelse(earnings$ethnicity == "Other", 1, 0)

这样命名新变量不是必须的，但有时这样的命名可以使追踪变得更加容易。无论如何，一旦我们创建了这些数值变量，就可以以通常的方式将它们包含在模型中，例如：

fit_6 <- stan_glm(weight ~ height + male + eth_Black + eth_Hispanic + eth_Other, data = earnings)

这种方法允许我们直接在回归模型中包含各个种族的指示变量，从而更清晰地评估不同种族相对于基线种族（本例中为“White”）的影响。

5.5.使用索引变量访问组级预测变量

当我们在个体层次上构建回归模型时，有时需要利用来自更宏观层面的预测变量。例如，我们可能需要根据学生的背景信息以及他们所在学校家长的平均收入水平来预测学生的考试成绩。假设我们收集了20所学校中1000名学生的数据。在这种情况下，我们有一个名为 school 的索引变量，它是一个长度为1000的向量，包含从1到20的整数值，用于标识每名学生所属的学校。此外，我们还有一个名为 income 的向量，长度为20，表示每所学校家长的平均收入水平。

为了在回归模型中包含这种组级（学校层次）的预测变量，我们可以创建一个新的学生层面的预测变量。具体做法是，根据每名学生所属学校的索引，从 income 向量中取出相应的家长平均收入值，形成一个新的长度为1000的向量，我们称之为 school_income。这个向量将作为回归模型中的一个预测变量，用于评估家长平均收入水平对学生考试成绩的影响。

在R语言中，这一过程可以通过以下代码实现：

# 假设 earnings 数据框中已经包含了 student 和 income 两列
# student 是一个长度为1000的向量，包含学生所属学校的索引（1到20）
# income 是一个长度为20的向量，包含每所学校家长的平均收入水平

# 创建一个新的预测变量 school_income，长度为1000
school_income <- income[earnings$student]

# 接下来，我们可以将 school_income 包含在回归模型中
# 例如，使用 stan_glm 函数从 rstan 包拟合一个模型
# 这里假设我们还考虑了其他学生个体层面的预测变量，如 pretest（预测试成绩）、age（年龄）和 male（性别）
fit_students <- stan_glm(score ~ pretest + age + male + school_income, data = earnings)

# 输出模型摘要
print(fit_students)

在上述代码中，stan_glm 函数用于拟合回归模型，score 是因变量（学生考试成绩），pretest、age、male 和新创建的 school_income 是自变量。通过这种方式，我们能够在个体学生层面上评估学校层面变量（如家长平均收入）对学生成绩的潜在影响。

6.配对和分块设计的回归表达

我们多次讨论了如何将回归系数解释为比较。相反，通常将比较表述为回归分析也是很有帮助的。

6.1.完全随机实验

在一项完全随机化的实验设计中，我们设想一个场景：共有 $n$ 名参与者，他们被随机且平均地分配到处理组和对照组，每组各有 $\frac{n}{2}$ 人。在这种设计下，我们对处理效果的直接评估可以通过计算处理组的平均值 $\bar{y}_T$ 与对照组的平均值 $\bar{y}_C$ 之差来实现，即 $\bar{y}_T - \bar{y}_C$ 。此估计量的标准误差可以通过合并两组的标准差来计算，公式为 $\sqrt{\frac{sd^2_T}{\frac{n}{2}} + \frac{sd^2_C}{\frac{n}{2}}}$ 。

在第7.3节中讨论过，我们可以将这种推断置于回归分析的框架中，通过引入一个指示变量来区分两组：处理组和对照组。在回归模型中，这个指示变量的系数估计值本质上就是我们所关心的处理效果差值 $\bar{y}_T - \bar{y}_C$ 。此外，回归模型提供的标准误差与直接使用均值差异公式得到的标准误差相近，尽管在从不混合（unpooled）到混合（pooled）方差估计时可能会有细微的差别。

对于没有预处理预测变量的情形，回归分析和差分分析的推断结果是基本一致的。然而，回归分析的优势在于其能够更广泛地适用于包含预处理预测变量的更复杂实验设计。通过这种方式，我们可以更灵活地处理实验数据，并为实验设计提供更深入的统计洞察。

6.2.配对实验

在配对设计实验中，我们首先将参与者分成n对，每对中的两位参与者随机地被分配到实验组和对照组。在这种设计下，传统上推荐的数据分析方法是对每对参与者的测量结果进行差分，差分值记为 $z_i$ ，其中 $i$ 从1到 $\frac{n}{2}$ 。然后，我们基于这些差分值来估计治疗效果和其标准误差，具体为差分值的平均值 $\bar{z}$ 和标准差 $\frac{sd(z)}{\sqrt{\frac{n}{2}}}$ 。

除了传统的差分方法外，我们还可以通过回归分析来处理配对设计的数据。具体来说，我们可以构建一个包含所有n个数据点的回归模型，该模型中包含一个表示是否接受处理的指示变量，以及用于标识每对参与者的指示变量。首先，我们创建一个索引变量pairs，其值从1到 $\frac{n}{2}$ 。然后，我们使用以下R代码来拟合回归模型：

fit <- stan_glm(y ~ treatment + factor(pairs), data = expt)

这里，factor(pairs)用于生成对应于每对的指示变量，以便在回归模型中使用。R会自动处理模型中的共线性问题，它通过选择一对作为基线，并为其他各对包含指示变量，从而避免了完全多重共线性。这样，回归模型将包含 $\frac{n}{2} + 1$ 个预测变量：一个常数项、一个处理指示变量，以及 $\frac{n}{2} - 1$ 个剩余的配对组指示变量。

在回归模型中，处理效应可以通过处理指示变量的系数来估计，其估计值将与直接使用配对差分的平均值 $\bar{z}$ 得到的估计大致相同，其标准误差为 $\frac{sd(z)}{\sqrt{\frac{n}{2}}}$ 。此外，回归分析框架的优势在于，它很容易扩展到包含预处理信息的情况，这些信息可以作为额外的预测变量加入到模型中，从而提供了一种灵活且强大的数据分析方法。

6.3.分块设计实验

在分块设计的实验中，我们通常将参与者分配到不同的组别中，以控制潜在的变异源，提高实验的统计效力。假设有 ( n ) 名参与者被分为 ( J ) 个组，这些组可以是教室中的学生群体。在这种设计下，每个组内成员会被随机分配到实验组或对照组。

为了分析这种类型的数据，我们可以构建一个回归模型，该模型将结果变量作为响应，将处理指示变量以及 ( J - 1 ) 个组指示变量作为预测因子。在估计处理效应时，选择哪一个组作为参照组并不重要，因为模型会考虑到所有组与处理效应的关系。

以下是使用回归分析进行分块设计估计的一般步骤：

创建组指示变量：为每个组创建一个指示变量，除了一个组作为基线外，每个组都有一个对应的指示变量。
拟合回归模型：

fit <- lm(outcome ~ treatment + factor(group), data = experiment_data)

这里，outcome 是我们想要预测的结果变量，treatment 是一个二元变量，表示个体是否接受了实验处理，group 是一个因子，表示个体所属的组别。

估计处理效应：在回归模型中，处理指示变量的系数将告诉我们，平均而言，接受处理的个体与未接受处理的个体在结果变量上的差异。
考虑其他预处理变量：如果存在其他可能影响结果的预处理变量，如个体的基线测量或已知的混杂因素，这些变量可以作为额外的预测因子加入回归模型中。

通过这种方式，分块设计允许我们更精确地估计处理效应，并且通过包含组指示变量和其他相关变量，可以控制实验误差和其他潜在的混杂因素。这种设计的优势在于其灵活性和对复杂实验条件的适应性。

7.示例分析

我们通过一个具体实例来阐释模拟预测的方法，该实例涉及对美国国会选举的预测模型。本节内容将引导我们完成以下步骤：

构建预测模型：首先，我们将创建一个基于1986年选举数据的预测模型，目的是预测1988年的选举结果。
模型应用：随后，我们将使用这个模型，基于1988年的选举数据，来预测1990年的选举结果。
结果验证：最终，我们将对这些预测结果与1990年实际的选举结果进行比对，以评估模型的准确性和预测中的不确定性。

这个过程不仅展示了如何使用历史数据来预测未来的事件，而且还揭示了在预测政治事件时可能遇到的复杂性和不确定性。通过对模型的迭代应用和结果的比较，我们可以更深入地理解预测模型在实际应用中的有效性和局限性。

7.1.背景

在美国，国会选举的预测模型构建是一个复杂的过程，需要考虑多种因素。全国共有435个选区，每个选区的选举结果可以用 $y_i$ 表示，这里的 $i$ 从1到435，代表1988年在第 $i$ 个选区民主党相对于两党总票数的比例，即不包括除民主党和共和党之外其他党派的选票。图5.5展示了这些选举数据 $y$ 的分布情况。

要理解数据的变异性，首先需要考虑的是两党是否都参与了选举。直方图的两端尖峰表明许多选区的选举是没有对手的。接下来，1986年的选举结果，即上一次选举的数据，对于预测当前选举结果非常有帮助。此外，还需考虑现任议员是否竞选连任这一信息。

我们的回归模型包括以下几个预测变量：

常数项，表示模型的基线预测值。
1986年选举中，第 $i$ 选区民主党在两党投票中的占比，作为历史数据的参考。
现任议员的指示变量，如果1988年的选举中，第 $i$ 选区的民主党议员竞选连任，则该指示变量为+1；如果共和党议员竞选连任，则为-1；如果选举是开放的，即两位主要候选人当时都没有占据席位，则该指示变量为0。

由于现任议员变量是分类变量，我们可以通过不同的符号在散点图中区分共和党现任议员、民主党现任议员和开放席位的数据，如图5.6a所示。

我们将拟合一个线性回归模型来分析这些数据。尽管每位候选人的得票数是离散的，但由于每个选区的票数足够多，使用连续模型拟合不会带来本质上的损失。这个模型可以帮助我们理解不同因素如何影响选举结果，并预测未来的选举趋势。

在这里插入图片描述
图5.5 1988年国会选举数据的直方图。左边和右边的尖峰分别代表无人竞争的共和党人和民主党人。

图5.6 (a) 展示了1986年和1988年的国会选举数据。十字符号代表1988年有共和党现任议员竞选的选举，点代表有民主党现任议员，而开放的圆圈代表开放席位。回归模型中的“现任议员”预测变量对于圆圈等于0，对于点等于+1，对于十字符号等于-1。无竞争的选举结果（在0和1处）已经轻微地加入了抖动处理。
(b) 展示了回归分析的数据，去除了1988年无竞争的选举，并将1986年无竞争的选举值替换为0.25和0.75。在两个图表中都包括了y=x线作为比较。

7.2.数据处理

在预测各个选区的选举结果时，我们依赖于前一次选举的数据和现任议员是否竞选连任的信息。初选一般安排在9月，因此我们可以预期，在11月的普选前大约两个月，就能得知这些相关信息。

在分析过程中，我们遇到了一些数据问题。在某些年份中，很多选举由于缺乏竞争对手而变得没有竞争性，导致得票比例达到0或1。虽然我们可以将这些数据直接纳入模型，但我们选择了另一种处理方式：对于没有竞争对手的共和党，我们估算其得票比例为0.25；对于没有竞争对手的民主党，我们估算其得票比例为0.75。这些估算值旨在近似模拟在有对手竞争的情况下民主党候选人可能获得的票数比例。此外，我们也可以从有竞争的选举结果的分布中随机抽取值来进行估算，但对于本研究而言，采用简单的估算方法已经足够。调整后的数据集在图10.6b中进行了展示。

在处理这些无竞争的选举数据时，我们采取了以下步骤：

对于共和党在无竞争选举中的得票比例，我们赋予了一个估算值 $P_{\text{Republican uncontested}} = 0.25$ 。
对于民主党在无竞争选举中的得票比例，我们赋予了一个估算值 $P_{\text{Democrat uncontested}} = 0.75$ 。
这些估算值是基于以往有竞争选举中民主党候选人获得的票数比例来确定的。

通过这种方式，我们能够更准确地反映选举的实际情况，并为模型提供一个更加均衡和代表性的数据集。这种处理方法有助于避免由于极端值（0或1）而导致的潜在偏差，从而提高模型预测的准确性和可靠性。

7.3.拟合模型

我们建立了一个回归模型，以预测选票（每个选区民主党在两党投票中的份额），根据过去的投票（上一次选举中民主党的份额）和现任议员情况（刚刚定义的-1/0/1变量）。首先，我们根据1986年的数据预测1988年的情况，这需要从数据集中提取适当的变量：

data88 <- data.frame(vote=congress$v88_adj, past_vote=congress$v86_adj,
inc=congress$inc88)
fit88 <- stan_glm(vote ~ past_vote + inc, data=data88)

得到的结果是：

中位数 MAD_SD
(Intercept) 0.24 0.02
past_vote 0.52 0.03
inc 0.10 0.01
辅助参数：
中位数 MAD_SD
sigma 0.07 0.00

这个模型存在一些问题，通过仔细检查图5.6b中的前后对比图可以看出：开放席位的选举（在该图中由开放的圆圈显示）大多偏离回归线，表明可能存在交互作用——即现任议员和开放席位的不同斜率。此外，x=0.5刚下方和刚上方平均y值之间的跳跃，没有被incumbency_88预测变量完全拟合。我们忽略了无竞争选举中潜在投票份额的不确定性，这些份额被简单地估算为0.25和0.75。可以对这些数据拟合更好的模型，但这里的简单回归拟合足以演示基于模拟的预测推断的原则。

在这里插入图片描述
图5.7 国会选举预测模型的模拟结果。预测值y~i对应于1990年的选举。

7.4.非线性预测新数据

在对新数据点进行非线性预测分析时，我们采用了模拟的方法来获取关于预测结果的不确定性信息。对于1990年民主党在国会选举中获胜席位数的推断，我们关注的是模拟结果中超过50%获胜阈值的选举次数。通过模拟，我们得到了一个预测矩阵，其中每一行代表一个随机模拟的结果。

我们首先创建了一个用于模拟预测的数据框data90，它包含了1988年的投票数据和1990年的现任议员信息。然后，我们利用posterior_predict函数从已有模型fit88中生成了1990年选举结果的预测模拟。

为了计算民主党预计获胜的选举数量，我们将预测结果中超过0.5的值进行求和，得到dems_pred。这个过程可以通过直接使用rowSums函数简化，也可以通过循环遍历每一行的模拟结果来实现。

dems_pred <- rowSums(pred90 > 0.5)
# 或者使用循环
dems_pred <- rep(NA, n_sims)
for (s in 1:n_sims) {
    dems_pred[s] <- sum(pred90[s,] > 0.5)
}

图5.7的底部几行展示了模型参数和预测结果的统计摘要，包括中位数、均值和标准差。参数σ和β的均值和中位数与点估计非常接近，差异较小，这反映了模拟抽样的变异性。模拟结果还揭示了参数估计的不确定性，这种不确定性通过标准误差（即回归输出中报告的mad sd）体现出来。

此外，预测结果的不确定性略高于参数σ̂的估计值，但由于原始数据集的规模较大，且预测的x值都在原始数据的范围内，这种不确定性的增加非常微小。这符合我们对模型预测能力的预期。

最后，图5.7的右下角给出了民主党预计赢得的席位数的预测均值和标准差，分别为260.1和2.5。这种估计考虑了预测结果的非线性特性，而不能简单地从单个选区的预测中得出。模拟成为了评估这种复杂预测不确定性的唯一实用方法。

值得注意的是，1990年民主党实际赢得的席位数为267，虽然落在了预测范围内，但这种准确性在一定程度上是偶然的。因为当前模型并未考虑到全国范围内党派支持度的波动，这种波动在不同选举周期之间是常见的。换句话说，1986年数据预测1988年选举的模型截距，与1988年数据预测1990年选举的模型截距是不同的。直接将第一个模型用于第二个数据集的预测是不恰当的。

7.5.模拟推断和数学分析

在特定情境下，结合模拟推断与数学分析的方法能够提供更深入的洞察。例如，我们再次审视了在前面的章节中讨论的问题：估计某个特定选区选举结果出现平局的概率，或者票数恰好接近平局的概率。这项估计对于计算个人选票可能具有决定性作用的概率至关重要，而在不同选区或州之间比较这些概率，对于如何分配竞选资源的决策过程可能具有指导意义。

我们通过正态分布来估计给定预测下的平局选举概率。现在，我们探讨如何利用预测模拟来计算这一概率。考虑一个有 $n_i$ 名选民的选区 $i$ ，我们已经将民主党候选人的得票比例 $\tilde{y}_i$ 的预测分布近似为连续分布，这是合理的，因为 $n_i$ 的数量通常达到数万或数十万。因此，平局可以近似地表示为 $\tilde{y}_i$ 连续变量落在区间 $\left[\frac{1}{2} - \frac{1}{2n_i}, \frac{1}{2} + \frac{1}{2n_i}\right]$ 内。

为了通过模拟来估计这个概率，最直接的方法是进行大量预测模拟，并统计 $\tilde{y}_i$ 落在0.5 ± $\frac{1}{2n_i}$ 范围内的模拟次数比例。然而，对于实际的 $n_i$ 值，这个范围可能非常小，以至于需要成千上万次的模拟才能准确估计这个概率。例如，如果在1000次模拟中没有一次结果落在该区间内，这样的模拟结果并不具有实际帮助。

为了提高估计的精确度，我们可以采取一种结合模拟和分析的方法：首先进行1000次 $\tilde{y}$ 的模拟，然后对每个选区计算模拟结果落在0.49和0.51之间比例，并将其除以0.02乘以 $n_i$ （即在0.49和0.51之间可以容纳多少个宽度为 $\frac{1}{n_i}$ 的区间）。或者，我们可以计算模拟结果落在0.45和0.55之间的比例，并将其除以0.1乘以 $n_i$ 。即使在1000次模拟后，对于一些选区，估计的概率可能仍然为零，但这样的零估计将更为精确。

此外，我们可以使用预测分布遵循自由度为 $n - k$ 的t分布这一事实来估计极不可能事件的概率。首先，利用1000次模拟来计算每个 $\tilde{y}_i$ 的预测均值和标准差，然后使用t分布的尾部概率来计算 $\tilde{y}_i$ 落在0.5 ± $\frac{1}{2n_i}$ 范围内的概率。这种方法利用了模拟结果和分析方法的优势，提供了一种更为精确和实用的预测不确定性评估。

8. 回归模型中的数学符号

在阐释具体例子时，使用描述性的变量名很有帮助。然而，为了讨论更一般的理论以及数据操作，我们将采用通用的数学符号。本节介绍这些符号，并讨论模型的随机性方面。

8.1.预测因子

在统计模型的构建中，我们通常将 $X$ 矩阵中除去常数项之外的列称为预测因子。此外，当我们希望强调用于预测结果的输入信息时，也会使用“预测因子”这一术语。

以一个包含母亲教育水平与母亲智商交互作用的模型为例，我们可以表示为：

$kid_score = 58 + 16 × mom_hs + 0.5 × mom_iq − 0.2 × mom_hs × mom_iq + error \text{kid\_score} = 58 + 16 \times \text{mom\_hs} + 0.5\\ \times \text{mom\_iq} - 0.2 \times \text{mom\_hs} \times \text{mom\_iq} + \text{error}$

在这个回归方程中，有三个关键的预测因子：

母亲的高中教育水平（ $mom_{hs}$ ）
母亲的智商（ $mom_{iq}$ ）
母亲教育水平与智商的交互作用（ $mom_{hs} \times{mom_{iq}}$ ）

根据具体的语境，有时也会将方程中的常数项视为一个预测因子。这种分类方式有助于我们更好地理解模型的结构以及各个组成部分对预测结果的影响。
在这里插入图片描述
图5.8 回归模型的符号表示。该模型是根据给定的预测因子 $X$ 来拟合观察到的结果 $y$ 的。正如文中所述，然后可以将该模型应用于预测新的数据上给定的预测因子X̃所对应的未观察到的结果 $\tilde{y}$ （用小问号表示）。

8.2.向量-矩阵表示法在回归分析中的应用

在回归分析中，我们通常使用特定的数学符号来表示数据和模型。对于第(i)个观测对象的结果，我们用 $y_i$ 来表示，而确定性的预测值则表示为 $X_i \beta$ ，其中 $X_i$ 是包含第(i)个观测对象特征的向量，(\beta)是相应的系数向量。数据集中的观测对象用 $\ldots, n$ 来索引。

以最近的例子为例， $y_i$ 代表第(i)个孩子的测试分数。在这个数据集中，共有 $n = 1378$ 个数据点，每个数据点的特征向量 $X_i$ 包含 $k = 4$ 个元素。具体来说， $X_{i1}$ 是一个常数项，对所有人来说都定义为1； $X_{i2}$ 代表母亲的高中毕业状态，用0或1编码； $X_{i3}$ 是母亲的测试分数； $X_{i4}$ 是母亲测试分数与高中毕业状态的交互项。系数向量(\beta)的长度也是4。

模型预测与实际观测结果之间的偏差称为误差，用 $\epsilon_i$ 表示，并假设它们遵循均值为0，标准差为(\sigma)的正态分布，记作$ \text{normal}(0, \sigma) $。

在模型估计中，观测结果与预测值之间的差异称为残差。因此， $\beta$ 和 $\hat{\beta}$ 分别表示误差向量和残差向量。我们用$ \hat{y} $来表示给定新数据$ X’$时，模型的预测结果。

在不同的学科领域，对于我们所说的预测因子和结果变量的术语使用习惯可能有所不同。有些领域使用“自变量”来指代预测因子，“因变量”来指代结果变量。这些术语源自于回归模型在实验中应用时，输入变量的操纵可能导致预测因子相互独立的历史背景。然而，在社会科学研究中，这种情况并不常见，因此我们倾向于避免使用这些术语。有时，预测因子和结果变量也被称为“左侧变量”和“右侧变量”。

8.3.模型的数学表达

经典线性回归模型可以用以下数学形式表示：

$y_i = \beta_1 X_{i1} + \ldots + \beta_k X_{ik} + \epsilon_i, \text{for:} i = 1, \ldots, n,$

其中误差项 $\epsilon_i$ 具有独立的正态分布，均值为0，标准差为(\sigma)。

模型的等价表达式为：

$y_i = X_i \beta + \epsilon_i, \quad \text{for： } i = 1, \ldots, n,$

这里， $X$ 是一个 $\times k$ 的矩阵，其第(i)行为 $X_i$ 。使用多元表示法，可以写作：

$y_i \sim \text{normal}(X_i \beta, \sigma), \quad \text{for： } i = 1, \ldots, n.$

为了得到更简洁的表示，我们可以使用：
$\sim \text{multivariate normal}(X \beta, \sigma^2 I),$

这里， $y$ 是一个长度为 $n$ 的向量， $X$ 是一个 $\times k$ 的预测因子矩阵，(\beta)是一个长度为 $k$ 的列向量， $I$ 是一个 $\times n$ 的单位矩阵。通过最小二乘法拟合模型，我们可以得到系数的估计值$ \hat{\beta} $和标准差的估计值$ \hat{\sigma} $，正如我们对于只有一个预测因子和常数项的简单回归所演示的那样。

8.4.最小二乘估计、最大似然估计与贝叶斯推断

在线性回归分析中，无论是使用一个还是多个预测变量，估计参数和进行统计推断的步骤是一致的。我们从最小二乘估计开始，即寻找一个系数向量 $\hat{\beta}$ ，它能够使残差平方和(RSS)达到最小：
$\text{RSS} = \sum_{i=1}^{n} (y_i - X_i \hat{\beta})^2$

对于标准线性回归模型，如果预测变量的测量是准确的，并且误差项相互独立、具有相同的方差且呈正态分布，最小二乘解同时也是最大似然估计。残差标准差的估计公式：
$\hat{\sigma} = \sqrt{\frac{\text{RSS}}{n - k}} \quad (5.5)$
这里的 $k$ 代表回归系数的总数。当回归模型中只有一个预测变量时，即 $k = 2$ ，公式(5.5)就简化为单变量的形式。

在贝叶斯推断中，虽然估计参数的方法与最小二乘和最大似然估计不同，但它们都旨在寻找最佳的参数估计，以便对模型进行解释和预测。贝叶斯方法特别强调先验信息的使用，并在参数估计中考虑不确定性。然而，在标准线性回归模型中，最小二乘估计因其简单性和在正态误差假设下的最优性质而被广泛采用。

8.4假设检验

我们建议避免使用传统的零假设显著性检验，这里我们简要回顾两种常见的检验方法。

t检验的目的是证明某个回归系数是否在统计上显著地不同于零。具体来说，它是一个针对特定系数 $\beta_j$ 等于零的零假设 $H_0: \beta_j = 0$ 的检验。在满足特定假设的情况下，标准化后的回归系数 $\frac{\hat{\beta}_j}{\text{s.e.}_j}$ 在零假设成立时，将近似地遵循自由度为 $n - k$ 的t分布。因此，如果该标准化系数的绝对值 $\left| \frac{\hat{\beta}_j}{\text{s.e.}_j} \right|$ 超过了相应t分布的临界值，我们可以在预设的显著性水平下拒绝零假设。

F检验旨在证明整个回归模型是否具有预测能力，而不仅仅是模型中的某个特定系数。它是一个针对模型中除了常数项外所有系数都等于零的零假设 $H_0$ 的检验。在一定假设下，总平方和与残差平方和的比值，经过适当缩放后，会遵循F分布。因此，如果这个比值超过了依赖于样本量 $n$ 和预测变量数量 $k$ 的临界值，我们可以拒绝零假设。

我们几乎不使用回归的假设检验，因为我们很少遇到需要将系数视为精确为零的情况。因此，零假设的拒绝对我们来说并不重要，因为我们从未真正认真考虑过系数为零的情况。在现实世界中，只要数据量足够，几乎所有的假设都可以被拒绝。

尽管如此，估计中的不确定性是实际存在的。我们尊重假设检验所关注的核心问题，即评估一个估计值何时可能被噪声所淹没，以至于从数据角度来看，某个特定系数或系数集合可能实际上为零。

我们建议通过查看标准误差以及参数估计来处理这些问题。当估计值较为嘈杂时，我们推荐使用贝叶斯推断，因为先验信息的使用有助于稳定估计和预测。如果目标是检验在不包含某些变量的情况下是否能够进行有效预测，我们建议通过交叉验证来比较不同模型的性能。

9.加权回归分析

当误差项满足独立同分布且具有常数方差的正态分布时，最小二乘回归等同于最大似然估计。在计算平方残差和时，每个数据点被赋予了相同的权重。

然而，在某些特定情境下，在拟合模型的过程中，对不同的数据点赋予不同的权重是合理的。这时，我们可以采用加权最小二乘回归方法。在这种方法中，我们寻找的估计量 $\hat{\beta}_{\text{wls}}$ 能够最小化以下加权平方残差和：
$\sum_{i=1}^{n} w_i (y_i - X_i \beta)^2$
这里的 $(w_1, \ldots, w_n)$ 是一个非负权重向量。权重较高的数据点在估计中占有更大的比重，这导致回归线更贴近这些点。

与普通最小二乘估计相似，我们可以通过矩阵代数来推导出加权最小二乘估计。最终得到的估计公式如下：
$\hat{\beta} = (X^t W^{-1} X)^{-1} X^t W^{-1} y \quad (5.6)$
在这个公式中， $W$ 是权重矩阵，定义为 $\text{Diag}(w)$ ，其中 $\text{Diag}(w)$ 表示以向量 $w$ 为对角线元素的对角矩阵。通过这种方式，加权最小二乘回归能够考虑到不同数据点的重要性，从而提供更为精确的模型估计。

9.1.导致加权回归的三种模型

加权最小二乘估计可以通过以下三种不同的模型来推导：

代表更大总体的观测数据：这是实际应用中回归权重最常用的方法。我们使用加权回归来拟合样本数据，目的是为了估计如果能够对整个总体进行拟合时会得到的未加权线性模型。例如，如果我们的数据来源于一项调查，该调查过度抽样了年长的白人女性，而我们想要估计整个人口的回归模型。在这种情况下，我们会根据该样本中个体所代表的人口类型的比例，给予相应的权重。在这个例子中，男性、年轻人和少数族裔成员将会被赋予更高的权重。将这些权重纳入回归分析中，是为了在总体层面上，而非仅仅在样本层面上，尽可能减少平方误差的总和。
重复观测的情况：更具体地说，假设每个数据点可以代表一个或多个实际观测值，即索引 $i$ 代表具有相同预测向量 $x_i$ 的 $w_i$ 个数据点的集合，并且 $y_i$ 是相应 $w_i$ 个结果变量的平均值。在这种情况下，对压缩后的数据集 $(x, y, w)$ 进行加权回归分析，实际上等同于对原始数据进行未加权回归分析。
不等方差的情况：从完全不同的角度来看，当误差项具有独立正态分布但方差不相等时，加权最小二乘是该回归模型的最大似然估计，其中误差项的标准差 $\text{sd}(\epsilon_i)$ 与 $1/\sqrt{w_i}$ 成正比。这意味着，方差较大的测量值在模型拟合过程中会被赋予较低的权重。正如第11.1节进一步讨论的，不等方差通常对于估计回归系数的目标来说并不是一个大问题，但当涉及到对个别案例进行预测时，它就变得更加重要了。

这三种模型虽然都会得到相同的点估计值，但它们隐含的标准误差和预测分布是不同的。在最常见的情况下，即权重用于调整样本和总体之间的差异，一旦我们得到了这些权重，我们首先会将权重向量重新标准化，使其均值为1（在R语言中，可以通过设置 $w/\text{mean}(w)$ 来实现），然后我们可以在回归分析中将这些权重作为一个参数（例如，使用 $stan_glm ( y ∼ x , d a t a = d a t a , w e i g h t s = w ) \text{stan\_glm}(y \sim x, data=data, weights=w)$ ）。

9.2.利用权重矩阵处理误差相关性

公式(10.6)同样适用于使用权重矩阵 $W$ 的情况，它能够为具有正态分布误差的模型提供最大似然估计，其中误差的协方差矩阵为 $W^{-1}$ 。这种方法有时被称为广义最小二乘法（Generalized Least Squares, GLS）。相关性模型在多种数据分析场景中都有应用，例如时间序列分析、空间统计、群集样本以及其他存在数据结构化特征的场合。

在数据点之间存在相关性的情况下，使用权重矩阵来调整误差的相关性变得尤为重要。这种相关性可能源于多种因素，如时间上的连续性、空间上的接近性或样本的群集特征。权重矩阵$W $能够捕捉这些误差项之间的相关结构，使我们在参数估计时能够考虑到它们。

广义最小二乘法是一种灵活的分析工具，它适用于以下等多种情况：

时间序列分析：在时间序列数据中，由于时间的连续性，不同时间点的观测值之间可能存在自相关性。
空间统计：在空间数据中，由于地理位置的接近性，不同位置的观测点可能表现出空间相关性。
群集样本：在群集样本设计中，由于群集内部的共同特征，同一群集内的个体之间可能存在相关性。

通过引入权重矩阵 $W$ ，我们能够更准确地估计模型参数，并且得到更为可靠的标准误差和预测分布。这种方法使我们能够深入理解数据的相关结构，并且在分析过程中予以适当的处理。

10.将同一模型应用于多个数据集

在数据分析中，对多个数据集或现有数据集的子集重复拟合同一回归模型是一种常见做法。例如，研究者可能会利用多年或多个国家的调查数据来探究身高与收入之间的关联，或者在美国的不同地区或州内进行类似的分析。

本书未涉及的是多层次建模技术，它允许我们通过部分地汇集不同拟合的信息来重复估计回归模型。而我们这里讨论的是一种更非正式的方法，即在不同年份或群体之间不进行信息汇集，分别估计回归模型，然后将这些估计结果一起展示，这可以看作是多层次建模的非正式先驱。这种复制方法易于在多种情况下应用。

重复建模后，通过时间序列图展示估计结果的方法有时被称为“秘密武器”，因为它简单、有效，但作为数据分析工具却很少被使用。我们认为，使用频率低的一个原因可能是，一旦数据集的时间序列特性被认识到，人们自然会想要迈出下一步，直接对其进行建模。然而，在实际操作中，存在许多问题，对于这些问题，进行跨截面分析是具有启发性的，并且通过时间序列展示来呈现趋势是恰当的。
在这里插入图片描述
图5.9 展示了从1972年到2000年每次总统选举期间的民意调查数据分别拟合得到的党派认同相对于政治意识形态、种族和其他预测变量的回归估计系数 ± 0.67 标准误差（因此，是50%的置信区间）。这些图表使用了不同的尺度，并且预测变量根据它们系数的大小大致按降序排列。这组图表说明了如何从一系列回归分析中展示推断结果。

10.1.党派认同的预测分析

本研究采用一系列横截面回归分析，阐释了所谓的“秘密武器”，即基于政治意识形态和人口统计变量来建立党派认同的模型。政治学者始终对党派认同及其演变趋势抱有浓厚兴趣。本研究利用全国选举研究的数据，该研究通过1至7的量表来衡量党派认同（1代表坚定的民主党支持者，7代表坚定的共和党支持者），我们将其作为连续变量进行分析。研究中纳入的预测变量包括：

政治意识形态（1代表坚定自由派，7代表坚定保守派）
种族（0代表白人，1代表黑人，0.5代表其他种族）
年龄（分为几个年龄段：18至29岁、30至44岁、45至64岁、65岁以上，以最低年龄段为参照组）
教育水平（1代表未完成高中，2代表高中毕业，3代表部分大学教育，4代表大学毕业）
性别（0代表男性，1代表女性）
收入水平（按百分位数分组：1代表0-16百分位，2代表17-33百分位，以此类推至5代表96-100百分位）

图5.9展示了这些估计系数随时间变化的趋势。意识形态和种族是影响党派认同的两个最为关键因素，而且它们的影响系数随时间呈现出变化趋势。由于意识形态是在一个7点量表上测量的，因此在比较自由派（意识形态=2）和保守派（意识形态=6）时，需要将系数乘以4来预测期望的变化。年龄和性别的预测差异在所研究的30年间发生了显著变化。总体而言，图5.9展示了在图表中并列展示多个模型拟合结果的强大之处，这有助于揭示平均模式和趋势，特别是与每个图表上显示的零线进行比较时。

11.参考文献说明

线性回归在社会和物理科学领域的应用可以追溯到数个世纪以前；参见Stigler（1986年）。关于回归分析，有许多深入的教科书资源，包括Neter等人（1996年）的作品；Ramsey和Schafer（2001年）的作品，专注于模型理解、图形展示和实验设计等议题；Woolridge（2001年）的作品，则从计量经济学的视角阐述回归建模。Fox（2002年）在应用回归分析的教学中融入了R语言的教学。Carlin和Forbes（2004年）提供了关于线性模型和回归概念的优秀概念性引导。

在母亲和儿童测试分数的案例中，我们采用了资格测试的数据，这些数据经过标准化处理，均值设为100，标准差为15，与智商测试分数的标准化方式相同。关于儿童测试分数和母亲就业的更深入讨论，请参阅Hill等人（2005年）的研究。身高和脚长的例子引自McElreath（2020年）的著作。