每日算法 — 使用java实现支持向量机:核函数与SMO序列最小优化

支持向量机(Support Vector Machine, SVM)是机器学习领域最具理论美感的算法之一。其核心思想是在特征空间中寻找一个最优超平面,使得不同类别样本之间的间隔(Margin)最大化。本文将用Java从零实现SVM,涵盖线性可分SVM、核技巧处理非线性问题,以及经典的SMO(Sequential Minimal Optimization)优化算法。

SVM基本原理

最大间隔分类器

给定训练集 ${(\mathbf{x}i, y_i)}{i=1}^n$,其中 $y_i \in {-1, +1}$。SVM寻找超平面 $\mathbf{w}^T \mathbf{x} + b = 0$ 使得分类间隔最大。通过拉格朗日乘子法,原问题转化为对偶问题:

$$\max_{lpha} \sum_{i=1}^n lpha_i – rac{1}{2}\sum_{i=1}^n\sum_{j=1}^n lpha_i lpha_j y_i y_j K(\mathbf{x}_i, \mathbf{x}_j)$$

约束条件:$\sum_{i=1}^n lpha_i y_i = 0$,$0 \le lpha_i \le C$

核函数

当数据线性不可分时,核函数 $K(\mathbf{x}_i, \mathbf{x}_j) = \phi(\mathbf{x}_i)^T \phi(\mathbf{x}_j)$ 将数据隐式映射到高维空间。常用核函数包括:

  • 线性核:$K(u,v) = u^T v$
  • 多项式核:$K(u,v) = (\gamma u^T v + r)^d$
  • RBF核:$K(u,v) = \exp(-\gamma |u-v|^2)$

SMO算法

SMO由Platt于1998年提出,每次仅优化两个拉格朗日乘子 $lpha_i$ 和 $lpha_j$,其余固定。这使得子问题有解析解,避免了调用数值优化库。

Java实现

import java.util.*;

/**
 * SVM支持向量机实现
 * 包含:线性核/RBF核、SMO优化、预测与评估
 */
public class SVM {

    // 训练样本
    static class Sample {
        double[] x; // 特征向量
        int y;      // 标签:+1 或 -1
        Sample(double[] x, int y) { this.x = x.clone(); this.y = y; }
    }

    // 核函数接口
    interface Kernel {
        double compute(double[] a, double[] b);
    }

    static class LinearKernel implements Kernel {
        public double compute(double[] a, double[] b) {
            double sum = 0;
            for (int i = 0; i < a.length; i++) sum += a[i] * b[i];
            return sum;
        }
    }

    static class RBFKernel implements Kernel {
        final double gamma;
        RBFKernel(double gamma) { this.gamma = gamma; }
        public double compute(double[] a, double[] b) {
            double sum = 0;
            for (int i = 0; i < a.length; i++) {
                double d = a[i] - b[i];
                sum += d * d;
            }
            return Math.exp(-gamma * sum);
        }
    }

    // SVM模型参数
    List<Sample> samples;
    double[] alpha;    // 拉格朗日乘子
    double b;          // 偏置
    Kernel kernel;
    double C;          // 惩罚参数
    double tol;        // 容差
    int maxIter;       // 最大迭代次数

    // 缓存误差
    double[] errorCache;

    public SVM(Kernel kernel, double C, double tol, int maxIter) {
        this.kernel = kernel;
        this.C = C;
        this.tol = tol;
        this.maxIter = maxIter;
    }

    /**
     * 计算模型输出 f(x) = sum(alpha_i * y_i * K(x_i, x)) + b
     */
    double evaluate(double[] x) {
        double sum = 0;
        for (int i = 0; i < samples.size(); i++) {
            if (alpha[i] > 0) {
                sum += alpha[i] * samples.get(i).y * kernel.compute(samples.get(i).x, x);
            }
        }
        return sum + b;
    }

    /**
     * 计算第k个样本的预测误差 E_k = f(x_k) - y_k
     */
    double computeError(int k) {
        double f = evaluate(samples.get(k).x);
        double e = f - samples.get(k).y;
        errorCache[k] = e;
        return e;
    }

