矩阵运算 - 网络统计学类函数(3)(2017-03-13 银河统计)

在数学中,矩阵(Matrix)是一个按照长方阵列排列的复数或实数集合,最早来自于方程组的系数及常数所构成的方阵。矩阵是高等代数学中的常见工具,也常用于统计分析等应用数学学科中。在Javascript中,矩阵和数组都可以用Array处理,但由数学定义可知,矩阵的行、列元素数量是规范的,数组则不一定。例如,一个\(2\times 3\)矩阵用JS数组可表示为“[[1,2,3],[4,5,6]]”,矩阵每行、每列元素必须相等。而数组则可以表示为“[[1,2,3],[4,5]]”。所以说“在JS中,矩阵式规范化的数组”。

在网络统计学中,矩阵类函数(webTJ.Matrix)是网络统计学大类(webTJ)下的一个重要子类,主要用于和线性代数密切相关的多元统计分析、多元回归末等方法的计算。

为方便运行本文中样例代码,可打开网络统计学代码调试窗口,复制、粘贴代码到数据处理代码窗口中运行即可。

矩阵函数(子类名称:webTJ.Matrix)一览表

序号 函数名称 参数1 参数2 参数3 功能 备注
1 getMEmpty(arrs) 二维数组 * * 生成空矩阵 *
2 getMCopy(arrs) 二维数组 * * 复制矩阵 *
3 getPlus(arrs1,arrs2) 二维数组1 二维数组2 * 矩阵加 数组1和数组2行列相同
4 getMinus(arrs1,arrs2) 二维数组1 二维数组2 * 矩阵减 同上
5 getMultiply(arrs1,arrs2) 二维数组1 二维数组2 * 矩阵乘 前列后行相同
6 getTranspose(arrs) 二维数组 * * 矩阵转置 *
7 getInverse(arrs) 二维数组 * * 矩阵求逆 方阵
8 getXTX(arrs) 二维数组 * * 转置矩阵乘原矩阵($X^TX$) *
9 getXTY(xarrs,yarr) 自变量数组 因变量数组 * 转置矩阵乘因变量矩阵($X^TY$) *
10 getInsertRow(arrs,row) 二维数组 插入行位置 * 矩阵添加行 *
11 getInsertRRow(arrs,rarr,row) 二维数组 给定行数组 插入行位置 矩阵添加给定行 *
12 getInsertCol(arrs,col) 二维数组 插入列位置 * 矩阵添加列 *
13 getInsertRCol(arrs,carr,col) 二维数组 给定列数组 插入列位置 矩阵添加给定列 *
14 getRemoveRow(arrs,row) 二维数组 删除行位置 * 矩阵删除行 *
15 getRemoveCol(arrs,col) 二维数组 删除列位置 * 矩阵删除列 *
16 getEig(arrs) 二维数组 * * 矩阵特征值和特征向量 方阵
17 getRandom(rows,cols) 矩阵行数 矩阵列数 * 二生成随机数矩阵 *
18 getDet(arrs) 二维数组 * * 计算矩阵行列式 方阵
19 getDiag(arr) 一维数组 * * 生成对角矩阵 *
20 getIdentity(rank) 矩阵阶数 * * 生成单位阵 *
21 getSVD(arrs) 二维数组 * * 矩阵SVD分解 *
22 getSparse(arrs) 二维数组 * * 稀疏矩阵压缩 *
23 getINVSparse(arrs) 二维数组 * * 稀疏矩阵解压 *
24 getCov(arrs,k) 二维数组 估计类型 * 样本协方差矩阵 k=0有偏、k=1无偏

注:本网页中所有数据管理类函数和代码样例都可以复制、粘贴到网页尾部“代码窗口”运行通过

一、矩阵运算函数手册###

1、生成空矩阵  [返回]

## 函数
    webTJ.Matrix.getMEmpty(arrs);
##参数
    【arrs】
    【二维数组】

代码样例

webTJ.clear();
var oArrs=webTJ.Matrix.getMEmpty(4,5);        //建立4行5列空数组
oArrs[3][3]=100;                              //矩阵第4行第4列复制为100
oArrs[0][0]=200;                              //矩阵第1行第1列复制为200
oArrs[2][1]=300;                              //矩阵第3行第2列复制为300
webTJ.display(oArrs,1);

2、复制矩阵  [返回]

## 函数
    webTJ.Matrix.getMCopy(arrs);
##参数
    【arrs】
    【二维数组】

代码样例

webTJ.clear();
var oArrs=[
    [3,2,5,1],
    [2,5,4,3],
    [3,1,2,4],
    [2,1,1,5],
    [4,1,3,1]];
var oArrs1=webTJ.Matrix.getMCopy(oArrs);    //复制数组oArrs并赋值给变量oArrs1
webTJ.display(oArrs1,1);

3、矩阵加  [返回]

## 函数
    webTJ.Matrix.getPlus(arrs1,arrs1);
