尧图网页设计 
简介这是一套基于Java的数值分析、线性代数与机器学习算法设计源码包定位为面向科研人员、工程师、学生及机器学习开发者的数学工具集。源码包覆盖数值积分、微分方程求解、插值与最优化等数值分析方法提供矩阵乘法、行列式计算、矩阵分解、特征值求解等线性代数工具并集成线性回归、逻辑回归、决策树、聚类及支持向量机等基础机器学习算法便于用户快速搭建实验环境、验证算法思路。包内共45个文件以36个Java源文件为主体另有文本说明文件与图片辅助示意压缩包约544KB目录按测试、矩阵、数值分析、机器学习等模块组织结构清晰便于检索。目前已有341人学习/使用适合需要在Java环境下高效完成数学计算和模型原型的开发者可大幅降低高级数学与机器学习技术的应用门槛提升工程效率。1. 这套源码不是课程设计凑数Java实现数值分析、线性代数与机器学习先回答“值不值”当Java后端要在没有Python环境、没有外部算法服务的进程里独立完成矩阵分解、线性回归或主成分降维时“数值分析、线性代数与机器学习”这三个词就从课本目录变成了实打实的代码逻辑。这个标题下的源码价值不在“写得多漂亮”而在于让你把算法吃透并且能交付一个可运行的jar。它不是教科书里那种只做增删改查的课程设计案例源码而是能在离线推荐、风控评分、报表拟合里当真工具的东西。适合谁想搞懂算法而不仅调API的Java工程师转机器学习方向的求职者以及需要一份能讲清楚源码的在校生。2. 从零搭算法库的工程骨架包划分、Maven依赖与矩阵抽象的三层设计动手写算法之前先决定目录长什么样。这决定了后面几千行代码会不会在两个月后变成一坨改不动的黑匣子。我的建议是把项目当成别人会来接手的东西来搭而不是当成自己的草稿纸。2.1 包划分的底层逻辑为什么算法模块必须和数学原语分离我用过的最顺手的结构是这样包名直接按职责隔离com.example.mathlib/ ├── linalg/ │ ├── Matrix.java │ ├── DenseMatrix.java │ └── LUDecomposition.java ├── numerical/ │ ├── JacobiSolver.java │ ├── PowerIteration.java │ └── OdeSolver.java ├── ml/ │ ├── LinearRegression.java │ ├── GradientDescent.java │ └── CrossValidator.java ├── util/ │ ├── MatrixIO.java │ └── Statistics.java └── benchmark/ └── MatrixBenchmark.java依赖方向必须是单向的ml调用numerical和linalgnumerical调用linalg谁也不许反向依赖。我第一次写这种库时把矩阵运算和回归逻辑放进同一个包后来又顺手塞了归一化、交叉验证进去整个包直接变成杂物间——代码都在但谁都不敢动。后来按依赖方向重写才敢在数值模块里加新算法而不担心砸了机器学习的实验。这种划分不是强迫症。机器学习算法要跑在矩阵分解和求解器之上数值算法要跑在矩阵存储之上底层存储一换上层不该跟着改。包边界就是这层缓冲。如果你只打算交课程作业这种结构也能让你在答辩时多讲两句“为什么这么设计”效果比贴几百行循环强。2.2 先写矩阵接口double[][]在算法库里有哪四个不可控问题用double[][]当然能写数值分析课程作业大多这么干。但做一个要持续生长的算法源码库二维数组有四个让你后半夜翻车的毛病它保存的是行对象的引用转置、切片、拷贝全得手工实现稍不注意就在原矩阵上改出副作用没有维度检查get(row, col)越界直接扔一个底层异常调试时看不出是谁传错了所有元素固定占满内存稀疏矩阵想省内存没门将来想接原生BLAS做加速数组头的内存布局不是JNI喜欢的样子。所以我一般先定一个极小的矩阵接口够用就好不贪多public interface Matrix { int rows(); int cols(); double get(int row, int col); void set(int row, int col, double value); void swapRows(int rowA, int rowB); Matrix times(Matrix other); double[] times(double[] vector); Matrix transpose(); Matrix copy(); }逻辑说明多数数值算法真正依赖的操作没那么多——读写元素、交换行、矩阵乘矩阵、矩阵乘向量、转置、拷贝。把这六个能力收进接口实现类可以自由选择存储方式上层算法只面向接口编程。swapRows非常重要高斯消元和LU分解都靠它换主元没有这个操作实现类就没法做选主元优化。参数说明接口里有意的没加getRow、getColumn这类便利方法。接口每多一个方法所有实现类都要跟着变。我宁愿在DenseMatrix里单独提供rowArray(int)返回拷贝也不把它推进接口保持抽象层的稳定。稠密实现我一般用一维数组下标从二维折成一维public class DenseMatrix implements Matrix { private final int rows; private final int cols; private final double[] data; public DenseMatrix(int rows, int cols) { this.rows rows; this.cols cols; this.data new double[rows * cols]; } Override public double get(int row, int col) { return data[row * cols col]; } Override public void set(int row, int col, double value) { data[row * cols col] value; } Override public Matrix times(Matrix other) { DenseMatrix result new DenseMatrix(rows, other.cols()); for (int i 0; i rows; i) { for (int k 0; k cols; k) { double aik get(i, k); if (aik 0.0) continue; for (int j 0; j other.cols(); j) { result.data[i * other.cols() j] aik * other.get(k, j); } } } return result; } Override public double[] times(double[] vector) { double[] out new double[rows]; for (int i 0; i rows; i) { for (int k 0; k cols; k) { out[i] get(i, k) * vector[k]; } } return out; } }逻辑说明乘法我特意写成i-k-j循环顺序而不是i-j-k。外层行循环固定后先遍历k让内层连续访问other的第k行缓存命中率比i-j-k高。if (aik 0.0) continue;对稀疏这种小矩阵没什么用但对含大量零元素的矩阵能省下不少乘法。实现类加了rowArray(int)做边界检查方法后上层算法就能更安全地批量读行。2.3 Maven依赖怎么选commons-math3做对照自家实现做核心这套源码的定位是“核心算法自己写第三方库只做对照和测试基准”。如果直接引入ND4J或Smile业务代码只剩一堆调用看起来包了一层很深的壳但面试或答辩时一问“LU分解怎么做”还是接不住。自研实现和成熟库做对比测试是确认自己没写错的底线手段。pom.xml我一般这么给project modelVersion4.0.0/modelVersion groupIdcom.example/groupId artifactIdmathlib/artifactId version1.0.0/version properties maven.compiler.release11/maven.compiler.release project.build.sourceEncodingUTF-8/project.build.sourceEncoding /properties dependencies dependency groupIdorg.apache.commons/groupId artifactIdcommons-math3/artifactId version3.6.1/version scopecompile/scope /dependency dependency groupIdorg.junit.jupiter/groupId artifactIdjunit-jupiter/artifactId version5.9.2/version scopetest/scope /dependency /dependencies /project参数说明maven.compiler.release11指定编译目标为Java 11而不是继承IDE默认的JDK版本这一步能省掉后面大量“在本机好好的换台机器就UnsupportedClassVersionError”的沟通成本。commons-math3用3.6.1而不是更新的版本是因为这版API足够稳定行为也够经典做对照测试时可信度高。这里不要求你有很深的Java基础但理解依赖scope的意义很重要test范围的JUnit不会进到交付jar里而commons-math3如果没有人为排除会在打包时被带进去。3. 把线性代数和数值分析写进JavaLU分解、迭代法与特征计算的实现与参数这一章是数值分析和线性代数交汇的地方也是整套源码里最容易写错的部分。我们做三件事LU分解解稠密方程组、雅可比迭代解稀疏对角线占优方程组、幂迭代求主特征值。它们覆盖了“直接法”“迭代法”“特征计算”三类最常用的数值分析场景。3.1 LU分解与高斯消元最小可运行实现与选主元阈值LU分解的目标是把方阵A拆成“置换矩阵P乘A等于L乘U”的形式。L是下三角U是上三角。常见做法是用部分选主元partial pivoting保证数值稳定也就是在每一列里挑绝对值最大的行来当主元行。public class LUDecomposition { private final Matrix lu; private final int[] pivot; private int pivotSign 1; public LUDecomposition(Matrix a) { int m a.rows(); int n a.cols(); lu a.copy(); pivot new int[m]; for (int i 0; i m; i) { pivot[i] i; } for (int col 0; col Math.min(m, n); col) { // 部分选主元找当前列绝对值最大的行抑制误差放大 int maxRow col; double maxAbs Math.abs(lu.get(col, col)); for (int row col 1; row m; row) { double abs Math.abs(lu.get(row, col)); if (abs maxAbs) { maxAbs abs; maxRow row; } } if (maxAbs 1e-12) { continue; // 接近奇异跳过当前列行列式近似为0 } lu.swapRows(col, maxRow); // 记录换行并翻转符号供行列式计算 int tmp pivot[col]; pivot[col] pivot[maxRow]; pivot[maxRow] tmp; pivotSign -pivotSign; double pivotVal lu.get(col, col); for (int row col 1; row m; row) { double factor lu.get(row, col) / pivotVal; lu.set(row, col, factor); // L部分直接存进下三角 for (int k col 1; k n; k) { lu.set(row, k, lu.get(row, k) - factor * lu.get(col, k)); } } } } public double[] solve(double[] b) { int n lu.rows(); double[] x b.clone(); // 先按pivot做行置换 for (int i 0; i n; i) { x[i] b[pivot[i]]; } // 前代解Ly Pb for (int i 0; i n; i) { double sum x[i]; for (int j 0; j i; j) { sum - lu.get(i, j) * x[j]; } x[i] sum; } // 回代解Ux y for (int i n - 1; i 0; i--) { double sum x[i]; for (int j i 1; j n; j) { sum - lu.get(i, j) * x[j]; } x[i] sum / lu.get(i, i); } return x; } public double determinant() { double det pivotSign; for (int i 0; i lu.rows(); i) { det * lu.get(i, i); } return det; } }逻辑说明分解时我们把L的乘数直接覆盖到原矩阵的下三角位置U覆盖到上三角位置对角线留给U。这种“原地分解”把两块矩阵合到同一块内存里省一半空间。solve里先做行置换再前代解LyPb最后回代解Uxy整个过程不需要显式求逆矩阵。参数说明选主元的阈值1e-12是从double机器精度推算的经验值。double能表示的精度大约在2.2e-16量级消元过程会让误差放大几个数量级取1e-12给中间误差留了余量。但这个阈值不能拍脑袋往大了调我试过1e-8结果把一个条件数在10^8附近的病态矩阵误判成奇异矩阵解出来的方程结果完全不能用。真正的工程做法是先看矩阵量纲再决定阈值。3.2 雅可比迭代解方程组收敛条件与容差参数怎么给直接法的问题在于矩阵一上万规模O(n^3)的复杂度直接把耗时拉到不可接受。这时候迭代法的价值就出来了。雅可比迭代是最容易实现的一种每次用上一轮的整个向量算这一轮的新值。public class JacobiSolver { public static double[] solve(double[][] a, double[] b, double tol, int maxIter) { int n b.length; double[] x new double[n]; double[] next new double[n]; for (int iter 0; iter maxIter; iter) { for (int i 0; i n; i) { double sum b[i]; for (int j 0; j n; j) { if (j ! i) { sum - a[i][j] * x[j]; } } next[i] sum / a[i][i]; } // 相对二范数误差判断是否收敛 double diff 0.0; double base 0.0; for (int i 0; i n; i) { double d next[i] - x[i]; diff d * d; base next[i] * next[i]; } System.arraycopy(next, 0, x, 0, n); if (base 0 Math.sqrt(diff / base) tol) { break; } } return x; } }逻辑说明雅可比迭代的收敛充分条件是矩阵严格对角占优也就是每一行对角元素绝对值大于该行其他元素绝对值之和。满足这个条件时x会逐步逼近真实解不满足时迭代可能发散或者震荡。next数组必须每次重新算完再整体替换这是和“高斯-赛德尔”最大的区别后者在循环内直接覆盖当前分量收敛更快但每一步都依赖最新值不好并行。参数说明tol1e-8是我常用的默认值对大多数业务场景足够而且迭代次数不会爆到离谱。tol1e-12理论上更精确但double在方程组里的累积误差本来就在这个量级容易让迭代停在“永远达不到”的阈值附近白耗大量算力。base 0的判断是防止解全为零向量时出现除零这个不加第一轮就会翻车。3.3 幂迭代求特征值机器学习里SVD的第一步特征值分解是线性代数的重头戏机器学习里主成分分析PCA本质就是对协方差矩阵做特征分解。虽然生产库用的是QR算法或Lanczos算法但从理解角度讲幂迭代是把“特征向量”这个概念落到代码里最直观的入口。public class PowerIteration { private static final java.util.Random RND new java.util.Random(42); private final Matrix matrix; public PowerIteration(Matrix matrix) { this.matrix matrix; } public double dominantEigenValue(int maxIter, double tol) { int n matrix.cols(); double[] v new double[n]; for (int i 0; i n; i) { v[i] RND.nextDouble(); } normalize(v); double lambda 0.0; for (int iter 0; iter maxIter; iter) { double[] w matrix.times(v); double newLambda dot(v, w); if (Math.abs(newLambda - lambda) tol) { break; } lambda newLambda; normalize(w); v w; } return lambda; } private static double dot(double[] a, double[] b) { double s 0.0; for (int i 0; i a.length; i) s a[i] * b[i]; return s; } private static void normalize(double[] v) { double norm 0.0; for (double value : v) norm value * value; norm Math.sqrt(norm); if (norm 0.0) return; for (int i 0; i v.length; i) { v[i] / norm; } } }逻辑说明反复用矩阵作用在向量上向量方向会不断向最大特征值对应的特征向量靠拢。每次乘完更新lambda时我用的是瑞利商vᵀAv / vᵀv在v已归一化的情况下就是vᵀAv。这比直接取w的最大分量要稳收敛判据也更平滑。参数说明初始向量用固定种子42的随机数是为了实验结果可复现。实际使用中如果矩阵对称且初始向量和主特征向量接近正交收敛会极慢甚至看起来不收敛。幂迭代只能求出绝对值最大的特征值如果你要做SVD降维更好的是直接对协方差矩阵做正交三角化或者用Java里成熟的commons-math3的SingularValueDecomposition幂迭代留在做特征值快照和教学演示里更合适。4. 机器学习算法落地线性回归的三条路径、梯度下降与交叉验证到了这一章线性代数和数值分析的工具才开始真正为机器学习服务。我把机器学习算法的落地路径拆成三段正规方程、梯度下降、交叉验证。很多刚接触的人以为机器学习就是“调个库”看完这一章你会明白核心算法在Java里落地每一步都是在跟矩阵和数值精度打交道。4.1 正规方程还是梯度下降特征量级决定你的选择线性回归是最简单的机器学习模型也是理解“算法设计源码”最合适的起点。常见做法是用最小二乘法目标是最小化||Xw - y||²。两条路一条是直接解正规方程XᵀX w Xᵀy一条是梯度下降。public class LinearRegression { private double[] w; public void fit(double[][] x, double[] y) { int m x.length; int n x[0].length; double[][] xt transpose(x); double[][] xtx multiply(xt, x); double[] xty multiplyVector(xt, y); // 解线性方程组而不是显式求逆 LUDecomposition lu new LUDecomposition(new DenseMatrix(xtx)); w lu.solve(xty); } public double predict(double[] sample) { double s 0.0; for (int j 0; j w.length; j) { s w[j] * sample[j]; } return s; } }逻辑说明正规方程把拟合问题转成一次矩阵求解不需要迭代也不需要调学习率。用LU分解解XᵀX w Xᵀy比先算(XᵀX)⁻¹再乘Xᵀy稳定得多。显式求逆会把条件数的影响放大——原本10的8次方的病态矩阵求逆后误差可能直接毁掉结果。参数说明正规方程适合特征维度几百以内的中小场景。特征上万时XᵀX的构造和分解直接吃掉O(n^3)的时间和内存这时候必须梯度下降。另一个前提是特征量级不能太悬殊如果某一列是0.001另一列是1000XᵀX的对角元素差出6个数量级LU分解也会皱眉头。遇到这种情况先做特征缩放再看用什么求解器。4.2 梯度下降的四个参数学习率、批量、归一化与收敛判据特征量级一大线性规划就切到梯度下降。写一套能用的梯度下降不难难的是选对四个参数学习率、批量、归一化方式、收敛判据。public class GradientDescent { public static double[] batchTrain(double[][] x, double[] y, double lr, int maxEpochs, double tol) { int n y.length; int d x[0].length; double[] w new double[d]; for (int epoch 0; epoch maxEpochs; epoch) { double[] grad new double[d]; for (int i 0; i n; i) { double err predict(w, x[i]) - y[i]; for (int j 0; j d; j) { grad[j] err * x[i][j]; } } double gradNorm 0.0; for (int j 0; j d; j) { grad[j] / n; // 取平均使学习率不依赖样本量 gradNorm grad[j] * grad[j]; } gradNorm Math.sqrt(gradNorm); if (gradNorm tol) { break; } for (int j 0; j d; j) { w[j] - lr * grad[j]; } } return w; } private static double predict(double[] w, double[] sample) { double s 0.0; for (int j 0; j w.length; j) { s w[j] * sample[j]; } return s; } }逻辑说明批量梯度下降每一轮用全量样本计算梯度方向再按学习率更新权重。将梯度除以样本量n是为了让学习率对数据集大小不敏感同样的lr1e-3在1000条样本和100万条样本上行为才不会差出天和地。参数说明如下表参数取值建议说明lr1e-3起步过大导致震荡不收敛过小收敛肉眼不可见maxEpochs1e4到1e6取决于数据量和学习率设太小容易“没跑完就停”tol1e-6梯度二范数低于它认为进了平缓区批量全量或32-128全量稳定但每轮贵小批量快但噪声大这组参数里最玄学的是学习率。我的血泪经验是第一次跑发现损失变成NaN先别查矩阵运算十有八九是lr大了发现损失下降慢得像蜗牛也别急着加轮次先调大一倍lr。特征归一化比调参更能根治问题因为归一化能让损失函数等高线接近圆形梯度下降就不会走“之字形”。z-score归一化我用得最多public static double[] zScore(double[] col) { double mean 0.0; for (double v : col) mean v; mean / col.length; double var 0.0; for (double v : col) var (v - mean) * (v - mean); var / col.length; double std Math.sqrt(var); if (std 1e-8) std 1.0; // 常数列不缩放 double[] out new double[col.length]; for (int i 0; i col.length; i) { out[i] (col[i] - mean) / std; } return out; }参数说明std 1e-8的判断是针对“这一列全是同一个数”的情况。如果某列是常量方差为0直接除会得到一堆NaN这是机器学习里最典型的翻车点之一。把它兜住后面才能安心训练。这也是Java面试里经常被追问的一个细节归一化、标准化用的到底是什么量。4.3 K折交叉验证评估模型泛化性的最小可跑框架同一个模型训练集上得分漂亮测试集上可能一塌糊涂。要判断这个模型到底是学到规律还是背下答案就得做交叉验证。K折交叉验证是套路最标准的做法把数据乱序切成K份轮流拿其中一份当测试集其余K-1份当训练集最后把K次误差取平均。public class CrossValidator { public static double kFoldMSE(double[][] x, double[] y, int k, boolean useNormalEquation) { int n y.length; Integer[] idx new Integer[n]; for (int i 0; i n; i) idx[i] i; java.util.Collections.shuffle( java.util.Arrays.asList(idx), new java.util.Random(2024)); int foldSize n / k; double totalMse 0.0; for (int fold 0; fold k; fold) { int testFrom fold * foldSize; int testTo (fold k - 1) ? n : testFrom foldSize; double[][] trainX new double[n - (testTo - testFrom)][]; double[] trainY new double[n - (testTo - testFrom)]; double[][] testX new double[testTo - testFrom][]; double[] testY new double[testTo - testFrom]; int tr 0, te 0; for (int pos 0; pos n; pos) { int sampleIdx idx[pos]; if (pos testFrom pos testTo) { testX[te] x[sampleIdx]; testY[te] y[sampleIdx]; te; } else { trainX[tr] x[sampleIdx]; trainY[tr] y[sampleIdx]; tr; } } LinearRegression model new LinearRegression(); model.fit(trainX, trainY); double foldMse 0.0; for (int i 0; i testX.length; i) { double err model.predict(testX[i]) - testY[i]; foldMse err * err; } totalMse foldMse / testX.length; } return totalMse / k; } }逻辑说明用Collections.shuffle对索引乱序保证每折里的样本分布接近整体。fold k - 1时把testTo设为n是因为n不一定能被k整除最后一折要吃掉剩余样本。模型在每一折里重新fit这个代价逃不掉所以交叉验证比单次训练要慢K倍K不要拍脑袋设太大。参数说明Random(2024)是固定种子。固定种子意味着“这组数据每次跑出的折划分一致”课程设计和实验复现强烈建议这么做。真实业务里我会换不同种子多跑几次报告均值和方差因为单次划分可能运气好也可能运气差均值方差才是模型真实水平的估计。K通常取5或10K10偏差更小但耗时更高数据量大到跑不动时K5是够用的妥协。5. 避坑指南数值稳定性、内存与调试的5条血泪教训写算法源码真正的门槛不在“算法没听过”而在“结果不对却找不到哪儿错”。这一章把我在Java里做数值分析和机器学习时踩过的坑总结成五条每条都是“现象、原因、解决”三段式照着排查能省一整晚。5.1 现象同样的数据Java算出来和Python差一截double精度与条件数现象Python里numpy.linalg.solve和Java里自研LU求解同一个方程组结果第四位小数就开始分叉有时甚至对不上前三位的数。原因两个问题叠加。第一浮点运算本身有舍入误差这是天生的第二方程组矩阵的条件数大会把输入数据和中间计算的微小误差放大。条件数是矩阵最大奇异值和最小奇异值的比值比值越大矩阵越病态解对误差越敏感。选主元策略不同也会带来差异但主因还是条件数。解决先算一下条件数心里有底再选算法。工程上最实用的兜底是“列缩放”把每列归一到同一量级降低条件数再做分解。public static double[] columnScale(double[][] a) { int n a[0].length; double[] scale new double[n]; for (int j 0; j n; j) { double max 0.0; for (double[] row : a) { max Math.max(max, Math.abs(row[j])); } scale[j] (max 0) ? 1.0 / max : 1.0; } return scale; }逻辑说明用每列绝对值的最大值倒数做缩放因子。求解前先把矩阵每列乘上scale[j]解出结果后再对应除回就能把某些列数值过大的影响压住。这不会改变方程组的解但能把条件数拉下来一个数量级甚至更多。5.2 现象1万乘1万的double[][]直接OOM从二维数组到CRS稀疏存储现象创建10000 * 10000的double[]数组内存眼睁睁看着涨到800MB然后OutOfMemoryError: Java heap space砸在脸上。原因double[][]并不是一块连续内存而是“一个包含10000个引用的一维数组每个引用又指向一个长度10000的一维数组”。每个一维数组都有对象头对象头开销巨大外层数组和内存碎片会额外吃几十MB。更关键的是真实业务里的矩阵大多是稀疏的——几万个元素里可能只有几十个非零值用二维数组存成满阵纯属浪费。解决有两种思路。对稠密矩阵改成一维数组连续存储已经在上文DenseMatrix里实现了对稀疏矩阵用CRS压缩行存储public class SparseMatrix implements Matrix { private final int rows; private final int cols; private final int[] rowPtr; // 每个起始行在colIdx/val中的位置 private final int[] colIdx; // 每行的列下标 private final double[] val; // 非零值 public SparseMatrix(int rows, int cols, int[] rowPtr, int[] colIdx, double[] val) { this.rows rows; this.cols cols; this.rowPtr rowPtr; this.colIdx colIdx; this.val val; } Override public double get(int row, int col) { int start rowPtr[row]; int end rowPtr[row 1]; for (int p start; p end; p) { if (colIdx[p] col) return val[p]; } return 0.0; } Override public void set(int row, int col, double value) { throw new UnsupportedOperationException( CRS不适合随机写请先构建好结构再创建矩阵); } }逻辑说明CRS用三个数组存储非零元素val存值colIdx存列编号rowPtr记录每一行的非零元素在val里从哪里开始。get方法需要线性扫描行内元素所以CRS适合“构建一次、反复做矩阵乘向量”的场景不适合频繁随机写。这也是我在接口设计里保留get而不强推set的原因——稀疏实现对写操作天然不友好。参数说明如果矩阵规模还要大get的线性扫描是性能瓶颈可以对每行里的colIdx做二分查找代价是构建时保持列有序。这些细节不用一开始就做但要在接口设计时留出空间。5.3 现象计算中途都是NaN却不知道哪一步翻车NaN传播与逐步验证现象训练完线性回归预测值全是NaN打印中间矩阵发现某些位置早就变成NaN程序却不报错、不中断。原因NaN的传播是静默的。某一步出现0.0 / 0.0或者Infinity - Infinity结果就是NaN之后所有涉及NaN的运算全被污染但不会抛出异常整个计算结果就成了一堆不可读的垃圾。最常见的来源有三个梯度下降学习率太大导致中间值超过double上限变成Infinity归一化时遇到std0矩阵分解里对角元素为0导致除零。解决在关键运算入口加一个checkFinite断言把“危险数据”提前暴露出来public static void checkFinite(double[] values, String where) { for (double v : values) { if (!Double.isFinite(v)) { throw new ArithmeticException( where : non-finite value v); } } }逻辑说明在每次矩阵乘向量、每次梯度更新之后调用一次比如checkFinite(w, gradient after epoch epoch)。这样程序一旦出现NaN或Inf立刻在第一个现场抛出异常而不是等到最后输出一堆问号。生产环境下这个检查要按需关闭否则每次迭代都扫描全数组会有性能开销但在开发和调试阶段它就是后悔药。5.4 现象JDK版本升级后符号找不到Maven依赖与字节码版本的坑现象项目一直在JDK 11上跑某天把IDE默认JDK切换成17编译正常一运行就报UnsupportedClassVersionError或NoClassDefFoundError看起来像代码被“幽灵改过”。原因UnsupportedClassVersionError是class文件的字节码版本和运行JDK不匹配编译出的class是55JDK 11用JDK 17能跑但反过来说如果编译目标升到61JDK 17而运行环境还是11就会炸。NoClassDefFoundError则通常是某个依赖的传递依赖用了更高版本字节码或者Maven在解析时拉到了不兼容的版本。解决pom里显式指定maven.compiler.release不要依赖IDE。mvn -version mvn dependency:tree -Dincludesorg.apache.commons mvn clean package参数说明maven.compiler.release11是编译目标和运行兼容性的双保险。dependency:tree能让你看清哪个依赖悄悄带进了不兼容的传递依赖如果发现commons-math3内部又拉了别的库可以在pom里用exclusions排除。这类问题靠眼睛看代码永远看不出来Maven依赖树才是第一排查工具。5.5 现象单元测试里断言浮点结果为什么assertEquals永远会挂现象用JUnit写assertEquals(2.0, det)明明手算行列式就是2测试却红得刺眼把期望值改成2.0000000001才通过。原因JUnit的assertEquals(double, double)在没有第三个参数时执行的是精确比较。浮点运算没有精确任何一次乘法和加法的舍入都可能让结果尾数差几个1e-16。尤其矩阵分解重排了计算顺序结果和手算路径不同差一点是正常的。解决用三参数的assertEquals(expected, actual, delta)并给一个合理的误差上限Test void testDeterminant() { double det new LUDecomposition(matrix).determinant(); assertEquals(2.0, det, 1e-8); }参数说明误差上限1e-8是经验值。矩阵规模小、运算次数少时误差在1e-12量级矩阵规模一大或条件数一高误差会涨到1e-6都可能。断言精度宁可给松一点让测试反映“数值上等价”而不是“bit级相等”否则每个commit都在跟浮点噪声较劲。6. 让这套源码真正值钱用JMH基准测试和真实数据集说服团队代码能跑只是起点能让团队愿意把这套源码放进去用得拿出数据说话。我的习惯是加两类东西JMH性能基准和一个能端到端运行的真实数据集验证入口。6.1 JMH基准测试自研矩阵乘法到底比Apache Commons慢多少自研算法最尴尬的时刻就是被问“凭什么不用现成的库”。与其嘴上解释不如把两段代码跑进JMH。JMH是JVM生态的微基准测试标准它解决的核心问题是JIT预热一段方法第一次跑很慢跑热之后会快几十倍用System.currentTimeMillis()手写计时测出来的全是噪声。BenchmarkMode(Mode.AverageTime) OutputTimeUnit(TimeUnit.MICROSECONDS) Warmup(iterations 3, time 1) Measurement(iterations 5, time 1) Fork(1) public class MatrixBenchmark { private Matrix myMatrix new DenseMatrix(128, 128); private RealMatrix commonsMatrix new Array2DRowRealMatrix(128, 128); Benchmark public Matrix myImpl() { return myMatrix.times(myMatrix); } Benchmark public RealMatrix commonsImpl() { return commonsMatrix.multiply(commonsMatrix); } }参数说明Warmup(iterations3, time1)是预热3轮每轮1秒让JIT把热点代码编译成机器码Measurement(iterations5, time1)是正式测量5轮每轮1秒取平均耗时Fork(1)让测试在独立的JVM进程里跑避免当前进程里的其他上下文干扰结果。矩阵维度定在128而不是10000是因为基准测试要测的是算法本身在同规模下的相对差距而不是比谁先OOM。跑完你会面对两种结果如果自研比对照库慢出10倍以上别自欺欺人老老实实换库如果只慢20%到50%说明问题多半出在数组拷贝或循环顺序上照着i-k-j缓存优化方向改一版再测值得。这个结论贴到代码评审里比任何说明书都有说服力。6.2 端到端验证用UCI数据集跑通线性回归基准测试验证性能端到端验证正确性。我一般会在util.MatrixIO里加一个读CSV的工具然后写一个入口方法读一份公开数据集把样本先做z-score归一化再分别用正规方程和梯度下降跑5折交叉验证对比两个MSE。public static void main(String[] args) { double[][] raw MatrixIO.readCsv(new File(args[0])); double[][] x new double[raw.length][raw[0].length - 1]; double[] y new double[raw.length]; for (int i 0; i raw.length; i) { y[i] raw[i][raw[i].length - 1]; System.arraycopy(raw[i], 0, x[i], 0, raw[i].length - 1); } x normalizeAllColumns(x); double normalMse CrossValidator.kFoldMSE(x, y, 5, true); double gdMse CrossValidator.kFoldMSE(x, y, 5, false); System.out.println(5-fold MSE, normal equation : normalMse); System.out.println(5-fold MSE, gradient descent: gdMse); }逻辑说明最后一行代码是最关键的习惯——在终端留下两组可复现的数字。它们证明三件事正则方程能解、梯度下降能收敛、交叉验证能跑通整套框架。这个main方法保持能独立运行的状态比堆十个单元测试更直观。参数说明normalizeAllColumns会对每一列做一次zScore处理并把均值、标准差存在对象里方便对未来新样本执行同样的缩放这个细节很容易漏。如果跳过这一步梯度下降那个结果几乎必然比正规方程差出好几个数量级你以为算法写错其实只是没归一化。这套源码走到这里已经不是课程作业的规模了。我自己后来的习惯是每加一个算法先写对比测试再写基准测试再跑一个公开数据集做端到端验证三步缺一不可。这样做确实比“写完就算收工”要多花半天但换来的是一套别人能接手、自己半年后还敢重构的源码。共勉希望帮到你。本文还有配套的精品资源点击获取 Java实现数值分析与机器学习：从LU分解到梯度下降的工程源码
简介这是一套基于Java的数值分析、线性代数与机器学习算法设计源码包定位为面向科研人员、工程师、学生及机器学习开发者的数学工具集。源码包覆盖数值积分、微分方程求解、插值与最优化等数值分析方法提供矩阵乘法、行列式计算、矩阵分解、特征值求解等线性代数工具并集成线性回归、逻辑回归、决策树、聚类及支持向量机等基础机器学习算法便于用户快速搭建实验环境、验证算法思路。包内共45个文件以36个Java源文件为主体另有文本说明文件与图片辅助示意压缩包约544KB目录按测试、矩阵、数值分析、机器学习等模块组织结构清晰便于检索。目前已有341人学习/使用适合需要在Java环境下高效完成数学计算和模型原型的开发者可大幅降低高级数学与机器学习技术的应用门槛提升工程效率。1. 这套源码不是课程设计凑数Java实现数值分析、线性代数与机器学习先回答“值不值”当Java后端要在没有Python环境、没有外部算法服务的进程里独立完成矩阵分解、线性回归或主成分降维时“数值分析、线性代数与机器学习”这三个词就从课本目录变成了实打实的代码逻辑。这个标题下的源码价值不在“写得多漂亮”而在于让你把算法吃透并且能交付一个可运行的jar。它不是教科书里那种只做增删改查的课程设计案例源码而是能在离线推荐、风控评分、报表拟合里当真工具的东西。适合谁想搞懂算法而不仅调API的Java工程师转机器学习方向的求职者以及需要一份能讲清楚源码的在校生。2. 从零搭算法库的工程骨架包划分、Maven依赖与矩阵抽象的三层设计动手写算法之前先决定目录长什么样。这决定了后面几千行代码会不会在两个月后变成一坨改不动的黑匣子。我的建议是把项目当成别人会来接手的东西来搭而不是当成自己的草稿纸。2.1 包划分的底层逻辑为什么算法模块必须和数学原语分离我用过的最顺手的结构是这样包名直接按职责隔离com.example.mathlib/ ├── linalg/ │ ├── Matrix.java │ ├── DenseMatrix.java │ └── LUDecomposition.java ├── numerical/ │ ├── JacobiSolver.java │ ├── PowerIteration.java │ └── OdeSolver.java ├── ml/ │ ├── LinearRegression.java │ ├── GradientDescent.java │ └── CrossValidator.java ├── util/ │ ├── MatrixIO.java │ └── Statistics.java └── benchmark/ └── MatrixBenchmark.java依赖方向必须是单向的ml调用numerical和linalgnumerical调用linalg谁也不许反向依赖。我第一次写这种库时把矩阵运算和回归逻辑放进同一个包后来又顺手塞了归一化、交叉验证进去整个包直接变成杂物间——代码都在但谁都不敢动。后来按依赖方向重写才敢在数值模块里加新算法而不担心砸了机器学习的实验。这种划分不是强迫症。机器学习算法要跑在矩阵分解和求解器之上数值算法要跑在矩阵存储之上底层存储一换上层不该跟着改。包边界就是这层缓冲。如果你只打算交课程作业这种结构也能让你在答辩时多讲两句“为什么这么设计”效果比贴几百行循环强。2.2 先写矩阵接口double[][]在算法库里有哪四个不可控问题用double[][]当然能写数值分析课程作业大多这么干。但做一个要持续生长的算法源码库二维数组有四个让你后半夜翻车的毛病它保存的是行对象的引用转置、切片、拷贝全得手工实现稍不注意就在原矩阵上改出副作用没有维度检查get(row, col)越界直接扔一个底层异常调试时看不出是谁传错了所有元素固定占满内存稀疏矩阵想省内存没门将来想接原生BLAS做加速数组头的内存布局不是JNI喜欢的样子。所以我一般先定一个极小的矩阵接口够用就好不贪多public interface Matrix { int rows(); int cols(); double get(int row, int col); void set(int row, int col, double value); void swapRows(int rowA, int rowB); Matrix times(Matrix other); double[] times(double[] vector); Matrix transpose(); Matrix copy(); }逻辑说明多数数值算法真正依赖的操作没那么多——读写元素、交换行、矩阵乘矩阵、矩阵乘向量、转置、拷贝。把这六个能力收进接口实现类可以自由选择存储方式上层算法只面向接口编程。swapRows非常重要高斯消元和LU分解都靠它换主元没有这个操作实现类就没法做选主元优化。参数说明接口里有意的没加getRow、getColumn这类便利方法。接口每多一个方法所有实现类都要跟着变。我宁愿在DenseMatrix里单独提供rowArray(int)返回拷贝也不把它推进接口保持抽象层的稳定。稠密实现我一般用一维数组下标从二维折成一维public class DenseMatrix implements Matrix { private final int rows; private final int cols; private final double[] data; public DenseMatrix(int rows, int cols) { this.rows rows; this.cols cols; this.data new double[rows * cols]; } Override public double get(int row, int col) { return data[row * cols col]; } Override public void set(int row, int col, double value) { data[row * cols col] value; } Override public Matrix times(Matrix other) { DenseMatrix result new DenseMatrix(rows, other.cols()); for (int i 0; i rows; i) { for (int k 0; k cols; k) { double aik get(i, k); if (aik 0.0) continue; for (int j 0; j other.cols(); j) { result.data[i * other.cols() j] aik * other.get(k, j); } } } return result; } Override public double[] times(double[] vector) { double[] out new double[rows]; for (int i 0; i rows; i) { for (int k 0; k cols; k) { out[i] get(i, k) * vector[k]; } } return out; } }逻辑说明乘法我特意写成i-k-j循环顺序而不是i-j-k。外层行循环固定后先遍历k让内层连续访问other的第k行缓存命中率比i-j-k高。if (aik 0.0) continue;对稀疏这种小矩阵没什么用但对含大量零元素的矩阵能省下不少乘法。实现类加了rowArray(int)做边界检查方法后上层算法就能更安全地批量读行。2.3 Maven依赖怎么选commons-math3做对照自家实现做核心这套源码的定位是“核心算法自己写第三方库只做对照和测试基准”。如果直接引入ND4J或Smile业务代码只剩一堆调用看起来包了一层很深的壳但面试或答辩时一问“LU分解怎么做”还是接不住。自研实现和成熟库做对比测试是确认自己没写错的底线手段。pom.xml我一般这么给project modelVersion4.0.0/modelVersion groupIdcom.example/groupId artifactIdmathlib/artifactId version1.0.0/version properties maven.compiler.release11/maven.compiler.release project.build.sourceEncodingUTF-8/project.build.sourceEncoding /properties dependencies dependency groupIdorg.apache.commons/groupId artifactIdcommons-math3/artifactId version3.6.1/version scopecompile/scope /dependency dependency groupIdorg.junit.jupiter/groupId artifactIdjunit-jupiter/artifactId version5.9.2/version scopetest/scope /dependency /dependencies /project参数说明maven.compiler.release11指定编译目标为Java 11而不是继承IDE默认的JDK版本这一步能省掉后面大量“在本机好好的换台机器就UnsupportedClassVersionError”的沟通成本。commons-math3用3.6.1而不是更新的版本是因为这版API足够稳定行为也够经典做对照测试时可信度高。这里不要求你有很深的Java基础但理解依赖scope的意义很重要test范围的JUnit不会进到交付jar里而commons-math3如果没有人为排除会在打包时被带进去。3. 把线性代数和数值分析写进JavaLU分解、迭代法与特征计算的实现与参数这一章是数值分析和线性代数交汇的地方也是整套源码里最容易写错的部分。我们做三件事LU分解解稠密方程组、雅可比迭代解稀疏对角线占优方程组、幂迭代求主特征值。它们覆盖了“直接法”“迭代法”“特征计算”三类最常用的数值分析场景。3.1 LU分解与高斯消元最小可运行实现与选主元阈值LU分解的目标是把方阵A拆成“置换矩阵P乘A等于L乘U”的形式。L是下三角U是上三角。常见做法是用部分选主元partial pivoting保证数值稳定也就是在每一列里挑绝对值最大的行来当主元行。public class LUDecomposition { private final Matrix lu; private final int[] pivot; private int pivotSign 1; public LUDecomposition(Matrix a) { int m a.rows(); int n a.cols(); lu a.copy(); pivot new int[m]; for (int i 0; i m; i) { pivot[i] i; } for (int col 0; col Math.min(m, n); col) { // 部分选主元找当前列绝对值最大的行抑制误差放大 int maxRow col; double maxAbs Math.abs(lu.get(col, col)); for (int row col 1; row m; row) { double abs Math.abs(lu.get(row, col)); if (abs maxAbs) { maxAbs abs; maxRow row; } } if (maxAbs 1e-12) { continue; // 接近奇异跳过当前列行列式近似为0 } lu.swapRows(col, maxRow); // 记录换行并翻转符号供行列式计算 int tmp pivot[col]; pivot[col] pivot[maxRow]; pivot[maxRow] tmp; pivotSign -pivotSign; double pivotVal lu.get(col, col); for (int row col 1; row m; row) { double factor lu.get(row, col) / pivotVal; lu.set(row, col, factor); // L部分直接存进下三角 for (int k col 1; k n; k) { lu.set(row, k, lu.get(row, k) - factor * lu.get(col, k)); } } } } public double[] solve(double[] b) { int n lu.rows(); double[] x b.clone(); // 先按pivot做行置换 for (int i 0; i n; i) { x[i] b[pivot[i]]; } // 前代解Ly Pb for (int i 0; i n; i) { double sum x[i]; for (int j 0; j i; j) { sum - lu.get(i, j) * x[j]; } x[i] sum; } // 回代解Ux y for (int i n - 1; i 0; i--) { double sum x[i]; for (int j i 1; j n; j) { sum - lu.get(i, j) * x[j]; } x[i] sum / lu.get(i, i); } return x; } public double determinant() { double det pivotSign; for (int i 0; i lu.rows(); i) { det * lu.get(i, i); } return det; } }逻辑说明分解时我们把L的乘数直接覆盖到原矩阵的下三角位置U覆盖到上三角位置对角线留给U。这种“原地分解”把两块矩阵合到同一块内存里省一半空间。solve里先做行置换再前代解LyPb最后回代解Uxy整个过程不需要显式求逆矩阵。参数说明选主元的阈值1e-12是从double机器精度推算的经验值。double能表示的精度大约在2.2e-16量级消元过程会让误差放大几个数量级取1e-12给中间误差留了余量。但这个阈值不能拍脑袋往大了调我试过1e-8结果把一个条件数在10^8附近的病态矩阵误判成奇异矩阵解出来的方程结果完全不能用。真正的工程做法是先看矩阵量纲再决定阈值。3.2 雅可比迭代解方程组收敛条件与容差参数怎么给直接法的问题在于矩阵一上万规模O(n^3)的复杂度直接把耗时拉到不可接受。这时候迭代法的价值就出来了。雅可比迭代是最容易实现的一种每次用上一轮的整个向量算这一轮的新值。public class JacobiSolver { public static double[] solve(double[][] a, double[] b, double tol, int maxIter) { int n b.length; double[] x new double[n]; double[] next new double[n]; for (int iter 0; iter maxIter; iter) { for (int i 0; i n; i) { double sum b[i]; for (int j 0; j n; j) { if (j ! i) { sum - a[i][j] * x[j]; } } next[i] sum / a[i][i]; } // 相对二范数误差判断是否收敛 double diff 0.0; double base 0.0; for (int i 0; i n; i) { double d next[i] - x[i]; diff d * d; base next[i] * next[i]; } System.arraycopy(next, 0, x, 0, n); if (base 0 Math.sqrt(diff / base) tol) { break; } } return x; } }逻辑说明雅可比迭代的收敛充分条件是矩阵严格对角占优也就是每一行对角元素绝对值大于该行其他元素绝对值之和。满足这个条件时x会逐步逼近真实解不满足时迭代可能发散或者震荡。next数组必须每次重新算完再整体替换这是和“高斯-赛德尔”最大的区别后者在循环内直接覆盖当前分量收敛更快但每一步都依赖最新值不好并行。参数说明tol1e-8是我常用的默认值对大多数业务场景足够而且迭代次数不会爆到离谱。tol1e-12理论上更精确但double在方程组里的累积误差本来就在这个量级容易让迭代停在“永远达不到”的阈值附近白耗大量算力。base 0的判断是防止解全为零向量时出现除零这个不加第一轮就会翻车。3.3 幂迭代求特征值机器学习里SVD的第一步特征值分解是线性代数的重头戏机器学习里主成分分析PCA本质就是对协方差矩阵做特征分解。虽然生产库用的是QR算法或Lanczos算法但从理解角度讲幂迭代是把“特征向量”这个概念落到代码里最直观的入口。public class PowerIteration { private static final java.util.Random RND new java.util.Random(42); private final Matrix matrix; public PowerIteration(Matrix matrix) { this.matrix matrix; } public double dominantEigenValue(int maxIter, double tol) { int n matrix.cols(); double[] v new double[n]; for (int i 0; i n; i) { v[i] RND.nextDouble(); } normalize(v); double lambda 0.0; for (int iter 0; iter maxIter; iter) { double[] w matrix.times(v); double newLambda dot(v, w); if (Math.abs(newLambda - lambda) tol) { break; } lambda newLambda; normalize(w); v w; } return lambda; } private static double dot(double[] a, double[] b) { double s 0.0; for (int i 0; i a.length; i) s a[i] * b[i]; return s; } private static void normalize(double[] v) { double norm 0.0; for (double value : v) norm value * value; norm Math.sqrt(norm); if (norm 0.0) return; for (int i 0; i v.length; i) { v[i] / norm; } } }逻辑说明反复用矩阵作用在向量上向量方向会不断向最大特征值对应的特征向量靠拢。每次乘完更新lambda时我用的是瑞利商vᵀAv / vᵀv在v已归一化的情况下就是vᵀAv。这比直接取w的最大分量要稳收敛判据也更平滑。参数说明初始向量用固定种子42的随机数是为了实验结果可复现。实际使用中如果矩阵对称且初始向量和主特征向量接近正交收敛会极慢甚至看起来不收敛。幂迭代只能求出绝对值最大的特征值如果你要做SVD降维更好的是直接对协方差矩阵做正交三角化或者用Java里成熟的commons-math3的SingularValueDecomposition幂迭代留在做特征值快照和教学演示里更合适。4. 机器学习算法落地线性回归的三条路径、梯度下降与交叉验证到了这一章线性代数和数值分析的工具才开始真正为机器学习服务。我把机器学习算法的落地路径拆成三段正规方程、梯度下降、交叉验证。很多刚接触的人以为机器学习就是“调个库”看完这一章你会明白核心算法在Java里落地每一步都是在跟矩阵和数值精度打交道。4.1 正规方程还是梯度下降特征量级决定你的选择线性回归是最简单的机器学习模型也是理解“算法设计源码”最合适的起点。常见做法是用最小二乘法目标是最小化||Xw - y||²。两条路一条是直接解正规方程XᵀX w Xᵀy一条是梯度下降。public class LinearRegression { private double[] w; public void fit(double[][] x, double[] y) { int m x.length; int n x[0].length; double[][] xt transpose(x); double[][] xtx multiply(xt, x); double[] xty multiplyVector(xt, y); // 解线性方程组而不是显式求逆 LUDecomposition lu new LUDecomposition(new DenseMatrix(xtx)); w lu.solve(xty); } public double predict(double[] sample) { double s 0.0; for (int j 0; j w.length; j) { s w[j] * sample[j]; } return s; } }逻辑说明正规方程把拟合问题转成一次矩阵求解不需要迭代也不需要调学习率。用LU分解解XᵀX w Xᵀy比先算(XᵀX)⁻¹再乘Xᵀy稳定得多。显式求逆会把条件数的影响放大——原本10的8次方的病态矩阵求逆后误差可能直接毁掉结果。参数说明正规方程适合特征维度几百以内的中小场景。特征上万时XᵀX的构造和分解直接吃掉O(n^3)的时间和内存这时候必须梯度下降。另一个前提是特征量级不能太悬殊如果某一列是0.001另一列是1000XᵀX的对角元素差出6个数量级LU分解也会皱眉头。遇到这种情况先做特征缩放再看用什么求解器。4.2 梯度下降的四个参数学习率、批量、归一化与收敛判据特征量级一大线性规划就切到梯度下降。写一套能用的梯度下降不难难的是选对四个参数学习率、批量、归一化方式、收敛判据。public class GradientDescent { public static double[] batchTrain(double[][] x, double[] y, double lr, int maxEpochs, double tol) { int n y.length; int d x[0].length; double[] w new double[d]; for (int epoch 0; epoch maxEpochs; epoch) { double[] grad new double[d]; for (int i 0; i n; i) { double err predict(w, x[i]) - y[i]; for (int j 0; j d; j) { grad[j] err * x[i][j]; } } double gradNorm 0.0; for (int j 0; j d; j) { grad[j] / n; // 取平均使学习率不依赖样本量 gradNorm grad[j] * grad[j]; } gradNorm Math.sqrt(gradNorm); if (gradNorm tol) { break; } for (int j 0; j d; j) { w[j] - lr * grad[j]; } } return w; } private static double predict(double[] w, double[] sample) { double s 0.0; for (int j 0; j w.length; j) { s w[j] * sample[j]; } return s; } }逻辑说明批量梯度下降每一轮用全量样本计算梯度方向再按学习率更新权重。将梯度除以样本量n是为了让学习率对数据集大小不敏感同样的lr1e-3在1000条样本和100万条样本上行为才不会差出天和地。参数说明如下表参数取值建议说明lr1e-3起步过大导致震荡不收敛过小收敛肉眼不可见maxEpochs1e4到1e6取决于数据量和学习率设太小容易“没跑完就停”tol1e-6梯度二范数低于它认为进了平缓区批量全量或32-128全量稳定但每轮贵小批量快但噪声大这组参数里最玄学的是学习率。我的血泪经验是第一次跑发现损失变成NaN先别查矩阵运算十有八九是lr大了发现损失下降慢得像蜗牛也别急着加轮次先调大一倍lr。特征归一化比调参更能根治问题因为归一化能让损失函数等高线接近圆形梯度下降就不会走“之字形”。z-score归一化我用得最多public static double[] zScore(double[] col) { double mean 0.0; for (double v : col) mean v; mean / col.length; double var 0.0; for (double v : col) var (v - mean) * (v - mean); var / col.length; double std Math.sqrt(var); if (std 1e-8) std 1.0; // 常数列不缩放 double[] out new double[col.length]; for (int i 0; i col.length; i) { out[i] (col[i] - mean) / std; } return out; }参数说明std 1e-8的判断是针对“这一列全是同一个数”的情况。如果某列是常量方差为0直接除会得到一堆NaN这是机器学习里最典型的翻车点之一。把它兜住后面才能安心训练。这也是Java面试里经常被追问的一个细节归一化、标准化用的到底是什么量。4.3 K折交叉验证评估模型泛化性的最小可跑框架同一个模型训练集上得分漂亮测试集上可能一塌糊涂。要判断这个模型到底是学到规律还是背下答案就得做交叉验证。K折交叉验证是套路最标准的做法把数据乱序切成K份轮流拿其中一份当测试集其余K-1份当训练集最后把K次误差取平均。public class CrossValidator { public static double kFoldMSE(double[][] x, double[] y, int k, boolean useNormalEquation) { int n y.length; Integer[] idx new Integer[n]; for (int i 0; i n; i) idx[i] i; java.util.Collections.shuffle( java.util.Arrays.asList(idx), new java.util.Random(2024)); int foldSize n / k; double totalMse 0.0; for (int fold 0; fold k; fold) { int testFrom fold * foldSize; int testTo (fold k - 1) ? n : testFrom foldSize; double[][] trainX new double[n - (testTo - testFrom)][]; double[] trainY new double[n - (testTo - testFrom)]; double[][] testX new double[testTo - testFrom][]; double[] testY new double[testTo - testFrom]; int tr 0, te 0; for (int pos 0; pos n; pos) { int sampleIdx idx[pos]; if (pos testFrom pos testTo) { testX[te] x[sampleIdx]; testY[te] y[sampleIdx]; te; } else { trainX[tr] x[sampleIdx]; trainY[tr] y[sampleIdx]; tr; } } LinearRegression model new LinearRegression(); model.fit(trainX, trainY); double foldMse 0.0; for (int i 0; i testX.length; i) { double err model.predict(testX[i]) - testY[i]; foldMse err * err; } totalMse foldMse / testX.length; } return totalMse / k; } }逻辑说明用Collections.shuffle对索引乱序保证每折里的样本分布接近整体。fold k - 1时把testTo设为n是因为n不一定能被k整除最后一折要吃掉剩余样本。模型在每一折里重新fit这个代价逃不掉所以交叉验证比单次训练要慢K倍K不要拍脑袋设太大。参数说明Random(2024)是固定种子。固定种子意味着“这组数据每次跑出的折划分一致”课程设计和实验复现强烈建议这么做。真实业务里我会换不同种子多跑几次报告均值和方差因为单次划分可能运气好也可能运气差均值方差才是模型真实水平的估计。K通常取5或10K10偏差更小但耗时更高数据量大到跑不动时K5是够用的妥协。5. 避坑指南数值稳定性、内存与调试的5条血泪教训写算法源码真正的门槛不在“算法没听过”而在“结果不对却找不到哪儿错”。这一章把我在Java里做数值分析和机器学习时踩过的坑总结成五条每条都是“现象、原因、解决”三段式照着排查能省一整晚。5.1 现象同样的数据Java算出来和Python差一截double精度与条件数现象Python里numpy.linalg.solve和Java里自研LU求解同一个方程组结果第四位小数就开始分叉有时甚至对不上前三位的数。原因两个问题叠加。第一浮点运算本身有舍入误差这是天生的第二方程组矩阵的条件数大会把输入数据和中间计算的微小误差放大。条件数是矩阵最大奇异值和最小奇异值的比值比值越大矩阵越病态解对误差越敏感。选主元策略不同也会带来差异但主因还是条件数。解决先算一下条件数心里有底再选算法。工程上最实用的兜底是“列缩放”把每列归一到同一量级降低条件数再做分解。public static double[] columnScale(double[][] a) { int n a[0].length; double[] scale new double[n]; for (int j 0; j n; j) { double max 0.0; for (double[] row : a) { max Math.max(max, Math.abs(row[j])); } scale[j] (max 0) ? 1.0 / max : 1.0; } return scale; }逻辑说明用每列绝对值的最大值倒数做缩放因子。求解前先把矩阵每列乘上scale[j]解出结果后再对应除回就能把某些列数值过大的影响压住。这不会改变方程组的解但能把条件数拉下来一个数量级甚至更多。5.2 现象1万乘1万的double[][]直接OOM从二维数组到CRS稀疏存储现象创建10000 * 10000的double[]数组内存眼睁睁看着涨到800MB然后OutOfMemoryError: Java heap space砸在脸上。原因double[][]并不是一块连续内存而是“一个包含10000个引用的一维数组每个引用又指向一个长度10000的一维数组”。每个一维数组都有对象头对象头开销巨大外层数组和内存碎片会额外吃几十MB。更关键的是真实业务里的矩阵大多是稀疏的——几万个元素里可能只有几十个非零值用二维数组存成满阵纯属浪费。解决有两种思路。对稠密矩阵改成一维数组连续存储已经在上文DenseMatrix里实现了对稀疏矩阵用CRS压缩行存储public class SparseMatrix implements Matrix { private final int rows; private final int cols; private final int[] rowPtr; // 每个起始行在colIdx/val中的位置 private final int[] colIdx; // 每行的列下标 private final double[] val; // 非零值 public SparseMatrix(int rows, int cols, int[] rowPtr, int[] colIdx, double[] val) { this.rows rows; this.cols cols; this.rowPtr rowPtr; this.colIdx colIdx; this.val val; } Override public double get(int row, int col) { int start rowPtr[row]; int end rowPtr[row 1]; for (int p start; p end; p) { if (colIdx[p] col) return val[p]; } return 0.0; } Override public void set(int row, int col, double value) { throw new UnsupportedOperationException( CRS不适合随机写请先构建好结构再创建矩阵); } }逻辑说明CRS用三个数组存储非零元素val存值colIdx存列编号rowPtr记录每一行的非零元素在val里从哪里开始。get方法需要线性扫描行内元素所以CRS适合“构建一次、反复做矩阵乘向量”的场景不适合频繁随机写。这也是我在接口设计里保留get而不强推set的原因——稀疏实现对写操作天然不友好。参数说明如果矩阵规模还要大get的线性扫描是性能瓶颈可以对每行里的colIdx做二分查找代价是构建时保持列有序。这些细节不用一开始就做但要在接口设计时留出空间。5.3 现象计算中途都是NaN却不知道哪一步翻车NaN传播与逐步验证现象训练完线性回归预测值全是NaN打印中间矩阵发现某些位置早就变成NaN程序却不报错、不中断。原因NaN的传播是静默的。某一步出现0.0 / 0.0或者Infinity - Infinity结果就是NaN之后所有涉及NaN的运算全被污染但不会抛出异常整个计算结果就成了一堆不可读的垃圾。最常见的来源有三个梯度下降学习率太大导致中间值超过double上限变成Infinity归一化时遇到std0矩阵分解里对角元素为0导致除零。解决在关键运算入口加一个checkFinite断言把“危险数据”提前暴露出来public static void checkFinite(double[] values, String where) { for (double v : values) { if (!Double.isFinite(v)) { throw new ArithmeticException( where : non-finite value v); } } }逻辑说明在每次矩阵乘向量、每次梯度更新之后调用一次比如checkFinite(w, gradient after epoch epoch)。这样程序一旦出现NaN或Inf立刻在第一个现场抛出异常而不是等到最后输出一堆问号。生产环境下这个检查要按需关闭否则每次迭代都扫描全数组会有性能开销但在开发和调试阶段它就是后悔药。5.4 现象JDK版本升级后符号找不到Maven依赖与字节码版本的坑现象项目一直在JDK 11上跑某天把IDE默认JDK切换成17编译正常一运行就报UnsupportedClassVersionError或NoClassDefFoundError看起来像代码被“幽灵改过”。原因UnsupportedClassVersionError是class文件的字节码版本和运行JDK不匹配编译出的class是55JDK 11用JDK 17能跑但反过来说如果编译目标升到61JDK 17而运行环境还是11就会炸。NoClassDefFoundError则通常是某个依赖的传递依赖用了更高版本字节码或者Maven在解析时拉到了不兼容的版本。解决pom里显式指定maven.compiler.release不要依赖IDE。mvn -version mvn dependency:tree -Dincludesorg.apache.commons mvn clean package参数说明maven.compiler.release11是编译目标和运行兼容性的双保险。dependency:tree能让你看清哪个依赖悄悄带进了不兼容的传递依赖如果发现commons-math3内部又拉了别的库可以在pom里用exclusions排除。这类问题靠眼睛看代码永远看不出来Maven依赖树才是第一排查工具。5.5 现象单元测试里断言浮点结果为什么assertEquals永远会挂现象用JUnit写assertEquals(2.0, det)明明手算行列式就是2测试却红得刺眼把期望值改成2.0000000001才通过。原因JUnit的assertEquals(double, double)在没有第三个参数时执行的是精确比较。浮点运算没有精确任何一次乘法和加法的舍入都可能让结果尾数差几个1e-16。尤其矩阵分解重排了计算顺序结果和手算路径不同差一点是正常的。解决用三参数的assertEquals(expected, actual, delta)并给一个合理的误差上限Test void testDeterminant() { double det new LUDecomposition(matrix).determinant(); assertEquals(2.0, det, 1e-8); }参数说明误差上限1e-8是经验值。矩阵规模小、运算次数少时误差在1e-12量级矩阵规模一大或条件数一高误差会涨到1e-6都可能。断言精度宁可给松一点让测试反映“数值上等价”而不是“bit级相等”否则每个commit都在跟浮点噪声较劲。6. 让这套源码真正值钱用JMH基准测试和真实数据集说服团队代码能跑只是起点能让团队愿意把这套源码放进去用得拿出数据说话。我的习惯是加两类东西JMH性能基准和一个能端到端运行的真实数据集验证入口。6.1 JMH基准测试自研矩阵乘法到底比Apache Commons慢多少自研算法最尴尬的时刻就是被问“凭什么不用现成的库”。与其嘴上解释不如把两段代码跑进JMH。JMH是JVM生态的微基准测试标准它解决的核心问题是JIT预热一段方法第一次跑很慢跑热之后会快几十倍用System.currentTimeMillis()手写计时测出来的全是噪声。BenchmarkMode(Mode.AverageTime) OutputTimeUnit(TimeUnit.MICROSECONDS) Warmup(iterations 3, time 1) Measurement(iterations 5, time 1) Fork(1) public class MatrixBenchmark { private Matrix myMatrix new DenseMatrix(128, 128); private RealMatrix commonsMatrix new Array2DRowRealMatrix(128, 128); Benchmark public Matrix myImpl() { return myMatrix.times(myMatrix); } Benchmark public RealMatrix commonsImpl() { return commonsMatrix.multiply(commonsMatrix); } }参数说明Warmup(iterations3, time1)是预热3轮每轮1秒让JIT把热点代码编译成机器码Measurement(iterations5, time1)是正式测量5轮每轮1秒取平均耗时Fork(1)让测试在独立的JVM进程里跑避免当前进程里的其他上下文干扰结果。矩阵维度定在128而不是10000是因为基准测试要测的是算法本身在同规模下的相对差距而不是比谁先OOM。跑完你会面对两种结果如果自研比对照库慢出10倍以上别自欺欺人老老实实换库如果只慢20%到50%说明问题多半出在数组拷贝或循环顺序上照着i-k-j缓存优化方向改一版再测值得。这个结论贴到代码评审里比任何说明书都有说服力。6.2 端到端验证用UCI数据集跑通线性回归基准测试验证性能端到端验证正确性。我一般会在util.MatrixIO里加一个读CSV的工具然后写一个入口方法读一份公开数据集把样本先做z-score归一化再分别用正规方程和梯度下降跑5折交叉验证对比两个MSE。public static void main(String[] args) { double[][] raw MatrixIO.readCsv(new File(args[0])); double[][] x new double[raw.length][raw[0].length - 1]; double[] y new double[raw.length]; for (int i 0; i raw.length; i) { y[i] raw[i][raw[i].length - 1]; System.arraycopy(raw[i], 0, x[i], 0, raw[i].length - 1); } x normalizeAllColumns(x); double normalMse CrossValidator.kFoldMSE(x, y, 5, true); double gdMse CrossValidator.kFoldMSE(x, y, 5, false); System.out.println(5-fold MSE, normal equation : normalMse); System.out.println(5-fold MSE, gradient descent: gdMse); }逻辑说明最后一行代码是最关键的习惯——在终端留下两组可复现的数字。它们证明三件事正则方程能解、梯度下降能收敛、交叉验证能跑通整套框架。这个main方法保持能独立运行的状态比堆十个单元测试更直观。参数说明normalizeAllColumns会对每一列做一次zScore处理并把均值、标准差存在对象里方便对未来新样本执行同样的缩放这个细节很容易漏。如果跳过这一步梯度下降那个结果几乎必然比正规方程差出好几个数量级你以为算法写错其实只是没归一化。这套源码走到这里已经不是课程作业的规模了。我自己后来的习惯是每加一个算法先写对比测试再写基准测试再跑一个公开数据集做端到端验证三步缺一不可。这样做确实比“写完就算收工”要多花半天但换来的是一套别人能接手、自己半年后还敢重构的源码。共勉希望帮到你。本文还有配套的精品资源点击获取