矩阵求导与逆矩阵

矩阵求导是反向传播算法的数学基础——全连接层权重梯度 ∂L/∂W = δ·x^T、卷积核梯度、注意力层梯度，全都归结为矩阵对矩阵的导数运算。理解矩阵求导才能真正看懂框架 autograd 在做什么。逆矩阵在线性回归的正规方程解析解 w = (X^TX)^-1X^Ty、岭回归 (X^TX + λI)^-1、卡尔曼滤波、高斯过程推断中均不可或缺；Fisher 信息矩阵的逆给出参数估计的最优方差下界（Cramér-Rao 不等式）。实际工程中，大矩阵逆通常用 LU/Cholesky 分解或迭代法替代直接求逆以保证数值稳定。本节用具体数字逐步推导 Softmax Jacobian、全连接层反向传播梯度、正规方程等核心公式。

一、逆矩阵的定义

数学定义：对于方阵 A（n×n），如果存在矩阵 B 使得 AB = BA = I（单位矩阵），则 B 就是 A 的逆矩阵，记为 A⁻¹。

二、2×2 矩阵逆的手动计算

逐步数值计算：

三、逆矩阵的重要性质

四、行列式

数学定义：行列式是方阵到标量的一个映射，记为 det(A) 或 |A|。它衡量矩阵代表的线性变换对空间体积的缩放因子。

逐步数值计算（3×3）：

五、矩阵求导——反向传播的数学核心

矩阵求导（矩阵微积分）描述的是"损失函数对模型参数矩阵的导数怎么算"。反向传播就是用链式法则逐层算矩阵导数的过程。

5.1 标量对向量求导

逐步数值计算：

5.2 标量对矩阵求导

5.3 常用矩阵求导公式

六、Jacobian 矩阵

定义：当函数的输入和输出都是向量时，Jacobian 矩阵包含所有偏导数。

七、代码验证（C# / Rust）

C#（.NET 10）

// dotnet run 即可执行
// 3×3 行列式 & 逆矩阵 (余子式展开)
double Det3(double[,] m) =>
    m[0,0]*(m[1,1]*m[2,2]-m[1,2]*m[2,1])
   -m[0,1]*(m[1,0]*m[2,2]-m[1,2]*m[2,0])
   +m[0,2]*(m[1,0]*m[2,1]-m[1,1]*m[2,0]);
double[,] Inv3(double[,] m){ double d=Det3(m); return new double[,]{
    {(m[1,1]*m[2,2]-m[1,2]*m[2,1])/d,(m[0,2]*m[2,1]-m[0,1]*m[2,2])/d,(m[0,1]*m[1,2]-m[0,2]*m[1,1])/d},
    {(m[1,2]*m[2,0]-m[1,0]*m[2,2])/d,(m[0,0]*m[2,2]-m[0,2]*m[2,0])/d,(m[0,2]*m[1,0]-m[0,0]*m[1,2])/d},
    {(m[1,0]*m[2,1]-m[1,1]*m[2,0])/d,(m[0,1]*m[2,0]-m[0,0]*m[2,1])/d,(m[0,0]*m[1,1]-m[0,1]*m[1,0])/d}
};}

// ===== 2×2 逆矩阵 =====
double[,] A = {{4,7},{2,6}};
double detA = A[0,0]*A[1,1]-A[0,1]*A[1,0];
Console.WriteLine($"行列式: {detA}");
double[,] Ai = {{A[1,1]/detA,-A[0,1]/detA},{-A[1,0]/detA,A[0,0]/detA}};
Console.WriteLine($"逆矩阵:\n  [{Ai[0,0]:F1}, {Ai[0,1]:F1}]\n  [{Ai[1,0]:F1}, {Ai[1,1]:F1}]");
Console.Write("验证 A@A⁻¹:\n");
for (int i = 0; i < 2; i++){
    for (int j = 0; j < 2; j++){ double v=0; for(int k=0;k<2;k++) v+=A[i,k]*Ai[k,j]; Console.Write($"  {v:F1}"); }
    Console.WriteLine();
}

