现代偏好获取中的回归和贝叶斯方法
现代偏好获取的回归和贝叶斯方法
应用于简单的冰沙制作
线性回归通常被认为是预测建模的工作马，但其应用不仅限于直接的预测任务。本文旨在通过介绍Probit线性回归作为建模偏好的有效工具，丰富回归技术的讨论。此外，我们采用贝叶斯框架，从经典线性回归过渡到贝叶斯线性回归，阐明了基于成本的优化（特别是二元交叉熵损失最小化）与最大似然估计之间的内在关系。
通过这样做，我们旨在证明正则化可以被视为贝叶斯先验选择的一种形式，从而将成本函数方法与概率推理联系起来。
最后，我们将讨论贝叶斯线性回归不仅允许点估计，而且为这些预测提供分布，提供了更丰富、具有不确定性意识的观点。
贝叶斯框架
贝叶斯框架识别了两个主要组成部分：数据D和模型w。通过指定似然函数P(D∣w)和模型的先验分布P(w)，我们旨在找到最大化后验概率P(w∣D)的模型，该概率通过贝叶斯定理导出，如下所示：
在偏好学习中，对w的分布的优势在于捕捉人类偏好中固有的不确定性，从而提供不仅是单一的“最佳猜测”，而是一系列合理的模型。
偏好征集问题
偏好征集是决策理论的一个关键组成部分，旨在根据现有数据确定决策者的选择。在这项研究中，我们通过将模型适配到部分偏好集来解决偏好征集问题。在我们的案例中，偏好以最直接的形式表达：成对比较。为了说明这个概念，考虑一组水果，用F表示，包括苹果、香蕉、橙子、荔枝和芒果。
在我们的背景下，备选集A包括使用F集合中的一个或多个成分制作的所有可能的冰沙。
用户通过一组有序对(A,B)来表达他们的偏好，其中A严格优于B。
本文的后续部分将介绍一系列函数族，特别选择用于捕捉用户偏好的函数：加法函数。这些数学构造为理解不同因素如何影响个体偏好提供了一个简单但强大的框架，从而实现了通过成对比较表达的选择的有效建模。
加法模型
线性加法模型是可以用于捕捉用户偏好的最直接的模型。
加法线性模型
加法效用模型是将特定权重分配给我们集合中的每个单独成分的模型。然后，通过对其组成成分的权重求和来计算冰沙的整体效用或“好感度”。形式上，给定一组权重向量
由成分子集A制作的冰沙的效用是：
其中I是一个检测我是否在A中的身份函数。
带二进制交互的加法模型
2-加法模型在1-加法模型的基础上引入了额外的复杂性层次。权重向量不仅包含每个单个成分的权重，还包括每对成分的权重。这使得模型能够捕捉到成对水果之间的协同效应，有效地识别两种成分的组合如何影响整体效用。形式上，权重向量 w 被扩展以包括除了单个成分之外的每对 (i,j) 的权重：
而使用2-加法线性模型，冰沙的效用由以下公式给出：
其中 F² 是单个成分和成对的集合。
n-加法模型
进一步扩展这个概念，n-加法模型提供了一个高度灵活的效用框架。在这个模型中，权重向量不仅考虑单个成分和成对，还扩展到包括最多 n 个成分的任意子集的权重。这种泛化允许模型同时捕捉多个成分之间的复杂关系和协同效应。
形式上，权重向量 w 扩展以包括最多 n 个成分的所有可能组合的权重：
这个 n-加法模型可以捕捉到成分之间的全部交互作用，使其成为理解复杂偏好结构的极其强大的工具。
为了本分析的目的，我们将限制在2-加法模型上，因为我们认为成分之间的偏好关系的复杂性不太可能超过成对交互作用。
通过解决Probit回归问题来学习偏好
传统回归模型输出实值预测，而我们的目标是预测二元偏好关系。
为了实现这一目标，我们修改了回归模型，将其输出为A相对于B的优先概率。然后我们导出一个适当的损失函数，以有效地将这个概率模型拟合到我们的数据上。
将一个值压缩到0和1之间的一种经典方法是使用Probit函数。Probit函数定义如下：
下图显示了它的形状：
通过将该函数应用于 f(A) 和 f(B) 之间的差异，如果 f(A) 显著超过 f(B)，我们的模型将产生接近1的概率。反之，如果 f(A) 大约等于 f(B)，它将产生接近0.5的概率。
因此，偏好引导问题可以重新表述为寻找一个最优权重向量 w，使得：
二元交叉熵（BCE）损失
二元交叉熵（BCE）损失，也称为对数损失，用作输出概率范围在0到1之间的分类模型的性能指标，通常用于二元分类任务。数学上，给定真实标签 y（为0或1）和预测概率 p，BCE 定义为：
合成数据生成
为了验证我们的方法，我们引入了一个生成合成数据的协议。
该过程首先随机抽样一个权重向量w。然后我们将其中的一些参数设为零，以引入一层现实和简单性。
在假设用户的偏好与这个抽样函数相一致的情况下，我们可以将其作为一个基准来评估我们预测模型的准确性。
def singletons_and_pairs(lst): singletons = [(x,) for x in lst] pairs = list(combinations(lst, 2)) return singletons + pairsingredients = ["o", "a", "b","l", "m"]model = singletons_and_pairs(ingredients)w = np.random.normal(0, 1, size = (len(model),))p = np.random.randint(0, 2, size = (len(model),))w = w * p然后我们使用以下函数将每个备选项编码为一个二进制向量，其中组成部分的顺序与模型参数中的顺序相同
def vectorize_smoothie(smoothie): arr = np.zeros(len(model)) print(list(smoothie)) for i in range(arr.shape[0]): if all(j in list(smoothie) for j in model[i]): arr[i] = 1 return arr然后，为了评估特定的冰沙，我们使用一个乘积
vectorize_smoothie("oa") @ w# 返回 w_a + w_o + w_oa为了构建我们的数据集，我们首先随机抽样一个权重向量w。接下来，我们生成一组冰沙，并根据抽样的权重对每个冰沙进行评估。对于每一对冰沙 A 和 B，其中 f(A)>f(B)，我们在数据集中添加一个相应的偏好。每个 A 和 B 之间的偏好都以以下向量的形式进行捕捉:
对于每对 A,B，其中 f(A) > f(B)，我们添加两行 v(A,B) 和 v(B,A)，第一行标记为类别 1，第二行标记为类别 0。
以下代码给出了一个包含 n 个冰沙的数据集。
def sample_dataset(n): ingredients = ["o", "a", "b","l", "m"] model = singletons_and_pairs(ingredients) X = [] y = [] w = sample_w(model) subsets = set() while len(subsets) != n: s = random_subset(ingredients) subsets.add(s) subsets = list(subsets) for i in range(len(subsets)-1): x_i = vectorize_smoothie(subsets[i]) for j in range(i+1, len(subsets)): x_j = vectorize_smoothie(subsets[j]) x1 = x_i - x_j x2 = x_j - x_i if f(subsets[i], w) == f(subsets[j], w): continue if f(subsets[i], w) > f(subsets[j], w): X.append(x1) X.append(x2) y.append(1) y.append(0) continue if f(subsets[i], w) < f(subsets[j], w): X.append(x1) X.append(x2) y.append(0) y.append(1) continue X = np.array(X) y = np.array(y) return X,y,w,model基于成本的解决方案
解决这个问题的一种方法是利用 BCE 损失的凸性以及像 Torch 这样的库。
我们首先将生成的数据封装到 PyTorch 提供的适当的数据集加载器中。
X,y,w,model = sample_dataset(30)X_tensor = torch.FloatTensor(X)y_tensor = torch.FloatTensor(y)dataset = TensorDataset(X_tensor, y_tensor)train_size = int(0.3 * len(dataset))test_size = len(dataset) - train_sizetrain_dataset, test_dataset = random_split(dataset, [train_size, test_size])train_loader = DataLoader(dataset=train_dataset, batch_size=32, shuffle=True)test_loader = DataLoader(dataset=test_dataset, batch_size=1, shuffle=False)现在，我们创建一个简单的线性模型
class BinaryClassifier(nn.Module): def __init__(self, input_dim): super(BinaryClassifier, self).__init__() self.fc1 = nn.Linear(input_dim, 1) def forward(self, x): x = torch.sigmoid(self.fc1(x)) return x然后，我们使用PyTorch的autograd功能来训练它。
input_dim = X.shape[1]model = BinaryClassifier(input_dim)# 损失函数和优化器criterion = nn.BCELoss()optimizer = optim.Adam(model.parameters(), lr=0.01)losses = []# 训练模型epochs = 200for epoch in range(epochs): for batch_idx, (data, target) in enumerate(train_loader): optimizer.zero_grad() output = model(data).squeeze() loss = criterion(output, target) loss.backward() optimizer.step()然后，我们使用测试数据集来测试所得到的模型
model.eval()with torch.no_grad(): correct = 0 total = 0 for data, target in test_loader: output = model(data).squeeze() predicted = (output > 0.5).float() total += target.size(0) correct += (predicted == target).sum().item() acc = correct / total accuracy.append(acc)if (epoch+1) % 50 == 0: print(f'Epoch [{epoch+1}/{epochs}], Loss: {loss.item():.4f}') print(f'Test Accuracy: {100 * correct / total:.2f}%')使用20%的数据进行训练，我们得到了约98.32%的准确率，这是相当不错的。
最大似然估计（MLE）
解决概率回归问题的另一种方法是显式地构建给定权重向量w的数据的似然。
我们假设模型产生一个概率 p，表示A相对于B的优先级。这种情况下的预测分布如下所示：
给定权重向量的一对（x，y）的似然表达如下：
数据集的概率为
似然值可以非常小，在将许多概率相乘时具有重要意义。这可能会导致数值下溢（其中非常小的浮点数被舍入为零）。取这些值的对数可以将它们转换为更易处理的数字，通常是负数且具有较大的幅度。
因此，对数似然为
您可能会注意到这个损失是BCE损失的负数，这就是为什么最大化似然等同于最小化BCE损失。
正则化技术
正则化是机器学习中的一项关键技术，用于解决过拟合问题，其中模型过度适应训练数据，包括其中的噪声，从而影响其在新数据上的性能。它通过向损失函数添加惩罚项来限制模型参数的复杂性。这促使模型更简单，既能拟合训练数据，又能保持模型的简洁性。
L1（Lasso）和L2（Ridge）是常见的正则化形式，每个形式都对模型的目标引入了独特的惩罚项。
L1根据参数的绝对值添加惩罚，导致具有一些权重为零的稀疏模型。
相比之下，L2对参数的平方幅值进行惩罚，使权重缩小但不为零。
L1（Lasso）和L2（Ridge）正则化技术在对模型参数进行惩罚的方式上有所不同。L1对绝对值进行惩罚，导致一些权重完全为零，便于特征选择。相比之下，L2对权重的平方幅值进行惩罚，确保它们保持小但通常不为零，保留了所有特征但降低了影响。
最大后验估计
如前所述，贝叶斯定理允许我们通过利用似然函数和选定的参数先验分布P(w)来估计模型参数的后验分布，表示为P(w∣X,y)。
本质上，先验封装了在观察任何数据之前对参数的初始信念或假设，而似然量化了参数对观察数据的解释程度。贝叶斯定理将这些元素结合起来，产生一个后验分布，代表了我们关于参数的更新信念，考虑了先验和数据的影响。
两个非常著名的先验分布是拉普拉斯（Laplace）先验和高斯（Gaussian）先验。
拉普拉斯先验基于这样的假设，即权重w服从参数μ=0和尺度参数b的拉普拉斯分布。
换句话说，它假设权重的分布以零为中心，并且随着值偏离这个中心点而呈指数衰减，反映出偏好许多权重设为零的更稀疏模型。
高斯先验基于这样的假设，即权重w服从均值μ=0和方差σ的高斯（或正态）分布。
本质上，它假设权重的分布在零的对称中心周围，呈钟形曲线，表明权重最有可能接近均值，随着移动到更远的值，不太可能发生。这导致对权重平滑正则化的偏好，确保它们在大小上保持较小，而不一定驱使它们中的任何一个完全为零。
对数后验如下：
通过优化我们的模型，我们发现最大化对数后验与最小化特定正则化损失的目标基本等价。
值得注意的是，L1和L2正则化之间的区别取决于所考虑的先验分布形式。
在MCMC方法中使用后验分布
在贝叶斯框架中，一切都以概率的方式对待。因此，贝叶斯线性回归不是像经典线性回归那样估计固定的回归系数值，而是估计可能的系数值的分布。
使用后验分布的一种方法是从分布P(w|X,y)中采样一组权重。
一种简单的方法是使用MCMC方法，理解MCMC方法的起点是Metropolis-Hasting方法。
Metropolis-Hasting方法
Metropolis-Hastings（M-H）算法是贝叶斯统计学中用于从复杂概率分布中采样的方法。
它使用一个更简单的“提议分布”来探索目标分布，根据计算的概率接受或拒绝样本。值得注意的是，M-H算法不需要准确目标分布的知识；具有与之成比例的分布就足够了。
我们不会使用它，因为其他方法更可靠和高效，但我们仍会简要解释它的工作原理，因为M-H算法是一种基础的MCMC方法。
- 选择一个初始猜测
- 设置一个提议分布，通常是以当前值w为中心的高斯分布。
然后对于每一次迭代，我们按照以下步骤进行：
- 从提案分布 P(w’|w) 中采样一个新的 w’。
- 计算接受概率
- 从均匀分布 [0,1] 中随机抽取一个数 u。如果 u ≤ α，则接受 w’ 作为新的样本；否则保留 w。
NUTS 采样器和 pyMC3
Metropolis-Hastings 方法涉及在参数空间中提出一个新的点，然后根据其似然性与当前点的似然性进行比较来决定是否接受这个新点。其效率在很大程度上取决于提案分布的选择，并且在高维空间中可能出现随机游走行为，导致收敛速度缓慢。
NUTS（No-U-Turn Sampler）是 Hamiltonian Monte Carlo（HMC）方法的扩展。NUTS 不是进行随机游走，而是利用目标分布的梯度信息提出“跳跃步”，从而更高效地遍历分布。它的主要优势之一是自动确定最佳的跳跃步数，从而避免了随机游走问题和手动调整的繁琐任务。
PyMC3 是一个流行的概率编程框架，可以无缝集成这两种方法（以及其他方法），使用户能够轻松拟合复杂的贝叶斯模型，而不必陷入底层算法的复杂性中。
在我们的案例中，这段代码将从后验分布 P(w|X,y) 中采样一系列权重。
import pymc3 as pmwith pm.Model() as probit_model: # 权重和偏差的先验分布 weights = pm.Normal('weights', mu=0, sd=4, shape=X.shape[1]) bias = pm.Normal('bias', mu=0, sd=4) # Probit 链接函数 mu = pm.math.dot(X, weights) + bias phi = pm.math.invprobit(mu) # Inverse probit 链接函数 # 似然函数 y_obs = pm.Bernoulli('y_obs', p=phi, observed=y) # 从后验分布中采样 trace = pm.sample(5000, tune=1000, chains=5, target_accept = 0.90)我们可以绘制每个权重的不同分布。
我们可以看到每个权重都收敛到一个高斯分布。因此，现在每个预测可以以概率的方式进行，预测的分布也将是一个高斯分布。
例如，我们虚构的决策者对于橙子冰沙、橙子-苹果冰沙和香蕉-苹果冰沙的喜好分别由以下高斯分布给出。
使用生成数据的模型，我们可以看到这三种冰沙的真实效用分别为-0.66、-0.24和0.79，因此高斯分布实际上很好地反映了偏好和它们之间的差距。
结论
在这篇博文中，我们从偏好获取的复杂性到贝叶斯线性回归模型的复杂性进行了探讨。我们的讨论始于探索二次可加模型，这是一种既逼真又易于计算的捕捉用户偏好的方法。通过从基本线性回归过渡到更高级的 probit 模型，我们提供了一个新的视角来理解偏好数据。
我们还深入探讨了基于成本的观点和概率观点之间的等价性，揭示了最小化二元交叉熵损失与最大化似然性的类比，以及正则化作为隐含选择先验的作用。
最后，我们讨论了贝叶斯框架在生成不仅仅是点估计而是整个预测分布方面的实用性。这种方法为我们的模型提供了更高的置信度和可解释性，特别适用于偏好学习这种微妙的任务。
有了这个基础，未来的研究可以更深入地探索将这些复杂模型应用于越来越复杂和大规模的偏好数据。
致谢
特别感谢我的同事/朋友Anissa Hacene对这项工作的贡献，以及TDS团队对我们的及时审核和富有洞察力的评论。