##参数
    【arrs1,arrs1】
    【二维数组1,二维数组1】

代码样例

webTJ.clear();
var oTxt1="3:2:6:5:2:7,5:3:6:5:9:6,5:2:1:5:2:6,1:2:6:5:2:0";
var oTxt2="2:1:2:3:2:2,5:3:5:5:2:6,5:2:1:3:2:2,3:2:2:5:1:1";
var oArr1=webTJ.getArrs(oTxt1,",",":");
var oArr2=webTJ.getArrs(oTxt2,",",":");
webTJ.display(oArr1,1);
webTJ.display(oArr2,1);
var oArr=webTJ.Matrix.getPlus(oArr1,oArr2);    //数组oArr1加oArr2并赋值给变量oArr
webTJ.display(oArr,1);

4、矩阵减  [返回]

## 函数
    webTJ.Matrix.getMinus(arrs1,arrs1);
##参数
    【arrs1,arrs1】
    【二维数组1,二维数组1】

代码样例

webTJ.clear();
var oTxt1="3:2:6:5:2:7,5:3:6:5:9:6,5:2:1:5:2:6,1:2:6:5:2:0";
var oTxt2="2:1:2:3:2:2,5:3:5:5:2:6,5:2:1:3:2:2,3:2:2:5:1:1";
var oArr1=webTJ.getArrs(oTxt1,",",":");
var oArr2=webTJ.getArrs(oTxt2,",",":");
webTJ.display(oArr1,1);
webTJ.display(oArr2,1);
var oArr=webTJ.Matrix.getMinus(oArr1,oArr2);    //数组oArr1减oArr2并赋值给变量oArr
webTJ.display(oArr,1);

5、矩阵乘  [返回]

## 函数
    webTJ.Matrix.getMultiply(arrs1,arrs1);
##参数
    【arrs1,arrs1】
    【二维数组1,二维数组1】

代码样例

webTJ.clear();
var oTxt1="3:2:6:5:2:7,5:3:6:5:9:6,5:2:1:5:2:6,1:2:6:5:2:0";
var oTxt2="2:1:2:3,5:5:2:6,5:2:2:2,3:2:2:5,3:4:2:1,5:2:6:4";
var oArr1=webTJ.getArrs(oTxt1,",",":");
var oArr2=webTJ.getArrs(oTxt2,",",":");
webTJ.display(oArr1,1);
webTJ.display(oArr2,1);
var oArr=webTJ.Matrix.getMultiply(oArr1,oArr2);    //数组oArr1乘oArr2并赋值给变量oArr
webTJ.display(oArr,1);

6、矩阵转置  [返回]

## 函数
    webTJ.Matrix.getTranspose(arrs);
##参数
    【arrs】
    【数组】

代码样例

webTJ.clear();
var oArr=[
    [3,2,6,5,2,7],
    [5,3,6,5,9,6],
    [5,2,1,5,2,6],
    [1,2,6,5,2,0]];
webTJ.display(oArr,1);
var oArr1=webTJ.Matrix.getTranspose(oArr);    //转置数组变量oArr并赋值给oArr1
webTJ.display(oArr1,1);

7、矩阵求逆  [返回]

## 函数
    webTJ.Matrix.getInverse(arrs);
##参数
    【arrs】
    【数组】

代码样例

webTJ.clear();
var oArrs=[[2,2,3],[2,1,2],[1,3,4]];
webTJ.display(oArrs,1);
var oArrs1=webTJ.Matrix.getInverse(oArrs);    //逆变换数组变量oArr并赋值给oArr1
webTJ.display(oArrs1,1);

8、转置矩阵乘原矩阵(\(X^TX\))  [返回]

## 函数
    webTJ.Matrix.getXTX(arrs);
##参数
    【arrs】
    【数组】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oArrs1=webTJ.Matrix.getXTX(oArrs);    //oArr转置乘原数组变量oArr并赋值给oXTX
webTJ.display(oArrs1,1);

9、转置矩阵乘因变量矩阵(\(X^TY\))  [返回]

## 函数
    webTJ.Matrix.getXTY(xarrs,yarr);
##参数
    【xarrs,yarr】
    【自变量数组,因变量数组】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oY=[2,5,4,6,7];
webTJ.display(oY,1);
var oArrs1=webTJ.Matrix.getXTY(oArrs,oY);    //oArr转置乘一维数组变量oY并赋值给oXTY
webTJ.display(oArrs1,1);

10、矩阵添加行  [返回]

## 函数
    webTJ.Matrix.getInsertRow(arrs,row);
##参数
    【arrs,row】
    【数组,插入行位置】

注:插入行为数字1

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oTs=webTJ.Matrix.getInsertRow(oArrs,1);    //在第一行处插入一行
webTJ.display(oTs,1);

11、矩阵添加给定行  [返回]

## 函数
    webTJ.Matrix.getInsertRRow(arrs,rarr,row);
