多元微分、矩阵微积分与自动微分

标量训练目标依赖许多参数，常经过重复计算和记录批次。导数需包括所有依赖，区分偏向观察与整个输入的近似，并保留数组形状。本课先发展全微分、雅可比、Hessian与矩阵梯度，再把链规则实现为自动微分。

检索检查：回忆模块10矩阵积、转置、批形状和模块16局部模型、链规则、不可导边界。NumPy实验使用本课固定版本，标量反向引擎只需标准库。

开域上f:Rⁿ→R对第i坐标偏导，是[f(a+h e_i)−f(a)]/h的极限，其余坐标固定。存在偏导只检查一条坐标线，不自动认证同时变化时连续或全部方向一阶近似。对数和分母等定义域仍可能排除邻近路径。

D_v f(a)=lim_{t→0}[f(a+tv)−f(a)]/t为向量v方向导数。单位v描述每单位距离变化率，非单位还缩放沿路径速度。此处t双侧趋零；边界或非光滑的单侧方向导数需另限定。存在这些极限不保证它们关于v构成一个线性映射。

等值集由输出相同输入组成。在全可导时，梯度指欧氏意义下一阶最陡上升，对可微等值路径非零切向正交。这些由微分与链规则证明，不只凭等高图。不同特征单位会改变欧氏几何，应声明坐标、内积或度量约定。

偏导能漏不连续：原点外f=xy/(x²+y²)，原点零。两坐标限制恒零，偏导均零；沿x=y=t≠0值1/2，故不连续、更不可全可导。实验采样仅说明已经精确推导的反例。

甚至连续且所有方向导数存在仍不足。原点外g=x³/(x²+y²)、原点零，幅度≤|x|≤‖(x,y)‖，所以连续。非零v方向速率v₁³/(v₁²+v₂²)，轴方向为一、零，(1,1)却为1/2而非线性候选一。方向率映射不线性，无统一全导数。检验多条直线不替代全部小向量的统一余项条件。

全可导要求线性泛函L，使f(a+h)=f(a)+L(h)+r(h)，r(h)/‖h‖→0，覆盖整个向量的全部接近路径。有限维欧氏线性泛函可写gᵀh：展开坐标基，g_i=L(e_i)，线性给表示。沿坐标限制余项则说明这些系数等于偏导。

本课用n×1列梯度∇f，所以df=∇fᵀdx是标量。输入到标量的行导数为1×n，是其转置。NumPy(n,)可存系数，却不是显式行列；应同时标数学约定与存储形状，避免方形例把转置错误隐藏。

F:Rⁿ→Rᵐ全可导指F(a+h)=F(a)+Jh+r(h)，‖r‖/‖h‖→0。J为m×n，行i是输出i对全部输入偏导，J_{ij}=∂F_i/∂x_j，因此dF=Jdx。有限维逐分量全可导等价，因为有限组消失余项的向量范数也消失。雅可比是输入到输出的映射，不是可任意摆向的坐标表。

邻域内连续偏导是充分条件。逐坐标改变输入并望远镜分解增量，对每一步用一维中值定理，系数为趋a的中间点偏导。减a处偏导，连续性将系数误差统一压为ε，余项≤εΣ|h_i|≤ε√n‖h‖，得到全可导。这证明C¹充分，而非必要定义。

全可导蕴含D_vf=∇fᵀv及连续：有界线性图和余项都趋零。单位v下柯西–施瓦茨给最大速率‖∇f‖，非零梯度方向取等；梯度零只表示一阶率都零。可微等值路径c(t)用链规则给∇f(c(t))ᵀc′(t)=0，得到法向解释。

微分不是Python存储的字面无穷小；用于有限扰动只是近似，余项另界定。改变坐标或度量会改变代表同一微分的梯度；正定加权内积M下梯度为M⁻¹乘欧氏梯度，标量作用不变，最陡方向依长度定义。

