模型估计

\(\)
最早对于模型估计的认识,来源于底下这组令我深恶痛绝的线性回归公式
$$\left\{
\begin{eqnarray}
\hat{b} &=& \frac{\sum_{i=1}^n(x_i-\bar{x})(y_I-\bar{y})}{\sum_{i=1}^n(x_i-\bar{x}^2)} = \frac{\sum_{i=1}^nx_iy_i-n\bar{x}\bar{y}}{\sum_{i=1}^nx_i^2-n\bar{x}^2} \\
\hat{a} &=& \bar{y}-\hat{b}\bar{x} \\
\end{eqnarray}
\right.$$
不是我记不住这个公式,而是在当时(现在也一样)落后的考试制度下,必须用计算器把一个个的数算出来。一般的考试时都会给5、6组数据,光摁计算器就要少说花两分钟。

另一个记忆深刻的是高中时期伏安法测电阻。这个的印象好多了,只要画个图,再算个斜率就得到电阻,整个求解过程满是轻松愉悦的气氛。

现在,也经常遇到的模型估计问题。比如,经常会需要做曲线拟合。虽然这类问题已经有无数现成的工具可以一步解决,但有些时候还是要自己动手,丰衣足食。下面开始正题,简单地介绍一下如何进行模型的估计。

模型估计问题看起来不难,但随便找本书都能写上百十来页。这篇文章就挑些简单的作简要的描述,难的我也不会。本文主要涉及两个方程,一个是$$\mathbf{A}\mathbf{x}=\mathbf{0}$$另一个是
$$\mathbf{A}\mathbf{x}=\mathbf{b}$$
或者换一个写法,就是
$$\arg\min_{\mathbf{x}} \|\mathbf{A}\mathbf{x}\|_2$$或$$\arg\min_{\mathbf{x}} \|\mathbf{A}\mathbf{x}-\mathbf{b}\|_2$$
其中,\(\mathbf{A}\)是矩阵,\(\mathbf{x}\),\(\mathbf{b}\),\(\mathbf{0}\)是维数与\(\mathbf{A}\)相适应的列向量,这里设\(\mathbf{x}\)是一个\(n\)维向量。

SVD分解

第一个方程一般用奇异值分解(SVD)求解。简单地说,就是可以把\(\mathbf{A}\)分解成
$$\mathbf{A} = \mathbf{U}\mathbf{D}\mathbf{V}^T$$
\(\mathbf{U}\)和\(\mathbf{V}\)都是正交矩阵(严格来说,是酉矩阵),\(\mathbf{D}\)为对角矩阵,而且对角线上的元素非负,从左上角到右下角递减。至于SVD为什么存在,移步维基百科。\(\mathbf{V}\)的最后一列就是使得\(\|\mathbf{A}\mathbf{x}\|_2\)取得最小值的单位向量。利用SVD的表达式,要验证这一点并不难。因为\(\mathbf{V}\)是正交矩阵,所以它的所有列\(\{\mathbf{v}_k\}\)构成\(\mathbf{x}\)所在空间的一组标准正交基,于是\(\mathbf{x} = \sum_{k}a_k\mathbf{v}_k\)。设\(\mathbf{D}=\mathrm{diag}\{\lambda_k\}\),利用SVD表达式计算\(\mathbf{Ax}\)的范数的平方
$$\begin{eqnarray}
\|\mathbf{Ax}\|_2^2 &=& \|\mathbf{UDV}^T\sum_ka_k\mathbf{v}_k\|_2^2 \\
&=& \|\mathbf{DV}^T\sum_ka_k\mathbf{v}_k\|_2^2 \quad\quad (\mathbf{U}\text{具有保范性})\\
&=& \|\mathbf{D}\sum_ka_k\mathbf{V}^T\mathbf{v}_k\|_2^2 \quad\quad (\text{利用}\mathbf{V}的正交性)\\
&=& \sum_k|\lambda_ka_k|_2^2 \quad\quad (\text{利用}\mathbf{V}的正交性和\mathbf{D}的对角性)\\
&\geq& \sum_k|\lambda_na_k|_2^2 \quad\quad (\lambda_{k+1}\geq\lambda_k)\\
&=& \lambda_n^2 \quad\quad (\mathbf{x}为单位向量)
\end{eqnarray}$$
上式说明在任何情况下,\(\mathbf{Ax}\)的范数都不可能小于\(\lambda_n\)。而当\(a_n = 1\)时,\(\|\mathbf{Ax}\|_2=\lambda_n\)。

最小二乘

一般称$$\arg\min_{\mathbf{x}} \|\mathbf{A}\mathbf{x}-\mathbf{b}\|_2$$的解为最小二乘解,其结果为$$\mathbf{x}^* = (\mathbf{A}^T\mathbf{A})^{-1}\mathbf{A}^T\mathbf{b}$$这个解可以直接通过对代价函数(即范数方程)求导,令其导数为零求得。另一个方案是利用几何结论。由几何关系,最优解应该满足$$\mathbf{A}^T(\mathbf{Ax}-\mathbf{b})=\mathbf{0}$$。这是一个很多人都知道的结论,但是有很多人反应不过来。

以三维空间为例,这个问题用人话描述就是已知空间中的一个点A,现要在平面上找到一个与A最近的点,怎么找?就像上图那样,过A点做平面的垂线,交点B就是满足条件的点。所谓与平面垂直,是指AB与平面上的所有向量都垂直,即积为零。由于这个平面是一个向量空间,所以只要AB与空间的一组基底正交即可。

这里的\(\mathbf{Ax}\)是\(\mathbf{A}\)的列向量所张成的空间中的一个向量。\(\mathbf{Ax}-\mathbf{b}\)就相当于三维空间中的AB。AB与空间中的所有向量正交,只需\(\mathbf{Ax}-\mathbf{b}\)与\(\mathbf{A}\)的每一列都正交,因此$$\mathbf{A}^T(\mathbf{Ax}-\mathbf{b})=\mathbf{0}$$从而可以得到最小二乘解为$$\mathbf{x}^* = (\mathbf{A}^T\mathbf{A})^{-1}\mathbf{A}^T\mathbf{b}$$

参数估计举例

接下来以最为简单的直线拟合来说明如何通过SVD分解或者是直接求最小二乘解进行模型的估计。我通过MATLAB生成了一组数据

x = [1:5];
y = 3*x + 4 + randn(1,5);

构造的数据如下表所示

\(x\) 1 2 3 4 5
\(y\) 6.7950 9.87592 14.4897 17.4090 20.4172

现在要通过模型构造出方程\(\mathbf{A}\mathbf{x}=\mathbf{0}\)或者是\(\mathbf{A}\mathbf{x}=\mathbf{b}\)。设直线方程为$$kx+d-y=0$$,写成矩阵形式为$$[x~~1~~-y][k~~d~~1]^T = 0$$或者是$$[x~~1][k~~d]^T = y$$
首先,通过SVD分解求解模型,利用MATLAB求解

n = length(x);
A = zeros(n,3);
for i=1:length(x)
    A(i,:) = [x(i) 1 -y(i)];
end
[~,~,v] = svd(A);
v_end = v(:,end);
v_end = v_end/v_end(3);
k = v_end(1)
d = v_end(2)

再用最小二乘模型,利用MATLAB求解

n = length(x);
A = zeros(n,2);
b = zeros(n,1);
for i=1:length(x)
    A(i,:) = [x(i) 1];
    b(i,:) = [y(i)];
end
result = (A'*A)\A'*b;
k = result(1)
d = result(2)

使用两种方法估计出的模型为

\(k\) \(d\)
SVD 3.4551 3.4549
最小二乘 3.4777 3.3641

这两种估计模型的方法不仅仅可用于直线拟合这样的简单问题,也能拟合其它的曲线,还可以拟合多项式、二次型等。经过一些变换,能胜任某些非常复杂的模型的估计。

最后,顺便说一下,SVD分解方法在很多情况下表现优于最小二乘解,SVD更不容易受到噪声的影响,应用也更多。

未经允许不得转载:Charlie小站 » 模型估计

赞 (6)
分享到:更多 ()