##参数
    【arrs,rarr,row】
    【数组,给定行数组,插入行位置】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oRowv=[3,4,1,7];
var oTs=webTJ.Matrix.getInsertRRow(oArrs,oRowv,1);
webTJ.display(oTs,1);

12、矩阵添加列  [返回]

## 函数
    webTJ.Matrix.getInsertCol(arrs,col);
##参数
    【arrs,col】
    【数组,插入列位置】

注:插入列为数字1

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oTs=webTJ.Matrix.getInsertCol(oArrs,1);
webTJ.display(oTs,1);

13、矩阵添加给定列  [返回]

## 函数
    webTJ.Matrix.getInsertRCol(arrs,carr,col);
##参数
    【arrs,carr,col】
    【数组,给定列数组,插入列位置】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oColv=[1,3,2,4,5];
var oTs=webTJ.Matrix.getInsertRCol(oArrs,oColv,1);
webTJ.display(oTs,1);

14、矩阵删除行  [返回]

## 函数
    webTJ.Matrix.getRemoveRow(arrs,row);
##参数
    【arrs,row】
    【数组,删除行位置】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oTs=webTJ.Matrix.getRemoveRow(oArrs,1);
webTJ.display(oTs,1);

15、矩阵删除列  [返回]

## 函数
    webTJ.Matrix.getRemoveCol(arrs,col);
##参数
    【arrs,col】
    【数组,删除列位置】

代码样例

webTJ.clear();
var oArrs=[[3,2,5,1],[2,5,4,3],[3,1,2,4],[2,1,1,5],[4,1,3,1]];
webTJ.display(oArrs,1);
var oTs=webTJ.Matrix.getRemoveCol(oArrs,1);
webTJ.display(oTs,1);

16、矩阵特征值和特征向量  [返回]

## 函数
    webTJ.Matrix.getEig(arrs);
##参数
    【arrs】
    【数组】

注:该函数返回复合数组,第一个数组为特征值数组,第二个数组为特征向量数组

代码样例

webTJ.clear();
var A = [[1,2,5],[3,5,-1],[7,-3,5]];
var oArrs=webTJ.Matrix.getEig(A);    //计算矩阵特征值和特征向量,返回数组
webTJ.display(oArrs[0],1);           //矩阵形式显示特征值数组oArrs[0]
webTJ.display(oArrs[1],1);           //矩阵形式显示特征向量数组oArrs[1]

17、生成随机数矩阵  [返回]

## 函数
    webTJ.Matrix.getRandom(rows,cols);
##参数
    【rows,cols】
    【矩阵行数,矩阵列数】

注:矩阵元素为0-1均匀分布随机数

代码样例

webTJ.clear();
var oArrs=webTJ.Matrix.getRandom(4,5);    //生成4行5列随机数矩阵
webTJ.display(oArrs,1);           

18、计算矩阵行列式  [返回]

## 函数
    webTJ.Matrix.getDet(arrs);
##参数
    【arrs】
    【数组】

代码样例

webTJ.clear();
var oArrs = [
  [6,8,4,2,8,5],
  [3,5,2,4,9,2],
  [7,6,8,3,4,5],
  [5,5,2,8,1,6],
  [3,2,2,4,2,2],
  [8,3,2,2,4,1]];
var oV=webTJ.Matrix.getDet(oArrs);
webTJ.display(oV,0);           

19、生成对角矩阵  [返回]

## 函数
    webTJ.Matrix.getDiag(arr);
##参数
    【arr】
    【一维数组】

代码样例

webTJ.clear();
var oArr = [2,5,2,3];
var oDiag = webTJ.Matrix.getDiag(oArr);//生成4x4对角矩阵
webTJ.display(oDiag,1);

20、生成单位阵  [返回]

## 函数
    webTJ.Matrix.getIdentity(rank);
##参数
    【rank】
    【矩阵阶数】

代码样例

webTJ.clear();
var oIdentity = webTJ.Matrix.getIdentity(5);    //生成5阶角矩阵
webTJ.display(oIdentity,1);

21、矩阵SVD分解  [返回]

## 函数
    webTJ.Matrix.getSVD(arrs);
##参数
    【arrs】
    【二维数组】

代码样例

webTJ.clear();
var oArrs = [
  [22,10,2,3,7],
  [14,7,10,0,8],
  [-1,13,-1,-11,3],
  [-3,-2,13,-2,4],
  [9,8,1,-2,4],
  [9,1,-7,5,-1],
  [2,-6,6,5,1],
  [4,5,0,-2,2]];
var oBrrs=webTJ.Matrix.getSVD(oArrs);
webTJ.display(oBrrs[0],1);            //显示U值
webTJ.display(oBrrs[1],1);            //显示S值
webTJ.display(oBrrs[2],1);            //显示V值

22、稀疏矩阵压缩  [返回]

## 函数
    webTJ.Matrix.getSparse(arrs);
##参数
    【arrs】
    【二维数组】