点上可导未必邻近C¹：非零x时s=x²sin(1/x)、零处零，商x sin(1/x)趋零，所以s′(0)=0；非零导数2x sin(1/x)−cos(1/x)振荡。f(x,y)=s(x)+y²幅度≤x²+y²，原点零线性余项除范数趋零，因此全可导，但导数不连续。把C¹充分条件当定义会删合法点可导函数。

可导F:Rⁿ→Rᵐ、G:Rᵐ→Rᵖ且域兼容，复合J=J_G(F(a))J_F(a)，形状(p×m)(m×n)=p×n。内先作用所以右乘，反序可能无定义或描述其他映射。

证明写ΔF=J_Fh+r_F，余项比范数趋零，于是‖ΔF‖/‖h‖有界。外增量J_GΔF+r_G(ΔF)，第一余J_G r_F消失，第二在ΔF非零时其比例拆为[‖r_G‖/‖ΔF‖][‖ΔF‖/‖h‖]，前趋零后有界；ΔF零时直接零。避免除消失内增量。

JVP是Jv，输入切向长度n映到输出m。VJP传统写输出行种子qᵀJ，本列约定转置为Jᵀq，长度n。不能误用J乘输出种子；种子规定微分哪个输出组合，向量输出没有未指定组合的单一梯度。

矩形例更暴露错误：F=(xy,x²,y²)，L=u+2v−3w，(1,2)处J行(2,1),(2,0),(0,4)，3×2。q=(1,2,−3)，Jᵀq=(6,−11)，与xy+2x²−3y²一致，Jq内维不匹配。长度二输入方向则Jv产生三个率。

前向逐局部JVP传输入切向，反向逐转置作用传输出伴随。标量输出多输入时一反向给种子一的全部输入系数。完整m×n雅可比通常要n次基输入前向或m次基输出反向，不自动以一次标量路径成本得到全部条目。

L(y)标量、y=Ax时dy=A dx、dL=qᵀA dx，拉回Aᵀq。非线性局部导数须在真实前向值求；替换输入导数或反向前修改缓存会失效，需要存储或重算中间值。

共享输入影响多路径，最终标量系数是影响之和。反向需等下游贡献全部到达共享节点再运行局部规则，唯一节点拓扑顺序给此日程；把图当树、无累积递归可能丢路径或重复传播已合计量。

Hessian是列梯度的雅可比，n×n，H_{ij}为梯度分量i对坐标j导数。邻域二阶偏导连续时混合可交换，H对称：考虑小坐标矩形四角增量，按两种顺序用中值定理，除面积、宽度趋零，连续使极限相同。不能仅因写了两种偏导符号就假定相等。

无连续性可不同：原点外f=xy(x²−y²)/(x²+y²)、原点零，f_x(0,y)=−y，f_y(x,0)=x，因此原点∂_y f_x=−1、∂_x f_y=1，邻域光滑假设失败。

C²下二阶模型f(a+h)=f(a)+∇f(a)ᵀh+hᵀH(a)h/2+o(‖h‖²)。沿g(t)=f(a+th)，g′=∇fᵀh、g″=hᵀHh。一维泰勒给线段中间H，连续使其与H(a)差在小h时统一消失，得到余项比例。仅列存在的坐标二阶项弱于梯度全可导；反例没有C²二阶模型。数值估计轻微非对称也可能是差分或程序错，应先判数学正则，再判估计器。

内部驻点梯度须零，可逐坐标用费马。C²下H正定严格局部最小、负定最大、不定给正负方向而鞍点；半正定有零方向则无结论。零处x⁴+y⁴最小，x⁴−y⁴鞍点，H都零，局部与后面的全局优化不同。

对称H下vᵀHv是方向二阶率，Hv为梯度沿v的一阶变化，长度n。可计算乘积而不存n²条目，但需足够光滑，折点经典H可不存在。梯度差分是另有步长误差的估计。