// ===== 正规方程求线性回归 =====
double[,] X = {{1,80,2},{1,120,3},{1,150,4}};
double[] y = {200,350,450};
double[,] XtX = new double[3,3]; double[] Xty = new double[3];
for(int i=0;i<3;i++) for(int j=0;j<3;j++) for(int k=0;k<3;k++) XtX[i,j]+=X[k,i]*X[k,j];
for(int i=0;i<3;i++) for(int k=0;k<3;k++) Xty[i]+=X[k,i]*y[k];
var inv = Inv3(XtX);
double[] w = new double[3];
for(int i=0;i<3;i++) for(int j=0;j<3;j++) w[i]+=inv[i,j]*Xty[j];
Console.Write("回归系数: "); foreach(var v in w) Console.Write($"{v:F2} ");
Console.Write("\n预测值:   ");
for(int i=0;i<3;i++){double p=0;for(int j=0;j<3;j++) p+=X[i,j]*w[j]; Console.Write($"{p:F1} ");}
Console.WriteLine();

// ===== 岭回归 =====
double lam=0.1;
var R=(double[,])XtX.Clone(); for(int i=0;i<3;i++) R[i,i]+=lam;
var ri=Inv3(R); double[] wr=new double[3];
for(int i=0;i<3;i++) for(int j=0;j<3;j++) wr[i]+=ri[i,j]*Xty[j];
Console.Write("岭回归系数: "); foreach(var v in wr) Console.Write($"{v:F2} ");
Console.WriteLine();

// ===== Softmax Jacobian =====
double[] z = {2.0,1.0,0.5};
double sumE=0; foreach(var v in z) sumE+=Math.Exp(v);
double[] sig = z.Select(v => Math.Exp(v)/sumE).ToArray();
Console.Write("Softmax: "); foreach(var s in sig) Console.Write($"{s:F4} ");
Console.Write("\nJacobian:\n");
for(int i=0;i<3;i++){
    for(int j=0;j<3;j++) Console.Write($"{(i==j?sig[i]:0)-sig[i]*sig[j],8:F4}");
    Console.WriteLine();
}
Console.Write("每行之和: ");
for(int i=0;i<3;i++){double s=0;for(int j=0;j<3;j++) s+=(i==j?sig[i]:0)-sig[i]*sig[j]; Console.Write($"{s:F4} ");}
Console.WriteLine();

// ===== 全连接层梯度 =====
double[,] W = {{0.15,0.20,-0.30},{-0.40,0.50,0.10}};
double[] x = {0.8,0.2,0.5};
double[] dLdz = {0.3,-0.5};
Console.Write("dL/dW:\n");
for(int i=0;i<2;i++){for(int j=0;j<3;j++) Console.Write($"  {dLdz[i]*x[j]:F2}");Console.WriteLine();}
Console.Write("dL/dx: ");
for(int j=0;j<3;j++){double v=0;for(int i=0;i<2;i++) v+=W[i,j]*dLdz[i]; Console.Write($"{v:F3} ");}
Console.WriteLine();

Rust

// 3×3 行列式
fn det3(m: &[[f64;3];3]) -> f64 {
    m[0][0]*(m[1][1]*m[2][2]-m[1][2]*m[2][1])
   -m[0][1]*(m[1][0]*m[2][2]-m[1][2]*m[2][0])
   +m[0][2]*(m[1][0]*m[2][1]-m[1][1]*m[2][0])
}
// 3×3 逆矩阵
fn inv3(m: &[[f64;3];3]) -> [[f64;3];3] {
    let d = det3(m);
    [[(m[1][1]*m[2][2]-m[1][2]*m[2][1])/d,(m[0][2]*m[2][1]-m[0][1]*m[2][2])/d,(m[0][1]*m[1][2]-m[0][2]*m[1][1])/d],
     [(m[1][2]*m[2][0]-m[1][0]*m[2][2])/d,(m[0][0]*m[2][2]-m[0][2]*m[2][0])/d,(m[0][2]*m[1][0]-m[0][0]*m[1][2])/d],
     [(m[1][0]*m[2][1]-m[1][1]*m[2][0])/d,(m[0][1]*m[2][0]-m[0][0]*m[2][1])/d,(m[0][0]*m[1][1]-m[0][1]*m[1][0])/d]]
}