代码样例

webTJ.clear();
var oArrs = [
  [ 3, 5, 8,10, 8],
  [ 7,10, 3, 5, 3],
  [ 6, 3, 5, 1, 8],
  [ 2, 6, 7, 1, 2],
  [ 1, 2, 9, 3, 9]];
var oBrrs=webTJ.Matrix.getSparse(oArrs);    //稀疏矩阵压缩并赋值
webTJ.display(oBrrs,1);           

23、稀疏矩阵解压  [返回]

## 函数
    webTJ.Matrix.getINVSparse(arrs);
##参数
    【arrs】
    【二维数组】

代码样例

var oArrs = [
  [ 3, 5, 8,10, 8],
  [ 7,10, 3, 5, 3],
  [ 6, 3, 5, 1, 8],
  [ 2, 6, 7, 1, 2],
  [ 1, 2, 9, 3, 9]];
var oArrs1=webTJ.Matrix.getSparse(oArrs);        //稀疏矩阵压缩
webTJ.display(oArrs1,1);
var oArrs2=webTJ.Matrix.getINVSparse(oArrs1);    //稀疏矩阵解压
webTJ.display(oArrs2,1);

24、样本协方差矩阵  [返回]

## 函数
    webTJ.Matrix.getCov(arrs,k)
##参数
    【arrs,k】
    【二维数组,估计类型】

注:k=0时取n、有偏估计;k=1时取n-1、无偏估计

代码样例

webTJ.clear();
var oArrs=[[1,2,4],[3,4,2],[4,6,5],[2,3,7]];    \\定义矩阵
webTJ.display(oArrs,1);
var oCov=webTJ.Matrix.getCov(oArrs,1);\\计算样本协方差矩阵(无偏估计)
webTJ.display(oCov,1);

二、矩阵类函数在数学和统计学中的运用##

目的:线性方程组和多元回归模型矩阵代码表示和矩阵类函数编程

1、根据克莱姆法则解线性方程组###

设有线性方程组一般形式如下:

\[\begin{equation*} \begin{cases} A_{11}X_1+A_{12}X_2+\dots+A_{1n}X_n=B_1\\ A_{21}X_1+A_{22}X_2+\dots+A_{2n}X_n=B_2\\ \dots\\ A_{n1}X_1+A_{n2}X_2+\dots+A_{nn}X_n=B_n\\ \end{cases} \end{equation*} \]

系数构成的行列式称为该方程组的系数行列式D,

\[\left| \begin{array}{ccc} A_{11} & A_{12} & \dots & A_{1n}\\ A_{21} & A_{22} & \dots & A_{2n}\\ \dots & \dots &\dots & \dots\\ A_{n1} & A_{n2} & \dots & A_{nn} \end{array} \right|\]

若线性方程组的系数矩阵可逆(非奇异),即系数行列式 D≠0。有唯一解,其解为,

\[X_j=\frac{D_j}{D}\hspace{1cm}(j=1,2,\dots,n) \]

其中\(D_j\)是把D中第j列元素对应地换成常数项而其余各列保持不变所得到的行列式。

现有线性方程组如下:

\[\begin{equation*} \begin{cases} 2X_1+X_2-5X_3+X_4=8\\ X_1-3X_2-6X_4=9\\ 2X_2-X_3+2X_4=-5\\ X_1+4X_2-7X_3+6X_4=0 \end{cases} \end{equation*} \]

该方程组系数行列式为,

\[D=\left| \begin{array}{ccc} 2 & 1 & -5 & 1\\ 1 & -3 & 0 & -6\\ 0 & 2 & -1 & 2\\ 1 & 4 & -7 & 6 \end{array} \right|=27\]

\(D\neq0\),故可用克莱姆法则求解,其中,

\[D_1=\left| \begin{array}{ccc} 8 & 1 & -5 & 1\\ 9 & -3 & 0 & -6\\ -5 & 2 & -1 & 2\\ 0 & 4 & -7 & 6 \end{array} \right|=81\]

\[D_2=\left| \begin{array}{ccc} 2 & 8 & -5 & 1\\ 1 & 9 & 0 & -6\\ 0 & -5 & -1 & 2\\ 1 & 0 & -7 & 6 \end{array} \right|=-108\]

\[D_3=\left| \begin{array}{ccc} 2 & 1 & 9 & 1\\ 1 & -3 & 9 & -6\\ 0 & 2 & -5 & 2\\ 1 & 4 & 0 & 6 \end{array} \right|=-27\]

\[D_4=\left| \begin{array}{ccc} 2 & 1 & -5 & 8\\ 1 & -3 & 0 & 9\\ 0 & 2 & -1 & -5\\ 1 & 4 & -7 & 0 \end{array} \right|=27\]

解得,\(X_1=\frac{81}{27}=3,X_2=\frac{-108}{27}=-4,X_3=\frac{-27}{27}=-1,X_4=\frac{27}{27}=1\)

