本讲导读 本讲（Toda 第 5 讲）介绍求解资产定价/动态经济问题所需的数值方法。§1 多项式插值：Lagrange 插值（Prop 1）、误差界（Prop 2）、Chebyshev 多项式（Thm 3：极小化 \(\max|\prod(x-x_k)|\) 的最优节点 \(x_k=\cos\frac{2k-1}{2n}\pi\)）。§2 数值积分 (quadrature)：Newton-Cotes（梯形法、Simpson 法）、复合法则、高斯求积（正交多项式 Prop 4、根 Lemma 5、Thm 6：\(n\) 个根作节点可精确积分至 \(2n-1\) 阶、Golub-Welsch Thm 7、Gauss-Legendre/Chebyshev/Hermite/Laguerre）。§3 离散化（把连续随机过程近似为有限状态马尔可夫链）：3.1 早期方法（Tauchen 1986 不准、Rouwenhorst 1995、Tauchen-Hussey 1991）；3.2 Farmer-Tanaka-Toda 最大熵法（KL 相对熵最小化精确匹配矩，原始 (P)/对偶 (D)、Thm 8、Markov 推广、Algorithm 1）。含每讲参考文献。

命题 1（Lagrange 插值）/ Proposition 1 (Lagrange interpolation) 多项式易于微分积分。\(n-1\) 次多项式由 \(n\) 个系数确定，故给定 \(n\) 个点至多有一个多项式经过它们。设 \(x_1<\dotsLagrange 多项式 \(L_k(x)=\dfrac{\prod_{l\neq k}(x-x_l)}{\prod_{l\neq k}(x_k-x_l)}\)，则Polynomials are easy to differentiate and integrate. A degree \(n-1\) polynomial is determined by \(n\) coefficients, so given \(n\) points there is at most one polynomial through them. Let \(x_1<\dotsLagrange polynomial \(L_k(x)=\dfrac{\prod_{l\neq k}(x-x_l)}{\prod_{l\neq k}(x_k-x_l)}\); then

$$p(x)=\sum_{k=1}^n y_k L_k(x)$$

是满足 \(p(x_k)=y_k\) 的唯一次数 \(\le n-1\) 的多项式。（证明：\(L_k(x_l)=\delta_{kl}\)，故 \(p(x_l)=\sum_k y_k\delta_{kl}=y_l\)。）is the unique polynomial of degree up to \(n-1\) with \(p(x_k)=y_k\). (Proof: \(L_k(x_l)=\delta_{kl}\), so \(p(x_l)=\sum_k y_k\delta_{kl}=y_l\).)

命题 2（插值误差界）/ Proposition 2 (interpolation error bound) 设 \(f\) 是 \(C^n\)、\(p_{n-1}\) 是 \(f\) 在 \(x_1,\dots,x_n\) 处的插值多项式。则对任意 \(x\)，存在 \(\xi\in\mathrm{co}\{x,x_1,\dots,x_n\}\) 使Let \(f\) be \(C^n\) and \(p_{n-1}\) the interpolating polynomial of \(f\) at \(x_1,\dots,x_n\). Then for any \(x\) there is \(\xi\in\mathrm{co}\{x,x_1,\dots,x_n\}\) with

$$f(x)-p_{n-1}(x)=\frac{f^{(n)}(\xi)}{n!}\prod_{k=1}^n(x-x_k).$$

命题 2 证明（Rolle 定理）/ Proof of Proposition 2 (Rolle's theorem) 记 \(I=\mathrm{co}\{x,x_1,\dots,x_n\}\)，\(R(t)=f(t)-p_{n-1}(t)\)，\(S(t)=\prod_{k=1}^n(t-x_k)\)。令 \(g(t)=R(t)S(x)-R(x)S(t)\)，则 \(g(x)=0\) 且 \(g(x_k)=0\)（因 \(R(x_k)=S(x_k)=0\)）。由 Rolle 定理反复应用，存在 \(\xi\in I\) 使 \(g^{(n)}(\xi)=0\)。因 \(S\) 是首项系数 1 的 \(n\) 次多项式，\(S^{(n)}=n!\)，故 \(0=g^{(n)}(\xi)=R^{(n)}(\xi)S(x)-R(x)n!\)。又 \(\deg p_{n-1}\le n-1\) 故 \(R^{(n)}(\xi)=f^{(n)}(\xi)\)，得 \(f(x)-p_{n-1}(x)=\frac{f^{(n)}(\xi)}{n!}\prod_k(x-x_k)\)。\(\blacksquare\)Let \(I=\mathrm{co}\{x,x_1,\dots,x_n\}\), \(R(t)=f(t)-p_{n-1}(t)\), \(S(t)=\prod_{k=1}^n(t-x_k)\). Set \(g(t)=R(t)S(x)-R(x)S(t)\); then \(g(x)=0\) and \(g(x_k)=0\) (since \(R(x_k)=S(x_k)=0\)). By repeated application of Rolle's theorem, there is \(\xi\in I\) with \(g^{(n)}(\xi)=0\). As \(S\) is a degree-\(n\) polynomial with leading coefficient 1, \(S^{(n)}=n!\), so \(0=g^{(n)}(\xi)=R^{(n)}(\xi)S(x)-R(x)n!\). Since \(\deg p_{n-1}\le n-1\), \(R^{(n)}(\xi)=f^{(n)}(\xi)\), giving \(f(x)-p_{n-1}(x)=\frac{f^{(n)}(\xi)}{n!}\prod_k(x-x_k)\). \(\blacksquare\)

1.2 Chebyshev 多项式（定理 3）/ 1.2 Chebyshev polynomials (Theorem 3) 在区间（不失一般性取 $[-1,1]$）上插值时如何选节点？因 \(f^{(n)}(\xi)\) 依赖 \(f\) 而 \(\prod(x-x_k)\) 不依赖，应选节点极小化 \(\max_{x\in[-1,1]}|\prod_{k=1}^n(x-x_k)|\)。Chebyshev 多项式 \(T_n(x)\)：把 \(\cos n\theta\) 展为 \(\cos\theta\) 的 \(n\) 次多项式并令 \(x=\cos\theta\)（\(T_0=1,T_1=x,T_2=2x^2-1\)，递推 \(T_{n+1}(x)=2xT_n(x)-T_{n-1}(x)\)）。定理 3：When interpolating on an interval (WLOG $[-1,1]\(), how to choose nodes? Since \)f^{(n)}(\xi)$ depends on \(f\) but \(\prod(x-x_k)\) does not, choose nodes to minimize \(\max_{x\in[-1,1]}|\prod_{k=1}^n(x-x_k)|\). The Chebyshev polynomial \(T_n(x)\): expand \(\cos n\theta\) as a degree-\(n\) polynomial of \(\cos\theta\) and set \(x=\cos\theta\) (\(T_0=1,T_1=x,T_2=2x^2-1\), recursion \(T_{n+1}(x)=2xT_n(x)-T_{n-1}(x)\)). Theorem 3:

$$\min_{x_1,\dots,x_n}\max_{x\in[-1,1]}\left|\prod_{k=1}^n(x-x_k)\right|\ \text{is attained at}\ x_k=\cos\frac{2k-1}{2n}\pi,\quad\text{where}\ \prod_{k=1}^n(x-x_k)=2^{1-n}T_n(x).$$

定理 3 证明（等振荡）/ Proof of Theorem 3 (equioscillation) 令 \(p(x)=2^{1-n}T_n(x)\)。由递推，\(T_n\) 首项系数 \(2^{n-1}\)，故 \(p\) 首项系数 1。\(p(\cos\theta)=2^{1-n}\cos n\theta\)，故 \(\sup_{[-1,1]}|p|=2^{1-n}\)。若存在首项系数 1 的 \(n\) 次多项式 \(q\) 使 \(\sup|q|<2^{1-n}\)：\(p\) 在 \(y_k=\cos(k\pi/n)\)（\(k=0,\dots,n\)）处取 \((-1)^k 2^{1-n}\)，由介值定理在 \(y_0,\dots,y_n\) 间有 \(z_1,\dots,z_n\) 使 \(p(z_k)-q(z_k)=0\)；但 \(r=p-q\) 次数 \(\le n-1\) 却有 \(n\) 个根，故 \(r\equiv0\) 即 \(p\equiv q\)，矛盾。故最优节点为 \(x_k=\cos\frac{2k-1}{2n}\pi\)。\(\blacksquare\)Let \(p(x)=2^{1-n}T_n(x)\). By the recursion, \(T_n\) has leading coefficient \(2^{n-1}\), so \(p\) has leading coefficient 1. \(p(\cos\theta)=2^{1-n}\cos n\theta\), so \(\sup_{[-1,1]}|p|=2^{1-n}\). Suppose a degree-\(n\) polynomial \(q\) with leading coefficient 1 had \(\sup|q|<2^{1-n}\): \(p\) takes \((-1)^k 2^{1-n}\) at \(y_k=\cos(k\pi/n)\) (\(k=0,\dots,n\)), so by the intermediate value theorem there are \(z_1,\dots,z_n\) between \(y_0,\dots,y_n\) with \(p(z_k)-q(z_k)=0\); but \(r=p-q\) has degree \(\le n-1\) yet \(n\) roots, so \(r\equiv0\), i.e. \(p\equiv q\), a contradiction. Hence the optimal nodes are \(x_k=\cos\frac{2k-1}{2n}\pi\). \(\blacksquare\)

Newton-Cotes：梯形法与 Simpson 法 / Newton-Cotes: trapezoidal and Simpson 许多经济问题涉及最大化期望（如 \(\max_\theta\mathbb E[u(R(\theta)w)]\)）。期望是积分，多数积分无显式解，需数值积分 (quadrature) \(\int_a^b f(x)\,dx\approx\sum_{n=1}^N w_n f(x_n)\)（\(\{x_n\}\) 节点、\(\{w_n\}\) 权重）。Newton-Cotes：取等距节点、选权重使能精确积分 \(N-1\) 次以下多项式。\(N=2\)（梯形法）：\(w_1=w_2=\frac12\)，Many economic problems maximize an expectation (e.g. \(\max_\theta\mathbb E[u(R(\theta)w)]\)). Expectations are integrals, most without closed form, requiring quadrature \(\int_a^b f(x)\,dx\approx\sum_{n=1}^N w_n f(x_n)\) (\(\{x_n\}\) nodes, \(\{w_n\}\) weights). Newton-Cotes: take evenly spaced nodes, choose weights to integrate polynomials up to degree \(N-1\) exactly. \(N=2\) (trapezoidal rule): \(w_1=w_2=\frac12\),

$$\int_a^b f(x)\,dx\approx\frac{b-a}2(f(a)+f(b)).$$

\(N=3\)（Simpson 法）：\(w_1=w_3=\frac16,w_2=\frac23\)，\(N=3\) (Simpson's rule): \(w_1=w_3=\frac16,w_2=\frac23\),

$$\int_a^b f(x)\,dx\approx\frac{b-a}6\left(f(a)+4f\!\left(\frac{a+b}2\right)+f(b)\right),$$

它实际能精确积分 3 次多项式，误差界 \(\left|\int_a^b f-\frac{b-a}6(\dots)\right|\le\frac{\|f^{(4)}\|}{2880}(b-a)^5\)。which actually integrates degree-3 polynomials exactly, with error bound \(\left|\int_a^b f-\frac{b-a}6(\dots)\right|\le\frac{\|f^{(4)}\|}{2880}(b-a)^5\).

复合法则 / Compound rule \(N\ge4\) 的 Newton-Cotes 几乎不用（高阶时部分权重变负、引入舍入误差）。复合法则把 \([a,b]\) 分成 \(N\) 个等距子区间，对每段用梯形/Simpson。\((N+1)\)-点梯形法误差 \(\frac{\|f''\|}{12}N^{-2}\)；\((2N+1)\)-点 Simpson 法误差 \(\frac{\|f^{(4)}\|}{2880}N^{-4}\)。被积函数越光滑、Simpson 比梯形越准。Newton-Cotes with \(N\ge4\) is almost never used (beyond some order some weights turn negative, introducing rounding errors). The compound rule divides \([a,b]\) into \(N\) even subintervals and applies trapezoidal/Simpson to each. The \((N+1)\)-point trapezoidal error is \(\frac{\|f''\|}{12}N^{-2}\); the \((2N+1)\)-point Simpson error is \(\frac{\|f^{(4)}\|}{2880}N^{-4}\). The smoother the integrand, the more accurate; Simpson beats trapezoidal.

高斯求积：正交多项式、根、定理 6 / Gaussian quadrature 设权函数 \(w\) 下的内积 \((p,q)=\int w(x)p(x)q(x)\,dx\) 与正交多项式 \(\{p_n^*\}\)。命题 4：正交多项式满足三项递推 (TTRR) \(p_{n+1}(x)=xp_n^*(x)-(xp_n^*,p_n^*)p_n^*(x)-(p_n,p_n)^{1/2}p_{n-1}^*(x)\)。引理 5：\(p_n^*(x)\) 在 \((a,b)\) 上恰有 \(n\) 个实根。定理 6（高斯求积）：取 \(p_n^*\) 的 \(n\) 个根 \(aLet the inner product \((p,q)=\int w(x)p(x)q(x)\,dx\) under weight \(w\) define orthogonal polynomials \(\{p_n^*\}\). Proposition 4: they satisfy a three-term recurrence (TTRR) \(p_{n+1}(x)=xp_n^*(x)-(xp_n^*,p_n^*)p_n^*(x)-(p_n,p_n)^{1/2}p_{n-1}^*(x)\). Lemma 5: \(p_n^*(x)\) has exactly \(n\) real roots on \((a,b)\). Theorem 6 (Gaussian quadrature): take the \(n\) roots \(a

$$\int w(x)p(x)\,dx=\sum_{k=1}^n w_k p(x_k).$$

（证明：写 \(p=p_n^* q+r\)，\(\deg q,r\le n-1\)，由正交性 \((p_n^*,q)=0\)，\(p(x_k)=r(x_k)\)，再对 \(r\) 用 \(L_k\) 表示。）(Proof: write \(p=p_n^* q+r\) with \(\deg q,r\le n-1\); by orthogonality \((p_n^*,q)=0\), \(p(x_k)=r(x_k)\), and represent \(r\) via the \(L_k\).)

定理 7（Golub-Welsch）与四类高斯求积 / Theorem 7 (Golub-Welsch) and four Gaussian quadratures 定理 7（Golub-Welsch）：设 \(p_n^*\) 首项系数 \(k_n>0\)、\(\alpha_n=k_{n-1}/k_n>0\)、\(\beta_n=(xp_n^*,p_n^*)\)，构造 \(N\times N\) 对称三对角 Jacobi 矩阵 \(J_N\)（对角 \(\beta_0,\dots,\beta_{N-1}\)、次对角 \(\alpha_1,\dots,\alpha_{N-1}\)）。则求积节点 \(\{x_n\}\) 是 \(J_N\) 的特征值，权重满足 \(\frac1{w_n}=\sum_{k=0}^{N-1}p_k^*(x_n)^2>0\)。常见四类：Gauss-Legendre（$(-1,1)$，\(w=1\)）；Gauss-Chebyshev（$(-1,1)$，\(w=1/\sqrt{1-x^2}\)，算 Fourier 系数）；Gauss-Hermite（\((-\infty,\infty)\)，\(w=e^{-x^2}\)，算正态期望）；Gauss-Laguerre（\((0,\infty)\)，\(w=e^{-x}\)，算指数分布期望）。Theorem 7 (Golub-Welsch): with \(p_n^*\) leading coefficient \(k_n>0\), \(\alpha_n=k_{n-1}/k_n>0\), \(\beta_n=(xp_n^*,p_n^*)\), build the \(N\times N\) symmetric tridiagonal Jacobi matrix \(J_N\) (diagonal \(\beta_0,\dots,\beta_{N-1}\), off-diagonal \(\alpha_1,\dots,\alpha_{N-1}\)). Then the quadrature nodes \(\{x_n\}\) are eigenvalues of \(J_N\), and the weights satisfy \(\frac1{w_n}=\sum_{k=0}^{N-1}p_k^*(x_n)^2>0\). Four common types: Gauss-Legendre ($(-1,1)$, \(w=1\)); Gauss-Chebyshev ($(-1,1)$, \(w=1/\sqrt{1-x^2}\), for Fourier coefficients); Gauss-Hermite (\((-\infty,\infty)\), \(w=e^{-x^2}\), for normal expectations); Gauss-Laguerre (\((0,\infty)\), \(w=e^{-x}\), for exponential expectations).

3.1 早期方法与 Tauchen-Hussey / 3.1 Earlier methods and Tauchen-Hussey 若仅解单个含期望的优化（如静态组合问题），高精度高斯求积自然合适；但动态问题需算条件期望，且希望求积节点预先固定（不依赖模型状态）以降复杂度。离散化（把连续随机过程近似为有限状态马尔可夫链）正是工具。以 Gaussian AR(1) \(x_t=\rho x_{t-1}+\varepsilon_t\)（\(\varepsilon_t\sim N(0,\sigma^2)\)）为例，\(x_t\mid x_{t-1}\sim N(\rho x_{t-1},\sigma^2)\)。Tauchen (1986) 经典但不准（勿用）；Rouwenhorst (1995) 好（条件矩精确至 2 阶、构造性、\(\rho\ge0.99\) 时尤佳）。Tauchen-Hussey (1991) 用 Gauss-Hermite：以节点 \(\{x_n\}\)、权重 \(\{w_n\}\) 离散 \(N(0,\sigma^2)\)，由 \(\mathbb E[g(X)]=\int g(x)\frac1{\sqrt{2\pi\sigma^2}}e^{-x^2/2\sigma^2}dx\approx\sum_n\frac{w_n}{\sqrt\pi}g(\sqrt2\sigma x_n)\)，取节点 \(x_n'=\sqrt2\sigma x_n\)、权重 \(w_n'=w_n/\sqrt\pi\)。AR(1) 转移矩阵If one only solves a single optimization with an expectation (e.g. a static portfolio problem), a high-accuracy Gaussian quadrature is natural; but dynamic problems need conditional expectations, and it is desirable for the quadrature nodes to be preassigned (independent of the model state) to cut complexity. Discretization (approximating a continuous process by a finite-state Markov chain) is the tool. Take the Gaussian AR(1) \(x_t=\rho x_{t-1}+\varepsilon_t\) (\(\varepsilon_t\sim N(0,\sigma^2)\)), so \(x_t\mid x_{t-1}\sim N(\rho x_{t-1},\sigma^2)\). Tauchen (1986) is classic but inaccurate (do not use); Rouwenhorst (1995) is good (conditional moments exact to order 2, constructive, especially for \(\rho\ge0.99\)). Tauchen-Hussey (1991) uses Gauss-Hermite: discretize \(N(0,\sigma^2)\) with nodes \(\{x_n\}\), weights \(\{w_n\}\) via \(\mathbb E[g(X)]=\int g(x)\frac1{\sqrt{2\pi\sigma^2}}e^{-x^2/2\sigma^2}dx\approx\sum_n\frac{w_n}{\sqrt\pi}g(\sqrt2\sigma x_n)\), taking nodes \(x_n'=\sqrt2\sigma x_n\), weights \(w_n'=w_n/\sqrt\pi\). The AR(1) transition matrix is

$$p_{mn}\propto w_n'\,e^{-\frac{\mu^2-2x_n'\mu}{2\sigma^2}},\qquad \mu=\rho x_m',\quad \sum_{n=1}^N p_{mn}=1.$$

\(\rho\le0.5\) 时不错，但假设高斯冲击，且 \(\rho\) 大时变差。It is decent for \(\rho\le0.5\), but assumes Gaussian shocks and deteriorates for large \(\rho\).

3.2 Farmer-Tanaka-Toda 最大熵法 / 3.2 Farmer-Tanaka-Toda maximum entropy Tanaka-Toda (2013, 2015)、Farmer-Toda (2017) 给出更准、更通用的离散化（应作首选）。3.2.1 离散概率分布：给定密度 \(f\)，先取求积公式 \(\mathbb E[g(X)]=\int g(x)f(x)dx\approx\sum_{n=1}^N w_n g(x_n)f(x_n)\) (1)，得粗近似 \(q_n=\frac{w_n f(x_n)}{\sum_n w_n f(x_n)}\) (2)。但 \(\{q_n\}\) 的矩一般不匹配 \(f\)。Tanaka-Toda (2013) 通过更新概率精确匹配一组矩：设矩函数 \(T:\mathbb R^K\to\mathbb R^L\)、精确矩 \(\bar T=\int T(x)f(x)dx\)，解Tanaka-Toda (2013, 2015), Farmer-Toda (2017) give a more accurate and general discretization (should be the first choice). 3.2.1 Discretizing a distribution: given density \(f\), first pick a quadrature \(\mathbb E[g(X)]=\int g(x)f(x)dx\approx\sum_{n=1}^N w_n g(x_n)f(x_n)\) (1), giving a coarse approximation \(q_n=\frac{w_n f(x_n)}{\sum_n w_n f(x_n)}\) (2). But the moments of \(\{q_n\}\) generally do not match \(f\). Tanaka-Toda (2013) match a finite set of moments exactly by updating the probabilities: with a moment function \(T:\mathbb R^K\to\mathbb R^L\) and exact moments \(\bar T=\int T(x)f(x)dx\), solve

$$\min_{\{p_n\}}\sum_{n=1}^N p_n\log\frac{p_n}{q_n}\quad\text{s.t.}\quad\sum_{n=1}^N p_n T(x_n)=\bar T,\ \sum_{n=1}^N p_n=1,\ p_n\ge0.\tag{P}$$

目标是 \(\{p_n\}\) 相对 \(\{q_n\}\) 的 Kullback-Leibler 信息（相对熵）：在精确匹配矩的同时尽量贴近初始近似。(P) 是凸最小化，解唯一。转为低维（\(L\) 个未知）无约束凹最大化对偶问题The objective is the Kullback-Leibler information (relative entropy) of \(\{p_n\}\) relative to \(\{q_n\}\): match the moments exactly while staying as close to the initial approximation as possible. (P) is convex, so the solution is unique. Convert to the low-dimensional (\(L\) unknowns) unconstrained concave maximization dual

$$\max_{\lambda\in\mathbb R^L}\left[\lambda'\bar T-\log\left(\sum_{n=1}^N q_n e^{\lambda'T(x_n)}\right)\right].\tag{D}$$

定理 8：(P) 有解当且仅当 \(\bar T\in\mathrm{co}\,T(D_N)\)；有解则唯一。可由对偶解 \(\lambda_N\) 算更新概率 \(p_n\propto q_n e^{\lambda_N'T(x_n)}\)；若无解则正则性 \(\bar T\in\mathrm{int}\,\mathrm{co}\,T(D_N)\) 不成立，须选更少的矩。3.2.2 离散一般马尔可夫过程：对时齐一阶马尔可夫过程 \(P(x_t\le x'\mid x_{t-1}=x)=F(x',x)\)，对每个条件分布分别施用上法——给定网格 \(D_N=\{x_n\}\) 与初始转移矩阵 \(Q=(q_{nn'})\)、条件矩 \(\bar T_n=\mathbb E[T(x_t)\mid x_n]\)，对每个 \(n\) 解 (P\(_n\))/(D\(_n'\))，正则性 \(\bar T_n\in\mathrm{int}\,\mathrm{co}\,T(D_N)\)。即 Algorithm 1（选网格与初始近似 → 选矩函数 → 逐行解对偶 → 得转移矩阵）。Theorem 8: (P) has a solution iff \(\bar T\in\mathrm{co}\,T(D_N)\); if so, it is unique. From the dual solution \(\lambda_N\) one computes the updated probabilities \(p_n\propto q_n e^{\lambda_N'T(x_n)}\); if no solution exists, the regularity condition \(\bar T\in\mathrm{int}\,\mathrm{co}\,T(D_N)\) fails and one must select fewer moments. 3.2.2 Discretizing a general Markov process: for a time-homogeneous first-order Markov process \(P(x_t\le x'\mid x_{t-1}=x)=F(x',x)\), apply the above to each conditional distribution separately — given a grid \(D_N=\{x_n\}\), an initial transition matrix \(Q=(q_{nn'})\), and conditional moments \(\bar T_n=\mathbb E[T(x_t)\mid x_n]\), solve (P\(_n\))/(D\(_n'\)) for each \(n\), with regularity \(\bar T_n\in\mathrm{int}\,\mathrm{co}\,T(D_N)\). This is Algorithm 1 (choose grid and initial approximation → choose moment function → solve the dual row by row → obtain the transition matrix).