矩阵运算 - 网络统计学类函数（3）

（2017-03-13 银河统计）

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

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

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

矩阵函数（子类名称：webTJ.Matrix）一览表

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

## 函数
    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、转置矩阵乘因变量矩阵（XTY）  [返回]

## 函数
    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、根据克莱姆法则解线性方程组###

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

{A11X1+A12X2+⋯+A1nXn=B1A21X1+A22X2+⋯+A2nXn=B2…An1X1+An2X2+⋯+AnnXn=Bn

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

|A11A12…A1nA21A22…A2n…………An1An2…Ann|

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

Xj=DjD(j=1,2,…,n)

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

现有线性方程组如下：

{2X1+X2−5X3+X4=8X1−3X2−6X4=92X2−X3+2X4=−5X1+4X2−7X3+6X4=0

该方程组系数行列式为，

D=|21−511−30−602−1214−76|=27

因D≠0，故可用克莱姆法则求解，其中，

D1=|81−519−30−6−52−1204−76|=81

D2=|28−51190−60−5−1210−76|=−108

D3=|21911−39−602−521406|=−27

D4=|21−581−30902−1−514−70|=27

解得，X1=8127=3,X2=−10827=−4,X3=−2727=−1,X4=2727=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，从而计算出替换系数矩阵D1,D2,D3,D4行列式的值

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

样例代码

复制代码
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、试将下面线性方程组写成矩阵形式，并求解###

{2X1+2X2+3X3=2X1−X2=2−X1+2X2+X4=4

矩阵形式为，AX=b，其中，

A=(2231−10−121),X=(x1x2x3),b=(224)

计算|A|=−1≠0，故A可逆. 因而有X=A−1b，即，

(x1x2x3)=(2231−10−121)−1(224)=(1−4−31−5−3−164)(224)=(−18−2026)

根据矩阵相等的定义，方程组的解为：x1=−18,x2=−20,x3=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=A0+A1X1+A2X2+⋯+AmXm

矩阵形式：XA=Y。其中，

X=(1X11X21…Xm11X12X21…Xm1……………1X1nX2n…Xmm),Y=(Y1Y2…Yn),A=(A1A2…Am)

多元线性回归模型解得矩阵形式为，(XTX)−1×(XTY)。其中，

(XTX)−1=(n∑i=1nX1i∑i=1nX2i…∑i=1nXmi∑i=1nX1i∑i=1nX1i2∑i=1nX1iX2i…∑i=1nX1iXmi∑i=1nX2i∑i=1nX1iX2i∑i=1nX2i2…∑i=1nX2iXmi……………∑i=1nXmi∑i=1nX1iXmi∑i=1nX2iXmi…∑i=1nXmi2)(XTY)=(∑i=1nYi∑i=1nX1iYi∑i=1nX2iYi…∑i=1nXmiYi)

现有某地区居民不同类型收入和总消费统计资料（单位：元）如下，

序号持久收入临时收入总消费
118007812248
2210013422531
324004872220
427003752469
530006622972
6330011453447
7360013033686
8390013224058
942008083726
1045006433955

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

样例代码

复制代码
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. 自变量矩阵转置(XT)
  4. 自变量转置矩阵乘自变量矩阵(XTX)
  5. 自变量转置矩阵乘因变量数组(XTY)
  6. 自变量转置矩阵乘自变量矩阵的逆矩阵((XTX)−1)
  7. 逆矩阵乘自变量转置矩阵乘因变量数组((XTX)−1×(XTY))

4、主成分分析###

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

省份国内生产居民消费固定资产职工工资货物周转消费价格商品零售工业产值
北京1394.892505519.018144373.9117.3112.6843.43
天津920.112720345.466501342.8115.2110.6582.51
河北2849.521258704.8748392033.3115.2115.81234.85
山西1092.481250290.94721717.3116.9115.6697.25
内蒙832.881387250.234134781.7117.5116.8419.39
辽宁2793.372397387.9949111371.7116.11141840.55
吉林1129.21872320.454430497.4115.2114.2762.47
黑龙江2014.532334435.734145824.8116.1114.31240.37
上海2462.575343996.489279207.4118.71131642.95
江苏5155.2519261434.9559431025.5115.8114.32026.64
浙江3524.7922491006.396619754.4116.6113.5916.59
安徽2003.5812544744609908.3114.8112.7824.14
福建2160.522320553.975857609.3115.2114.4433.67
江西1205.111182282.844211411.7116.9115.9571.84
山东5002.3415271229.5551451196.6117.6114.22207.69
河南3002.741034670.3543441574.4116.5114.91367.92
湖北2391.421527571.684685849120116.61220.72
湖南2195.71408422.6147971011.8119115.5843.83
广东5381.7226991639.838250656.5114111.61396.35
广西1606.151314382.595105556118.4116.4554.97
海南364.171814198.355340232.1113.5111.364.33
四川35341261822.544645902.3118.51171431.81
贵州55.98111017.8773824.2117.3114.95.57
陕西1000.031208300.274396500.9119117600.98
甘肃553.351007114.815493507119.8116.5468.79
青海165.31144547.76575361.6118116.3105.8
宁夏169.75135561.985079121.8117.1115.3114.4
新疆834.571469376.965348339119.7116.7428.76

基本计算步骤：

I、样本协方差矩阵

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

Cov(Xi,Xj)=∑k=1n(Xik−X¯i)(Xjk−X¯j)n−1

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

A=(Cov(X1,X1)Cov(X2,X1)…Cov(Xm,X1)Cov(X1,X2)Cov(X2,X2)…Cov(Xm,X2)…………Cov(X1,Xm)Cov(X2,Xm)…Cov(Xm,Xm))

=1n−1(∑k=1n(X1k−X¯1)2∑k=1n(X2k−X¯2)(X1k−X¯1)…∑k=1n(Xmk−X¯m)(X1k−X¯1)∑k=1n(X1k−X¯1)(X2k−X¯2)∑k=1n(Xi2−X¯2)2…∑k=1n(Xmk−X¯m)(X2k−X¯2)…………∑k=1n(X1k−X¯1)(Xmk−X¯m)∑k=1n(Xi2−X¯2)(Xmk−X¯m)…∑k=1n(Xmk−X¯m)2)

II、协方差矩阵特征值

已知E单位矩阵，A为样本协方差矩阵（m阶方阵），求m阶矩阵A的特征值的基本方法即求齐次线性方程组(λE−A)X=0有非零解的值λ。即要求行列式(λE−A)=0。计算行列式获得的λ值即为矩阵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、将指标数组插入输出数组第一行

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

代码窗口

注：可将例题实例代码复制、粘贴到“代码窗口”，点击“运行代码”获得计算结果（鼠标选择实例代码→Ctrl+C:复制→鼠标点击“代码窗口”使其获得焦点→Ctrl+V:粘贴）

运行效果

posted @   银河统计  阅读(1625)  评论(0)  编辑  收藏  举报
编辑推荐：
· 10年+ .NET Coder 心语，封装的思维：从隐藏、稳定开始理解其本质意义
· .NET Core 中如何实现缓存的预热？
· 从 HTTP 原因短语缺失研究 HTTP/2 和 HTTP/3 的设计差异
· AI与.NET技术实操系列：向量存储与相似性搜索在 .NET 中的实现
· 基于Microsoft.Extensions.AI核心库实现RAG应用
阅读排行：
· 10年+ .NET Coder 心语 ── 封装的思维：从隐藏、稳定开始理解其本质意义
· 地球OL攻略 —— 某应届生求职总结
· 提示词工程——AI应用必不可少的技术
· Open-Sora 2.0 重磅开源！
· 周边上新：园子的第一款马克杯温暖上架
点击右上角即可分享
微信分享提示
AI IDE Trae