样例代码

1.  webTJ.clear();
2.  var oDs=[[2,1,-5,1],[1,-3,0,-6],[0,2,-1,2],[1,4,-7,6]];
3.  var oY=[8,9,-5,0];
4.  webTJ.display(oDs,1);
5.  var oDet=webTJ.Matrix.getDet(oDs);
6.  oDet=webTJ.getDecimal(oDet,8);
7.  webTJ.display("D="+oDet,0);
8.  var oDD=webTJ.Matrix.getRemoveCol(oDs,0);
9.  var oD1=webTJ.Matrix.getInsertRCol(oDD,oY,0);
10. webTJ.display(oD1,1);
11. var oDet1=webTJ.Matrix.getDet(oD1);
12. oDet1=webTJ.getDecimal(oDet1,8);
13. webTJ.display("D1="+oDet1,0);
14. oDD=webTJ.Matrix.getRemoveCol(oDs,1);
15. var oD2=webTJ.Matrix.getInsertRCol(oDD,oY,1);
16. webTJ.display(oD2,1);
17. var oDet2=webTJ.Matrix.getDet(oD2);
18. oDet2=webTJ.getDecimal(oDet2,8);
19. webTJ.display("D2="+oDet2,0);
20. oDD=webTJ.Matrix.getRemoveCol(oDs,2);
21. var oD3=webTJ.Matrix.getInsertRCol(oDD,oY,2);
22. webTJ.display(oD3,1);
23. var oDet3=webTJ.Matrix.getDet(oD3);
24. oDet3=webTJ.getDecimal(oDet3,8);
25. webTJ.display("D3="+oDet3,0);
26. oDD=webTJ.Matrix.getRemoveCol(oDs,3);
27. var oD4=webTJ.Matrix.getInsertRCol(oDD,oY,3);
28. webTJ.display(oD4,1);
29. var oDet4=webTJ.Matrix.getDet(oD4);
30. oDet4=webTJ.getDecimal(oDet4,8);
31. webTJ.display("D4="+oDet4,0);
32. var oX=[oDet1/oDet,oDet2/oDet,oDet3/oDet,oDet4/oDet];
33. webTJ.display("方程的解为:",0);
34. webTJ.display(oX,1);

代码功能

2.  设置系数矩阵oDs
3.  设置常数项数组oY
5.  计算系数矩阵行列式
6.  保留8位小数(防止出现过长小数)
8.  删除系数矩阵第1列
9.  将常数项数组oY作为第1列插入系数矩阵
11. 计算替换系数矩阵D1行列式
32. 计算得各自变量的解

注:8-13、14-19、20-25、26-31功能类似,都是在系数矩阵D中分别替换各列为常数项数组oY,从而计算出替换系数矩阵\(D_1,D_2,D_3,D_4\)行列式的值

为了练习使用矩阵类函数和展示计算过程,使得上面的代码显得较多。如果只是为了解线性方程,可以简化代码如下:

样例代码

1.  webTJ.clear();
2.  var oD=[[2,1,-5,1],[1,-3,0,-6],[0,2,-1,2],[1,4,-7,6]];
3.  var oD1=[[8,1,-5,1],[9,-3,0,-6],[-5,2,-1,2],[0,4,-7,6]];
4.  var oD2=[[2,8,-5,1],[1,9,0,-6],[0,-5,-1,2],[1,0,-7,6]];
5.  var oD3=[[2,1,8,1],[1,-3,9,-6],[0,2,-5,2],[1,4,0,6]];
6.  var oD4=[[2,1,-5,8],[1,-3,0,9],[0,2,-1,-5],[1,4,-7,0]];
7.  oD=webTJ.Matrix.getDet(oD);
8.  oD1=webTJ.Matrix.getDet(oD1);
9.  oD2=webTJ.Matrix.getDet(oD2);
10. oD3=webTJ.Matrix.getDet(oD3);
11. oD4=webTJ.Matrix.getDet(oD4);
12. var oX=[oD1/oD,oD2/oD,oD3/oD,oD4/oD];
13. webTJ.display(oX,0);

注:显示的计算结果中的长小数是由于JS内部进制转换产生的。代码中2-6定义系数行列式和替换系数行列式,7-11计算所有行列式

2、试将下面线性方程组写成矩阵形式,并求解###

\[\begin{equation*} \begin{cases} 2X_1+2X_2+3X_3=2\\ X_1-X_2=2\\ -X_1+2X_2+X_4=4 \end{cases} \end{equation*} \]

矩阵形式为,\(AX=b\),其中,

\[A=\left(\begin{array}{ccc} 2 & 2 & 3 \\ 1 &-1 & 0 \\ -1 & 2 & 1 \end{array}\right), \hspace{1cm} X=\left(\begin{array}{ccc} x_{_1}\\ x_{_2}\\ x_{_3} \end{array}\right), \hspace{1cm} b=\left(\begin{array}{ccc} 2\\ 2\\ 4 \end{array}\right)\]

