从txt中读取矩阵和三元线性回归的C++实现.

++i; j = 0; } } // 打印输出数组 printArray(arr[0],Row,Col); fin.close(); cout << endl; // 计算线性方程组系数 double labArray[4][4]={0}; for(i=0;i<3;i++) { for(j=0;j<4;j++) { labArray[i][j] = Lab(arr[0],Row,Col,i,j); } } labArray[3][3] = Lab(arr[0],Row,Col,3,3); // print printArray(labArray[0],4,4); cout << endl; // 计算线性方程组的解 double fenmu = 0; double fenzi[3]; double b[3][3]={{labArray[0][0],labArray[0][1],labArray[0][2]}, {labArray[1][0],labArray[1][1],labArray[1][2]}, {labArray[2][0],labArray[2][1],labArray[2][2]}}; double tmp[3][3]; memcpy(tmp,b,sizeof(double)*9); fenmu = det(b[0],3); cout << fenmu << endl; for(i = 0;i < 3; ++i) { memcpy(tmp,b ,sizeof(double)*9); for(j = 0; j < 3 ;++j) {
三元线性回归主要是研究因变量y与三个自变量之间的线性关系。 其数学模型为：
具体参考线性回归相关理论计算公式。 （理论部分链接由于审核问题不便给出，请自行百度） 本次代码数据均来源文库中的数据，并经过验证可行。 下面给出其C++实现代码： /*功能：从.txt中读取一组n*4（n>3）的矩阵，（或者自己在代码中直 接输入数组矩阵）计算其线性回归系数并计算其相关系数*/ #include <iostream> #include <fstream> #include <iomanip> #include <string> using namespace std; // 打印二元数组，arr表示数组首地址，row表示行（1表示一元数 组），col表示列 void printArray(double* arr,int row,int col); // 获取.txt文件中矩阵的列数 int getColumns(const char * fileName); // 获取.txt文件中矩阵的行数 int getRows(const char * fileName); // 计算方程组的系数 row，col表示行列，a，b表示角标 double Lab(double *arr,int row,int col,int a,int b); // 计算三元数组行列式的值 double det(double *arr,int n=3); // 计算x的n次方 int power(int x,unsigned int n); int main() { double data; int Col,Row; int i,j; char c; // 从arr.txt中读取矩阵
} // 获取.txt文件中矩阵的行数 int getRows(const char *fileName) { ifstream ifs; int row=0; string tmp; ifs.open(fileName,std::ios::in); if(ifs) { while(getline(ifs,tmp,'\n')) { ++row; } ifs.close(); } return row ; } // 计算方程组的系数 row，col表示行列，a，b表示角标 double Lab(double *arr,int row,int col,int a,int b) { int i=0; double sum1=0,sum2=0,asum=0,bsum=0; for(i = 0;i < row; ++i) { sum1 += (*(arr+i*col+a)) * (*(arr+i*col+b)); asum += *(arr+i*col+a); bsum += *(arr+i*col+b); } sum2 = (asum*bsum)/row; return sum1-sum2; } // 计算三元数组行列式的值
-
return 0; } // 打印二元数组，arr表示数组首地址，row表示行（1表示一元数 组），col表示列 void printArray(double* arr,int row=1,int col=1) { int i,j; for(i=0;i<row;i++) { for(j=0;j<col;j++) { cout << setw(10) << *(arr+i*col+j); } cout << endl; } cout << endl; } // 获取.txt文件中矩阵的列数 int getColumns(const char * fileName) { ifstream ifs ; ifs.open(fileName,std::ios::in); double tmp; char c = ' '; int col=0; while(c != '\n') { ifs >> tmp; ++col; c = ifs.peek(); } ifs.close(); return col;
tmp[j][i] = labArray[j][3]; } fenzi[i] = det(tmp[0],3); } printArray(fenzi,1,3); double a[4]={0}; for(i = 0; i < 3; ++i) { a[i+1] = fenzi[i] / fenmu; } // 计算a[0] double sumX1=0,sumX2=0,sumX3=0,sumY=0; for(i = 0; i < Row; ++i) { sumX1 += arr[i][0]; sumX2 += arr[i][1]; sumX3 += arr[i][2]; sumY += arr[i][3]; } a[0] = sumY/Row - sumX1*a[1]/Row - sumX2*a[2]/Row sumX3*a[3]/Row; printArray(a,1,4); cout << endl; cout << "线性关系为： "; cout << " y = " << a[0] << " + " << a[1]; cout << "x1 + " << a[2] << "x2 + " << a[3] << "x3 " << endl; cout << endl; // 计算相关系数 double r=0,S=0; r = sqrt( (a[1]*labArray[0][0]+a[2]*labArray[1][0] +a[3]*labArray[2][2]) / labArray[3][3] ); S=sqrt((labArray[3][3] (a[1]*labArray[0][0]+a[2]*labArray[1][0] +a[3]*labArray[2][2])) / (Row - Col) ); cout << "相关系数 r ：" << setw(10) << r << endl; cout << "剩余标准差 S ：" << setw(10) << S << endl; cout << endl;
/*矩阵存储格式如下： 1111 1111 .... .... .... 1111 一二三四列分别代表x1,x2,x3,y */ ifstream fin("arr.txt",ios::in); // 获取矩阵的行和列 Col = getColumns("arr.txt"); Row = getRows("arr.txt"); // 输出 cout << "矩阵的列数是： " << Col << " " << "矩阵的行数是：" << Row; cout << endl; double arr[Row][Col]; // 数组初始化为0 for(i=0;i<Row;i++) { memset(arr[i],0,sizeof(double)*Col); } i = 0; j = 0; // 从文件中读取矩阵 while(!fin.eof()) { fin >> data; arr[i][j] = data; ++j; // 判断是否换行 c=fin.peek(); if('\n'==c) {
来自百度文库
double det(double *arr,int n=3) { double result = 0; double det1,det2,det3; int i; det1 = (*arr) * ((*(arr+1*n+1)) * (*(arr+2*n+2)) - (*(arr+1*n+2)) * (*(arr+2*n+1))); det2 = (*(arr+1)) * ((*(arr+1*n+0)) * (*(arr+2*n+2)) - (*(arr+1*n+2)) * (*(arr+2*n+0))); det3 = (*(arr+2)) * ((*(arr+1*n+0)) * (*(arr+2*n+1)) - (*(arr+1*n+1)) * (*(arr+2*n+0))); result = det1 - det2 + det3; return result; } // 计算x的n次方 int power(int x,unsigned int n) { int result = 1; if(0 == n) return 1; while(n) { result *= x; --n; } return result; }