正定局部判别需统一二次裕量。单位球面vᵀHv连续正，紧致给最小λ>0，故hᵀHh≥λ‖h‖²。足够小邻域余项幅度≤λ‖h‖²/4，驻点非零增量≥λ‖h‖²/4，完成严格最小证明。采样少数正曲率不能替代该全方向裕量。

向量df=gᵀdx；矩阵用弗罗贝尼乌斯内积df=Σ(∇A f){ij}dA_{ij}=tr((∇A f)ᵀdA)，梯度与A同形状，每条目一个系数。tr为方阵对角和，兼容矩形B,C下tr(BC)=tr(CB)由同一个ΣB{ij}C_{ji}得，循环可用来收集合法维数微分。

固定b的bᵀx梯度b。一般方阵二次型df=dxᵀAx+xᵀA dx=[(A+Aᵀ)x]ᵀdx，梯度(A+Aᵀ)x，只有对称A才2Ax。反对称部分标量贡献零。输出形状适合不能修复错误简化。

归约是数学定义。单输出均方梯度(2/B)Xᵀr；多输出B×o所有条目平均除B·o，而先每行输出和再行平均只除B，目标差因子o。需准确报告，否则改变梯度与更新有效尺度。

X:B×d、W:d×o、Y:B×o，R=XW−Y，半平方Frobenius损失微分tr(RᵀdR)，dR=X dW+dX W−dY，收系数得W梯度XᵀR:d×o，X梯度R Wᵀ:B×d，Y梯度−R:B×o。先通过索引或微分确认，不混平目标与列预测造成矩阵广播。

共享数学列偏置b:o×1，预测XW+1bᵀ，db在B行重复，偏置梯度Rᵀ1，即行方向求和。NumPy常存(o,)，sum(axis=0)恢复同形状。前向广播对应反向求和，因为一个参数影响所有复制项；保留最后一行会丢其他贡献。

两层Z=XW₁+1b₁ᵀ、H=tanh Z、P=HW₂+1b₂ᵀ，输出伴随G。W₂梯度HᵀG，b₂行求和，隐藏伴随G W₂ᵀ，逐元素乘1−H²得D，W₁梯度XᵀD、b₁行求和。矩阵乘与逐元素乘作用不同。

实验B=3、输入二、隐藏三、输出一：W₁2×3，Z,H3×3，W₂3×1，P,G3×1，G W₂ᵀ3×3，XᵀD2×3。偏置存(3,)、(1,)。缺转置特意会失败，不用方形阵掩盖。均方三标量所以G=2R/3，形状对仍可能归约错。

AD按执行计算图使用局部导数规则，不是扰动取商的有限差分，也不是构造代数表达的符号微分。前向带值与切向，反向先记前向值再逆依赖传播伴随。导数规则按数学原语计算，求值仍有数值舍入。

字面浮点输入输出图有离散舍入阶梯，AD通常微分在存储值处求的目标实原语公式，而非每条舍入指令阶梯。“自动”不等于精确实算术，不消除条件敏感性，不证明跨分支可导。分支结果描述实际路径与原语约定，零处ReLU反向仍是约定。

实验引擎存父节点与局部闭包，先建每可达节点一次的拓扑列表，新反向清零所有伴随、设种子、逆序遍历。每父边用+=，两操作数同对象也加两次，确保共享节点已收到下游总和。无重置重复反向会累积旧请求，只有明确有意时才合理。

种子不是一则给缩放VJP，x=2种子二返回74、零返回零，普通导数仍37。多输出向量种子指定加权组合，梯度回答定义好的微分量，不是变量永久属性。状态、图重用、归约均在规格内。

梯度检查在选定光滑点、多个合适步长比较解析AD与中心差分，先检查形状，再按幅度设置绝对相对容差。错误公式可在对称或零残差点巧合，需非平凡值和逐参数或独立方向扰动。网络实验检验全部小参数并展示转置、丢批贡献、均值因子故障，支持这些例的实现，通式依据微分推导。