计算\(|A|=-1\neq0\),故A可逆. 因而有\(X=A^{-1}b\),即,

\[\left(\begin{array}{ccc} x_{_1}\\ x_{_2}\\ x_{_3} \end{array}\right)= \left(\begin{array}{ccc} 2 & 2 & 3 \\ 1 &-1 & 0 \\ -1 & 2 & 1 \end{array}\right)^{-1} \left(\begin{array}{ccc} 2\\ 2\\ 4 \end{array}\right)= \left(\begin{array}{ccc} 1 &-4 &-3 \\ 1 &-5 &-3 \\ -1 & 6 & 4 \end{array}\right) \left(\begin{array}{ccc} 2\\ 2\\ 4 \end{array}\right)= \left(\begin{array}{ccc} -18\\ -20\\ 26 \end{array}\right) \]

根据矩阵相等的定义,方程组的解为:\(x_{_1}=-18,x_{_2}=-20,x_{_3}=26\)

样例代码

1.  webTJ.clear();
2.  var oD=[[2,2,3],[1,-1,0],[-1,2,1]];
3.  var oY=[2,2,4];
4.  var oDet=webTJ.Matrix.getDet(oD);
5.  oDet=webTJ.getDecimal(oDet,8);
6.  webTJ.display("行列式="+oDet,0);
7.  var invD=webTJ.Matrix.getInverse(oD);
8.  invD=webTJ.getArrDecimal(invD,8);
9.  webTJ.display(invD,1);
10. var oArr=webTJ.Matrix.getMultiply(invD,oY);
11. webTJ.display(oArr,1);

代码功能

2.  定义系数矩阵
3.  定义常数项矩阵
4.  计算系数矩阵的行列式值
5.  保留标量有效小数8位(去除长小数)
7.  计算系数矩阵的逆矩阵
8.  保留向量(矩阵或数组)有效小数8位
10. 系数矩阵逆矩阵乘常数项矩阵(主要数项矩阵写法)

3、多元线性回归模型求解###

多元线性回归模型:\(Y=A_0+A_1X_1+A_2X_2+\dots+A_mX_m\)

矩阵形式:\(XA=Y\)。其中,

\[X=\left(\begin{array}{ccc} 1 & X_{11} & X_{21} & \dots & X_{m1}\\ 1 & X_{12} & X_{21} & \dots & X_{m1}\\ \dots & \dots & \dots & \dots & \dots\\ 1 & X_{1n} & X_{2n} & \dots & X_{mm}\\ \end{array}\right), \hspace{1cm} Y=\left(\begin{array}{ccc} Y_1\\ Y_2\\ \dots\\ Y_n \end{array}\right), \hspace{1cm} A=\left(\begin{array}{ccc} A_1\\ A_2\\ \dots\\ A_m \end{array}\right)\]

多元线性回归模型解得矩阵形式为,\((X^TX)^{-1}\times(X^TY)\)。其中,

\[\small{(X^TX)^{-1}=\left(\begin{array}{ccc} n & \sum\limits_{i=1}^n X_{1i} & \sum\limits_{i=1}^n X_{2i} & \dots & \sum\limits_{i=1}^n X_{mi}\\ \sum\limits_{i=1}^n X_{1i} & \sum\limits_{i=1}^n X_{1i}^2 & \sum\limits_{i=1}^n X_{1i}X_{2i} & \dots & \sum\limits_{i=1}^n X_{1i}X_{mi}\\ \sum\limits_{i=1}^n X_{2i} & \sum\limits_{i=1}^n X_{1i}X_{2i} & \sum\limits_{i=1}^n X_{2i}^2 & \dots & \sum\limits_{i=1}^n X_{2i}X_{mi}\\ \dots & \dots & \dots & \dots & \dots\\ \sum\limits_{i=1}^n X_{mi} & \sum\limits_{i=1}^n X_{1i}X_{mi} & \sum\limits_{i=1}^n X_{2i}X_{mi} & \dots & \sum\limits_{i=1}^n X_{mi}^2 \end{array}\right) \hspace{1cm} (X^TY)=\left(\begin{array}{ccc} \sum\limits_{i=1}^n Y_i\\ \sum\limits_{i=1}^n X_{1i}Y_i\\ \sum\limits_{i=1}^n X_{2i}Y_i\\ \dots\\ \sum\limits_{i=1}^n X_{mi}Y_i \end{array}\right)} \]

现有某地区居民不同类型收入和总消费统计资料(单位:元)如下,

序号 持久收入 临时收入 总消费
1 1800 781 2248
2 2100 1342 2531
3 2400 487 2220
4 2700 375 2469
5 3000 662 2972
6 3300 1145 3447
7 3600 1303 3686
8 3900 1322 4058
9 4200 808 3726
10 4500 643 3955

根据表中样本数据,运用多元线性回归模型进行参数估计。

样例代码