    /**
     * 选择第二个优化变量:启发式选择使得 |E_i - E_j| 最大的j
     */
    int selectJ(int i, double Ei) {
        int j = -1;
        double maxDelta = 0;
        // 优先从非边界样本中选择
        for (int k = 0; k < samples.size(); k++) {
            if (alpha[k] > 0 && alpha[k] < C) {
                double Ej = errorCache[k];
                double delta = Math.abs(Ei - Ej);
                if (delta > maxDelta && k != i) {
                    maxDelta = delta;
                    j = k;
                }
            }
        }
        if (j >= 0) return j;
        // 随机选择一个不等于i的j
        do { j = (int)(Math.random() * samples.size()); } while (j == i);
        return j;
    }

    /**
     * 裁剪alpha到[L, H]区间
     */
    double clip(double a, double L, double H) {
        return Math.max(L, Math.min(H, a));
    }

    /**
     * 单步SMO:优化alpha[i]和alpha[j]
     * 返回是否进行了有效更新
     */
    boolean takeStep(int i, int j) {
        if (i == j) return false;
        Sample si = samples.get(i);
        Sample sj = samples.get(j);
        double alphaI = alpha[i], alphaJ = alpha[j];
        int yi = si.y, yj = sj.y;

        double Ei = computeError(i);
        double Ej = computeError(j);

        // 计算边界L和H
        double L, H;
        if (yi != yj) {
            L = Math.max(0, alphaJ - alphaI);
            H = Math.min(C, C + alphaJ - alphaI);
        } else {
            L = Math.max(0, alphaI + alphaJ - C);
            H = Math.min(C, alphaI + alphaJ);
        }
        if (L == H) return false;

        double kii = kernel.compute(si.x, si.x);
        double kjj = kernel.compute(sj.x, sj.x);
        double kij = kernel.compute(si.x, sj.x);
        double eta = 2 * kij - kii - kjj;

        double alphaJNew;
        if (eta >= 0) {
            // 目标函数沿该方向非凸,尝试边界点
            double f1 = yi * Ei - alphaI * kii * yi * yi - alphaJ * kij * yi * yj;
            double f2 = yj * Ej - alphaI * kij * yi * yi - alphaJ * kjj * yi * yj;
            double Lobj = f1 * L + f2 * L + 0.5 * kii * L * L + 0.5 * kjj * L * L + eta * L * L;
            double Hobj = f1 * H + f2 * H + 0.5 * kii * H * H + 0.5 * kjj * H * H + eta * H * H;
            alphaJNew = (Lobj > Hobj + tol) ? L : ((Hobj > Lobj + tol) ? H : alphaJ);
        } else {
            alphaJNew = alphaJ - yj * (Ei - Ej) / eta;
            alphaJNew = clip(alphaJNew, L, H);
        }

        if (Math.abs(alphaJNew - alphaJ) < 1e-5) return false;

        double alphaINew = alphaI + yi * yj * (alphaJ - alphaJNew);

        // 更新b
        double b1 = b - Ei - yi * (alphaINew - alphaI) * kii - yj * (alphaJNew - alphaJ) * kij;
        double b2 = b - Ej - yi * (alphaINew - alphaI) * kij - yj * (alphaJNew - alphaJ) * kjj;

        if (alphaINew > 0 && alphaINew < C) b = b1;
        else if (alphaJNew > 0 && alphaJNew < C) b = b2;
        else b = (b1 + b2) / 2;

        alpha[i] = alphaINew;
        alpha[j] = alphaJNew;
        errorCache[i] = errorCache[j] = 0;
        return true;
    }