反向需存储或重算中间值及正确图日程，成本依计算而不只参数数。大系统检查点、Hv可改变工程代价，链与累积规格仍同。这个基础足够解释手动和自动反向，再进入优化。

实验1 · 向量导数与偏导反例

下载 lab1_vector_derivatives.py

"""Finite differences check stated analytic vector derivatives at smooth points."""
import numpy as np
np.set_printoptions(precision=6, suppress=True)

def F(z):
    x, y = z
    return np.array([x*x*y, np.sin(x)+y*y])

point = np.array([0.7, -0.4])
x, y = point
J = np.array([[2*x*y, x*x], [np.cos(x), 2*y]])
h = 1e-5
basis = np.eye(2)
numerical_J = np.column_stack([(F(point+h*e)-F(point-h*e))/(2*h) for e in basis])
print("Point / function value:", point, F(point))
print("Analytic Jacobian (outputs by inputs):\n", J)
print("Finite-difference Jacobian:\n", numerical_J)
print("Jacobian maximum error:", np.max(np.abs(J-numerical_J)))
assert J.shape == (2, 2) and np.allclose(J, numerical_J, atol=1e-9, rtol=0)

direction = np.array([0.6, 0.8])
directional = (F(point+h*direction)-F(point-h*direction))/(2*h)
print("JVP / directional difference:", J@direction, directional)
assert np.allclose(J@direction, directional, atol=1e-9, rtol=0)

def loss(z):
    return z[0]**2*z[1]+np.exp(z[1])
def grad(z):
    return np.array([2*z[0]*z[1], z[0]**2+np.exp(z[1])])
H = np.array([[2*y, 2*x], [2*x, np.exp(y)]])
numerical_grad = np.array([(loss(point+h*e)-loss(point-h*e))/(2*h) for e in basis])
numerical_H = np.column_stack([(grad(point+h*e)-grad(point-h*e))/(2*h) for e in basis])
print("Scalar gradient:", grad(point))
print("Hessian:\n", H)
print("Gradient / Hessian maximum errors:", np.max(np.abs(grad(point)-numerical_grad)), np.max(np.abs(H-numerical_H)))
assert np.allclose(grad(point), numerical_grad, atol=1e-9, rtol=0)
assert np.allclose(H, numerical_H, atol=1e-9, rtol=0)

# Axis partials alone do not imply continuity, hence do not imply differentiability.
def bad(z):
    a, b = z
    return a*b/(a*a+b*b) if a*a+b*b else 0.0
print("Origin axis samples / diagonal samples:")
for t in [1e-1, 1e-3, 1e-6]:
    print(t, bad([t, 0]), bad([0, t]), bad([t, t]))
    assert bad([t, 0]) == 0 and bad([0, t]) == 0 and bad([t, t]) == 0.5
print("Smooth derivative checks illustrate the analytic claims; axis samples cannot certify total differentiability.")

Point / function value: [ 0.7 -0.4] [-0.196     0.804218]
Analytic Jacobian (outputs by inputs):
 [[-0.56      0.49    ]
 [ 0.764842 -0.8     ]]
Finite-difference Jacobian:
 [[-0.56      0.49    ]
 [ 0.764842 -0.8     ]]
Jacobian maximum error: 1.702149532434305e-11
JVP / directional difference: [ 0.056    -0.181095] [ 0.056    -0.181095]
Scalar gradient: [-0.56     1.16032]
Hessian:
 [[-0.8      1.4    ]
 [ 1.4      0.67032]]
Gradient / Hessian maximum errors: 1.24074084340009e-11 1.5253354135325026e-11
Origin axis samples / diagonal samples:
0.1 0.0 0.0 0.5
0.001 0.0 0.0 0.5
1e-06 0.0 0.0 0.5
Smooth derivative checks illustrate the analytic claims; axis samples cannot certify total differentiability.

先推导J、梯度、H，再比较形状、误差和光滑假设。最后轴对角采样说明精确反例，不是用采样证明连续。