fn main() {
    // ===== 2×2 逆矩阵 =====
    let a = [[4.0,7.0],[2.0,6.0_f64]];
    let det = a[0][0].mul_add(a[1][1], -(a[0][1]*a[1][0]));
    println!("行列式: {det}");
    let ai = [[a[1][1]/det,-a[0][1]/det],[-a[1][0]/det,a[0][0]/det]];
    println!("逆矩阵:\n  [{:.1}, {:.1}]\n  [{:.1}, {:.1}]", ai[0][0],ai[0][1],ai[1][0],ai[1][1]);
    println!("验证 A@A⁻¹:");
    for i in 0..2 {
        for j in 0..2 {
            let v = a[i][0].mul_add(ai[0][j], a[i][1]*ai[1][j]);
            print!("  {v:.1}");
        }
        println!();
    }

    // ===== 正规方程 =====
    let x = [[1.0,80.0,2.0],[1.0,120.0,3.0],[1.0,150.0,4.0_f64]];
    let y = [200.0,350.0,450.0_f64];
    let mut xtx = [[0.0_f64;3];3];
    let mut xty = [0.0_f64;3];
    for i in 0..3 { for j in 0..3 { for k in 0..3 { xtx[i][j]+=x[k][i]*x[k][j]; } } }
    for i in 0..3 { for k in 0..3 { xty[i]+=x[k][i]*y[k]; } }
    let inv = inv3(&xtx);
    let mut w = [0.0_f64;3];
    for i in 0..3 { for j in 0..3 { w[i]+=inv[i][j]*xty[j]; } }
    print!("回归系数: "); for v in w { print!("{v:.2} "); }
    print!("\n预测值:   ");
    for i in 0..3 { let p: f64=(0..3).map(|j|x[i][j]*w[j]).sum(); print!("{p:.1} "); }
    println!();

    // ===== 岭回归 =====
    let mut r = xtx;
    for i in 0..3 { r[i][i] += 0.1; }
    let ri = inv3(&r);
    let mut wr = [0.0_f64;3];
    for i in 0..3 { for j in 0..3 { wr[i]+=ri[i][j]*xty[j]; } }
    print!("岭回归系数: "); for v in wr { print!("{v:.2} "); }
    println!();

    // ===== Softmax Jacobian =====
    let z = [2.0, 1.0, 0.5_f64];
    let sum_e: f64 = z.iter().map(|v| v.exp()).sum();
    let sig: [f64;3] = core::array::from_fn(|i| z[i].exp()/sum_e);
    print!("Softmax: "); for s in sig { print!("{s:.4} "); }
    println!("\nJacobian:");
    for i in 0..3 {
        for j in 0..3 {
            let jac = if i==j { sig[i] } else { 0.0 } - sig[i]*sig[j];
            print!("{jac:8.4}");
        }
        println!();
    }
    print!("每行之和: ");
    for i in 0..3 {
        let s: f64=(0..3).map(|j|(if i==j{sig[i]}else{0.0})-sig[i]*sig[j]).sum();
        print!("{s:.4} ");
    }
    println!();

    // ===== 全连接层梯度 =====
    let ww = [[0.15,0.20,-0.30_f64],[-0.40,0.50,0.10]];
    let xv = [0.8, 0.2, 0.5_f64];
    let dl = [0.3, -0.5_f64];
    println!("dL/dW:");
    for i in 0..2 { for j in 0..3 { print!("  {:.2}",dl[i]*xv[j]); } println!(); }
    print!("dL/dx: ");
    for j in 0..3 { let v=ww[0][j]*dl[0]+ww[1][j]*dl[1]; print!("{v:.3} "); }
    println!();
}