    /**
     * 外层循环:遍历所有alpha,选择违反KKT条件的进行优化
     */
    boolean examineExample(int i) {
        Sample si = samples.get(i);
        double Ei = computeError(i);
        double r = Ei * si.y;
        // KKT条件检查
        if ((r < -tol && alpha[i] < C) || (r > tol && alpha[i] > 0)) {
            // 优先从非边界样本中选择j
            int[] nonBound = new int[samples.size()];
            int nbCount = 0;
            for (int k = 0; k < samples.size(); k++)
                if (alpha[k] > 0 && alpha[k] < C) nonBound[nbCount++] = k;

            if (nbCount > 0) {
                int j = selectJ(i, Ei);
                if (takeStep(i, j)) return true;
            }

            // 随机遍历非边界样本
            int start = (int)(Math.random() * nbCount);
            for (int k = 0; k < nbCount; k++) {
                int idx = nonBound[(start + k) % nbCount];
                if (takeStep(i, idx)) return true;
            }

            // 遍历所有样本
            start = (int)(Math.random() * samples.size());
            for (int k = 0; k < samples.size(); k++) {
                int idx = (start + k) % samples.size();
                if (idx != i && takeStep(i, idx)) return true;
            }
        }
        return false;
    }

    public void fit(List<Sample> data) {
        this.samples = data;
        int n = data.size();
        alpha = new double[n];
        errorCache = new double[n];
        b = 0;

        int numChanged = 0;
        boolean examineAll = true;
        int iter = 0;

        while (iter < maxIter && (numChanged > 0 || examineAll)) {
            numChanged = 0;
            if (examineAll) {
                for (int i = 0; i < n; i++)
                    if (examineExample(i)) numChanged++;
            } else {
                for (int i = 0; i < n; i++)
                    if (alpha[i] > 0 && alpha[i] < C)
                        if (examineExample(i)) numChanged++;
            }
            iter++;
            examineAll = !examineAll && numChanged == 0;
        }
        System.out.println("SMO收敛,迭代次数:" + iter);
    }

    public int predict(double[] x) {
        return evaluate(x) >= 0 ? 1 : -1;
    }

    // 主函数:在简单的二维数据集上测试
    public static void main(String[] args) {
        List<Sample> data = new ArrayList<>();
        // 类别 +1
        data.add(new Sample(new double[]{2, 3}, 1));
        data.add(new Sample(new double[]{3, 3}, 1));
        data.add(new Sample(new double[]{3, 4}, 1));
        data.add(new Sample(new double[]{4, 4}, 1));
        // 类别 -1
        data.add(new Sample(new double[]{0, 1}, -1));
        data.add(new Sample(new double[]{1, 0}, -1));
        data.add(new Sample(new double[]{1, 1}, -1));
        data.add(new Sample(new double[]{0, 0}, -1));

        SVM svm = new SVM(new RBFKernel(0.5), 1.0, 1e-3, 1000);
        svm.fit(data);

        // 测试
        int correct = 0;
        for (Sample s : data) {
            int pred = svm.predict(s.x);
            if (pred == s.y) correct++;
            System.out.println("真实:" + s.y + " 预测:" + pred);
        }
        System.out.println("准确率:" + (100.0 * correct / data.size()) + "%");

        // 新样本预测
        System.out.println("预测 (2.5, 2.5): " + svm.predict(new double[]{2.5, 2.5}));
    }
}

算法复杂度分析

  • SMO子问题:每次优化两个变量有解析解,时间复杂度 $O(1)$(忽略核函数计算)。
  • 核函数计算:RBF核计算两个样本需要 $O(d)$,$d$ 为特征维度。
  • 总体训练复杂度:最坏情况下每次外层循环需要 $O(n^2)$ 次核计算,总体约为 $O(n^2 \cdot d \cdot ext{iter})$。对于大规模数据,需要采用核缓存或近似方法。
  • 预测复杂度:$O(s \cdot d)$,$s$ 为支持向量数量,通常 $s \ll n$。

扩展方向

  1. 核缓存:将频繁访问的核矩阵元素缓存到内存中,避免重复计算。
  2. 多分类扩展:通过One-vs-One或One-vs-Rest策略扩展到多分类场景。
  3. 概率输出:通过Platt缩放将SVM输出转化为概率估计。

总结

SVM通过最大间隔原则和核技巧,在有限样本下具有出色的泛化能力。SMO算法将复杂的二次规划问题分解为一系列可解析求解的二变量子问题,是SVM高效训练的关键。Java实现中通过启发式变量选择和KKT条件检查,确保了算法的收敛速度和数值稳定性。