1.  webTJ.clear();
2.  var oX=[[1,1800,781],[1,2100,1342],[1,2400,487],[1,2700,375],[1,3000,662],
           [1,3300,1145],[1,3600,1303],[1,3900,1322],[1,4200,808],[1,4500,643]];
3.  var oY=[2248,2531,2220,2469,2972,3447,3686,4058,3726,3955];
4.  webTJ.display(oX,1);
5.  var oXT=webTJ.Matrix.getTranspose(oX);
6.  webTJ.display(oXT,1);
7.  var oXTX=webTJ.Matrix.getMultiply(oXT,oX);
8.  webTJ.display(oXTX,1);
9.  var oXTY=webTJ.Matrix.getMultiply(oXT,oY);
10. webTJ.display(oXTY,1);
11. var invXTX=webTJ.Matrix.getInverse(oXTX);
12. webTJ.display(invXTX,1);
13. var oA=webTJ.Matrix.getMultiply(invXTX,oXTY);
14. webTJ.display(oA,1);

代码功能

  1. 定义自变量样本矩阵(X)
  2. 定义因变量样本数组(Y)
  3. 自变量矩阵转置(\(X^T\))
  4. 自变量转置矩阵乘自变量矩阵(\(X^TX\))
  5. 自变量转置矩阵乘因变量数组(\(X^TY\))
  6. 自变量转置矩阵乘自变量矩阵的逆矩阵(\((X^TX)^{-1}\))
  7. 逆矩阵乘自变量转置矩阵乘因变量数组(\((X^TX)^{-1}\times(X^TY)\))

4、主成分分析###

现有全国30个省、直辖市的8项经济指标如下表,试进行主成分分析。

省份 国内生产 居民消费 固定资产 职工工资 货物周转 消费价格 商品零售 工业产值
北京 1394.89 2505 519.01 8144 373.9 117.3 112.6 843.43
天津 920.11 2720 345.46 6501 342.8 115.2 110.6 582.51
河北 2849.52 1258 704.87 4839 2033.3 115.2 115.8 1234.85
山西 1092.48 1250 290.9 4721 717.3 116.9 115.6 697.25
内蒙 832.88 1387 250.23 4134 781.7 117.5 116.8 419.39
辽宁 2793.37 2397 387.99 4911 1371.7 116.1 114 1840.55
吉林 1129.2 1872 320.45 4430 497.4 115.2 114.2 762.47
黑龙江 2014.53 2334 435.73 4145 824.8 116.1 114.3 1240.37
上海 2462.57 5343 996.48 9279 207.4 118.7 113 1642.95
江苏 5155.25 1926 1434.95 5943 1025.5 115.8 114.3 2026.64
浙江 3524.79 2249 1006.39 6619 754.4 116.6 113.5 916.59
安徽 2003.58 1254 474 4609 908.3 114.8 112.7 824.14
福建 2160.52 2320 553.97 5857 609.3 115.2 114.4 433.67
江西 1205.11 1182 282.84 4211 411.7 116.9 115.9 571.84
山东 5002.34 1527 1229.55 5145 1196.6 117.6 114.2 2207.69
河南 3002.74 1034 670.35 4344 1574.4 116.5 114.9 1367.92
湖北 2391.42 1527 571.68 4685 849 120 116.6 1220.72
湖南 2195.7 1408 422.61 4797 1011.8 119 115.5 843.83
广东 5381.72 2699 1639.83 8250 656.5 114 111.6 1396.35
广西 1606.15 1314 382.59 5105 556 118.4 116.4 554.97
海南 364.17 1814 198.35 5340 232.1 113.5 111.3 64.33
四川 3534 1261 822.54 4645 902.3 118.5 117 1431.81
贵州 55.98 1110 17.87 7382 4.2 117.3 114.9 5.57
陕西 1000.03 1208 300.27 4396 500.9 119 117 600.98
甘肃 553.35 1007 114.81 5493 507 119.8 116.5 468.79
青海 165.31 1445 47.76 5753 61.6 118 116.3 105.8
宁夏 169.75 1355 61.98 5079 121.8 117.1 115.3 114.4
新疆 834.57 1469 376.96 5348 339 119.7 116.7 428.76

基本计算步骤:

I、样本协方差矩阵

两个样本总体间的样本协方差公式:

\[Cov(X_i,X_j)=\frac{\sum\limits_{k=1}^n (X_{ik}-{\overline X_i})(X_{jk}-{\overline X_j})}{n-1} \]

m个样本总体间的样本协方差矩阵:

\[A=\left(\begin{array}{ccc} Cov(X_1,X_1) & Cov(X_2,X_1) & \dots & Cov(X_m,X_1)\\ Cov(X_1,X_2) & Cov(X_2,X_2) & \dots & Cov(X_m,X_2)\\ \dots & \dots & \dots & \dots\\ Cov(X_1,X_m) & Cov(X_2,X_m) & \dots & Cov(X_m,X_m) \end{array}\right)\]

