利用C++vector实现OLS多元回归

使用最小二乘法实现多元回归，只需将数据格式化存入txt文本(格式化方法见代码内部)然后调整read_data("your_path")运行即可打印出模型中参数的估计值；
```cpp
#include<iostream>
#include<vector>
#include<fstream>
#include<sstream>
using namespace std;
typedef vector<double> vd;
typedef vector<vector<double>> vdd;
void matrixInversion(vdd a,vdd &b)  
{
	int n = a.size();
	int* is = new int[n];
	int* js = new int[n];
	int i, j, k;
	double d, p;
	for (k = 0; k < n; k++)
	{
		d = 0.0;
		for (i = k; i <= n - 1; i++)
			for (j = k; j <= n - 1; j++)
			{
				p = fabs(a[i][j]);
				if (p > d) { d = p; is[k] = i; js[k] = j; }
			}
		if (0.0 == d)
		{
			free(is); 
			free(js);
			printf("err**not inv\n");//奇异矩阵错误程序异常终止
			system("pause");
		}
		if (is[k] != k)
			for (j = 0; j <= n - 1; j++)
			{
				p = a[k][j];
				a[k][j] = a[is[k]][j];
				a[is[k]][j] = p;
			}
		if (js[k] != k)
			for (i = 0; i <= n - 1; i++)
			{
				p = a[i][k];
				a[i][k] = a[i][js[k]];
				a[i][js[k]] = p;
			}
		a[k][k] = 1.0 / a[k][k];
		for (j = 0; j <= n - 1; j++)
			if (j != k)
			{
				a[k][j] *= a[k][k];
			}
		for (i = 0; i <= n - 1; i++)
			if (i != k)
				for (j = 0; j <= n - 1; j++)
					if (j != k)
					{
						a[i][j] -= a[i][k] * a[k][j];
					}
		for (i = 0; i <= n - 1; i++)
			if (i != k)
			{
				a[i][k] = -a[i][k] * a[k][k];
			}
	}
	for (k = n - 1; k >= 0; k--)
	{
		if (js[k] != k)
			for (j = 0; j <= n - 1; j++)
			{
				p = a[k][j];
				a[k][j] = a[js[k]][j];
				a[js[k]][j] = p;
			}
		if (is[k] != k)
			for (i = 0; i <= n - 1; i++)
			{
				p = a[i][k];
				a[i][k] = a[i][is[k]];
				a[i][is[k]] = p;
			}
	}
	b = a;
	free(is); 
	free(js);
	return;
	//矩阵求逆有点生，这里参考了某位大神的代码
	//但是回头又找不到在哪了，但是绝非原创！！！对大神表示歉意

}
vdd read_data(string path)//设置一个函数来读取文件里面的数据
{
/*
		格式规范：
				1.因变量必须放在第一列
				2.只能是数字不带任何数据标签
				3.path必须是C++可以识别的格式
          
*/

	ifstream f;
	f.open(path);
	string str;
	vector<vector<double> > num;
	while (getline(f, str))
	{
		istringstream input(str);
		vector<double> tmp;
		int a;
		while (input >> a)
			tmp.push_back(a);
		num.push_back(tmp);
	}
	return num;
}
void add_constant(vdd& X)//向设计矩阵的第一行加上一来估计出常量
{
	int m = X[0].size();
	int n = X.size();//m 为列，n 为行
	for (int i = 0; i < n; i++)
		X[i].insert(X[i].begin(), 1.0);//向vector的开始处添加一个向量；
}
vdd mat_muti(vdd X1, vdd X2)//实现矩阵乘法运算
{
	int m1 = X1[0].size();
	int m2 = X2[0].size();
	int n1 = X1.size();
	int n2 = X2.size();//n为行，m为列  n1xm1  n2xm2
	vector<vector<double>> res;

	if (m1 != n2)
	{
		cout << "矩阵发生数据错误" << endl;
		system("pause");
	}
	else
	{
		for (int i = 0; i < n1; i++)
		{
			vd ress;
			for (int j = 0; j < m2; j++)
			{
				double v_res = 0.0;//储存向量乘法的结果
				for (int k = 0; k < m1; k++)
					for (int l = 0; l < n2; l++)
						v_res = v_res + X1[i][k] * X2[l][j];
				cout << "before this is okkk" << endl;
				ress.push_back(v_res);
			}
			res.push_back(ress);
		}
	}
	return res;
}
vdd  reverse(vdd x)
{
	vdd temp;
	for (int i = 0; i < x[0].size(); i++)
	{
		vd tp;
		for (int j = 0; j < x.size(); j++)
			tp.push_back(x[j][i]);
		temp.push_back(tp);
	}
	return temp;
}
vdd regression(vdd Y, vector<vector<double>> X)//Y初始化为外层只有一个元素的双层vector
{
	vdd reg_res;//最后回归的结果
	int m = X[0].size();//获取列
	int n = Y.size();//获取行
	int n1 = X.size();
	if (n1 != n)
	{
		cout << "发生输入错误请检查数据的完整性" << endl;
		system("pause");
	}
	add_constant(X);
	m++;//列数在原来矩阵上加1,得到设计矩阵
	vdd XT = reverse(X);//求转置矩阵
	vdd XI;
	matrixInversion(mat_muti(X, XT), XI);
	reg_res =mat_muti(mat_muti(XI, XT), Y);           //这里返回的还是双层vector但是依然可以查看结果
	return reg_res;
}
void print(vdd x)//测试函数用于打印VDD
{
	for (int i = 0; i < x.size(); i++)
	{
		for (int j = 0; j < x[i].size(); j++)
			cout << x[i][j] << "  ";
		cout << endl;
	}
}
void main()
{
//测试程序开始
	//测试add_constant函数
	vdd data,y;
	vd yy;
	//print(x);
	//add_constant测试成功!
	//vdd  z;
	//z = mat_muti(x, x);
	//print(z);
	//矩阵乘法测试成功
	//开始测试逆矩阵功能
	/*vdd inz;
	matrixInversion(x, inz);
	print(inz);//逆矩阵测试成功
	*/
	//vdd f;
	//f = reverse(x);
	//print(f);
	//完善regression函数的代码
	data = read_data("your_path");
	for (int i = 0; i < data.size(); i++)//提取x，y；
	{
		yy.push_back(data[i][0]);
		data[i].erase(data[i].begin());
	}
	//此时x就是data；
	y.push_back(yy);//必须初始化为单元素二维vector;
	vdd reg_result=regression(data, y);
	print(reg_result);//打印出估计出来的最小二乘参数；
}