实验2 · 最小反向引擎

下载 lab2_reverse_engine.py

"""A tiny scalar reverse engine: unique-node order, reset, seed, and += adjoints."""
import math

class Value:
    def __init__(self, data, parents=(), label=""):
        self.data, self.grad = float(data), 0.0
        self.parents, self.label = parents, label
        self._backward = lambda: None

    def __add__(self, other):
        other = other if isinstance(other, Value) else Value(other)
        out = Value(self.data+other.data, (self, other), "+")
        def backward():
            self.grad += out.grad
            other.grad += out.grad
        out._backward = backward
        return out
    __radd__ = __add__

    def __mul__(self, other):
        other = other if isinstance(other, Value) else Value(other)
        out = Value(self.data*other.data, (self, other), "*")
        def backward():
            self.grad += other.data*out.grad
            other.grad += self.data*out.grad
        out._backward = backward
        return out
    __rmul__ = __mul__

    def tanh(self):
        out = Value(math.tanh(self.data), (self,), "tanh")
        def backward():
            self.grad += (1-out.data*out.data)*out.grad
        out._backward = backward
        return out

    def backward(self, seed=1.0):
        visited, order = set(), []
        def visit(node):
            if node in visited:
                return
            visited.add(node)
            for parent in node.parents:
                visit(parent)
            order.append(node)
        visit(self)
        for node in order:
            node.grad = 0.0
        self.grad = float(seed)
        for node in reversed(order):
            node._backward()
        return order

x = Value(2, label="x")
u = x*x
u.label = "u=x*x"
square = u*u
square.label = "u*u"
z = square+u+x
z.label = "z"
order = z.backward()
print("Shared graph z=(x*x)^2+x*x+x at x=2")
for node in order:
    print(node.label, "value", node.data, "adjoint", node.grad)
assert z.data == 22 and x.grad == 37 and u.grad == 9
z.backward(seed=2)
print("Seed two input adjoint:", x.grad)
assert x.grad == 74
z.backward()
print("Fresh seeded pass resets rather than doubles:", x.grad)
assert x.grad == 37

a, b = Value(2), Value(-1)
mixed = a*b+a*a
mixed.backward()
print("Two-input value / gradient:", mixed.data, a.grad, b.grad)
assert mixed.data == 2 and (a.grad, b.grad) == (3, 2)

v = Value(0.4)
nonlinear = (v*v+0.5*v).tanh()
nonlinear.backward()
def reference(t):
    return math.tanh(t*t+0.5*t)
h = 1e-6
finite = (reference(0.4+h)-reference(0.4-h))/(2*h)
print("Nonlinear AD / finite difference:", v.grad, finite)
assert abs(v.grad-finite) < 1e-8
print("Overwriting a shared adjoint would discard a path; addition is part of the chain rule.")

Shared graph z=(x*x)^2+x*x+x at x=2
x value 2.0 adjoint 37.0
u=x*x value 4.0 adjoint 9.0
u*u value 16.0 adjoint 1.0
+ value 20.0 adjoint 1.0
z value 22.0 adjoint 1.0
Seed two input adjoint: 74.0
Fresh seeded pass resets rather than doubles: 37.0
Two-input value / gradient: 2.0 3.0 2.0
Nonlinear AD / finite difference: 1.1450754518266433 1.145075451791655
Overwriting a shared adjoint would discard a path; addition is part of the chain rule.

追拓扑唯一节点与父边贡献的区别，解释同对象x*x仍加两次。先预测种子二和重置，再比较光滑非线性差分。

实验3 · 全部两层参数梯度

下载 lab3_network_backward.py

"""Check every two-layer tanh-network parameter and expose shape/scaling faults."""
import numpy as np
np.set_printoptions(precision=6, suppress=True)
X = np.array([[0.2,-0.3], [0.5,0.7], [-0.4,0.1]])
Y = np.array([[0.1], [0.8], [-0.2]])
parameters = {
    "W1": np.array([[0.3,-0.2,0.1], [0.4,0.2,-0.5]]),
    "b1": np.array([0.1,-0.1,0.05]),
    "W2": np.array([[0.5], [-0.3], [0.2]]),
    "b2": np.array([0.02]),
}