\[\small{=\frac{1}{n-1}\left(\begin{array}{ccc} \sum\limits_{k=1}^n (X_{1k}-{\overline X_1})^2 & \sum\limits_{k=1}^n (X_{2k}-{\overline X_2})(X_{1k}-{\overline X_1}) & \dots & \sum\limits_{k=1}^n (X_{mk}-{\overline X_m})(X_{1k}-{\overline X_1})\\ \sum\limits_{k=1}^n (X_{1k}-{\overline X_1})(X_{2k}-{\overline X_2}) & \sum\limits_{k=1}^n (X_{i2}-{\overline X_2})^2 & \dots & \sum\limits_{k=1}^n (X_{mk}-{\overline X_m})(X_{2k}-{\overline X_2})\\ \dots & \dots & \dots & \dots\\ \sum\limits_{k=1}^n (X_{1k}-{\overline X_1})(X_{mk}-{\overline X_m}) & \sum\limits_{k=1}^n (X_{i2}-{\overline X_2})(X_{mk}-{\overline X_m}) & \dots & \sum\limits_{k=1}^n (X_{mk}-{\overline X_m})^2 \end{array}\right)} \]

II、协方差矩阵特征值

已知E单位矩阵,A为样本协方差矩阵(m阶方阵),求m阶矩阵A的特征值的基本方法即求齐次线性方程组\((\lambda E-A)X=0\)有非零解的值\(\lambda\)。即要求行列式\((\lambda E-A)=0\)。计算行列式获得的\(\lambda\)值即为矩阵A的特征值。

III、特征值取绝对值、排序(倒序)、累计、计算贡献率

注:这里着重计算,主成分分析具体方法参见后续文章。

样例代码

1.  webTJ.clear();
2.  var oStr = webTJ.getGSData("sysData");
3.  var oArrs=webTJ.getArrs(oStr,"|",":");
4.  oArrs=webTJ.Matrix.getRemoveCol(oArrs,0);
5.  var oTrr=oArrs[0];
6.  oArrs=webTJ.Matrix.getRemoveRow(oArrs,0);
7.  oArrs=webTJ.Array.getQuantify(oArrs);
8.  var oCov=webTJ.Matrix.getCov(oArrs,1);
9.  var oErrs=webTJ.Matrix.getEig(oCov);
10. var oVrrs=[];
11. var oLen=oTrr.length;
12. for (var i=0; i<oLen; i++) {
13.   oVrrs[i]=[]; oVrrs[i][0]=oTrr[i];
14.   oVrrs[i][1]=Math.abs(webTJ.getDecimal(oErrs[0][i],4));
15.   }
16. oVrrs=webTJ.Array.getNArrsSort(oVrrs,1,1);
17. var oS=0;
18. for (var i=0; i<oLen; i++) {
19.   oS+=oVrrs[i][1];
20.   oVrrs[i][2]=webTJ.getDecimal(oS,4); oVrrs[i][3]=0;
21.   }
22. for (var i=0; i<oLen; i++) {
23.   oVrrs[i][3]=webTJ.getDecimal(100*oVrrs[i][2]/oS,6)+"%";
24.   }
25. var oT=["指标","特征值","累计贡献额","累计贡献率"];
26. oVrrs=webTJ.Matrix.getInsertRRow(oVrrs,oT,0);
27 webTJ.show(oVrrs,2);

代码功能

2.  获得系统数据(表中30省份经济数据,字符)
3.  根据字符格式转换为矩阵
4.  删除矩阵第一列(省份名称列)
5.  将矩阵第一行赋值给指标数组oTrr
6.  删除矩阵第一行(指标名称行)
7.  量化矩阵,赋值存在数量字符
8.  计算样本协方差矩阵
9.  根据样本协方差矩阵计算特征值和特征向量
10. 定义二维输出数组
11. 定义指标数量
12-15. 循环中依列向输出数组中写入指标名称和特征值绝对值
16. 根据特征值对输出数组进行倒序排序
17. 定义累加变量及初始值
18-21. 循环中向输出数组中写入特征值累计值
22-24. 循环中向输出数组中写入特征值累计百分比(贡献率) 
25. 定义指标数组
26、将指标数组插入输出数组第一行

三、在线数据操作练习###


代码窗口

注:可将例题实例代码复制、粘贴到“代码窗口”,点击“运行代码”获得计算结果(鼠标选择实例代码$\rightarrow$Ctrl+C:复制$\rightarrow$鼠标点击“代码窗口”使其获得焦点$\rightarrow$Ctrl+V:粘贴)

运行效果

相关文章:

  • 2022-02-15
  • 2022-02-08
  • 2022-12-23
  • 2022-12-23
  • 2021-12-03
猜你喜欢
  • 2021-04-10
  • 2021-12-15
  • 2021-11-18
  • 2021-05-23
  • 2022-12-23
  • 2021-11-22
相关资源
相似解决方案