def forward(p):
    H = np.tanh(X@p["W1"]+p["b1"])
    prediction = H@p["W2"]+p["b2"]
    residual = prediction-Y
    return np.mean(residual*residual), H, prediction, residual

loss, H, prediction, residual = forward(parameters)
# Mean of B*o entries, here B=3 and o=1; the exact denominator is explicit.
G = 2*residual/residual.size
D = (G@parameters["W2"].T)*(1-H*H)
gradients = {"W2": H.T@G, "b2": G.sum(axis=0), "W1": X.T@D, "b1": D.sum(axis=0)}
print("Loss / prediction:", loss, prediction.ravel())
print("Batch / input / hidden / output shapes:", X.shape, parameters["W1"].shape, H.shape, prediction.shape)
h = 1e-6
for name, array in parameters.items():
    numerical = np.zeros_like(array)
    for index in np.ndindex(array.shape):
        saved = array[index]
        array[index] = saved+h
        plus = forward(parameters)[0]
        array[index] = saved-h
        minus = forward(parameters)[0]
        array[index] = saved
        numerical[index] = (plus-minus)/(2*h)
    error = np.max(np.abs(numerical-gradients[name]))
    print(name, "shape", gradients[name].shape, "maximum derivative error", error)
    assert gradients[name].shape == array.shape and error < 1e-9

try:
    G@parameters["W2"]  # Missing transpose: (3,1) cannot multiply (3,1).
except ValueError:
    print("Missing transpose rejected by incompatible derivative shapes.")
wrong_bias = D[-1]  # Overwrite rather than accumulate contributions from all rows.
wrong_sum_scaling = 2*residual  # This differentiates a sum, not the stated mean.
print("Bias overwrite maximum error:", np.max(np.abs(wrong_bias-gradients["b1"])))
print("Missing mean factor scales output adjoint by:", residual.size)
assert not np.allclose(wrong_bias, gradients["b1"])
assert np.allclose(wrong_sum_scaling, residual.size*G)
print("Shape, shared contributions, and loss reduction are separate parts of the backward contract.")

Loss / prediction: 0.12469122832324815 [0.142506 0.231685 0.022003]
Batch / input / hidden / output shapes: (3, 2) (2, 3) (3, 3) (3, 1)
W1 shape (2, 3) maximum derivative error 2.7340560371236222e-11
b1 shape (3,) maximum derivative error 8.312184274217316e-12
W2 shape (3, 1) maximum derivative error 4.48385772955362e-12
b2 shape (1,) maximum derivative error 5.846850781310309e-12
Missing transpose rejected by incompatible derivative shapes.
Bias overwrite maximum error: 0.13066171364054682
Missing mean factor scales output adjoint by: 3
Shape, shared contributions, and loss reduction are separate parts of the backward contract.

标前向反向形状、精确均方定义，逐条目检查。解释转置拒绝、覆盖偏置、丢均值因子为何各违反规格，即使数组看似合理。

1–12必做，13–14选做额外35分钟，不计入十二小时。

用作者Mathematics for Machine Learning伴随站与一手MIT OpenCourseWare多元微积分梯度、J、H与链材料。本课微分推导和引擎给具体形状与累积约定。

给偏导反例和C¹证明，写列梯度与输出乘输入J，推链拉回、二阶模型、二次型与最小二乘微分。

结业任务：X:B×d、W:d×o，半平方R=XW−Y，推XᵀR、RWᵀ及偏置行和，延伸tanh层，写逐元素因子和精确均值。解释唯一节点日程仍累加x*x两边。

可以继续：你能推导验证向量矩阵敏感性、解释自动反向。下一块优化研究导数对极小值与更新收敛能证明什么。见总览。