数据建模与分析线性回归小论文VIP

- 1 -上海住房面积和房价的线性回归分析王明黔（上海大学 机电工程与自动化学院，上海 200444） 摘要：在数据构建统计模型的学习中，统计学习是其一种基础的学习方法。本文针对城市人口数目与饮品连锁店利润的关系，就已有的数据进行线性回归分析，利用 Matlab 工具进行数据的线性回归模拟，进而得出城市人口数目与饮品连锁店利润的散点图、拟合直线图和三维等高线图。为了分析上海地区的住房面积和房价的关系，收集最近的售房成交数据，将数据导入到 Matlab 进行分析，得出上海房价与住房面积的线性关系。关键词：Matlab；线性回归；目标函数；梯度下降；统计学习基于数据的机器学习是现代智能技术中十分重要的一个方面，主要研究如何从一些观测数据(样本) 出发，得出目前尚不能通过原理性分析得到的规律，并用以对未来数据或无法观测的数据进行预测。现实生活中大量存在我们尚无法准确认识但却可以进行观测的事件。因此，这种机器学习在从现代科学、技术到社会、经济等各领域中都有着十分重要的应用 1。使用线性回归方法可以对一些观测数据进行分析，把预测事件中一些因素作为自变量，另一些随自变量变化而变化的变量作为因变量，研究它们之间的非确定性因果关系，以便预测因变量的未来发展趋势。根据若干观测数据寻找描述变量之间的函数或统计相关关系的最佳数学表达式，或者匹配数据之间相关关系的最佳拟合曲线，来表达随机性变量间的规律 2。利用线性回归通过多变量机器学习的方法，可以建立上海住房面积和价格的线性关系，建立数学模型并评估其中的未知参数。1 案例分析1.1 目标函数的建立根据已知给出的城市人口数目与饮品连锁店利润的一些数据，可以得到一个样本集，如图 1，为样本在 Matlab 软件加载数据图，第一列表示城市人口数目，第二列表示饮品连锁店利润。图 1 城市人口数目与饮品连锁店利润的样本集Fig 1 Urban population and beverage chain profits of sample set- 2 -针对样本集，我们可以假设一个线性模型：（1）01()hx式中： 假设的线性模型；()hx样本/特征；参数。01、其中， 为模型参数，因此问题就可以转换为，求出 的值。i 01和为了得到较为准确的回归线，应该满足一个前提，即各样本点尽可能分布在所建立的线性模型周围，因此我们建立目标函数：(2)()()2011(,)2miiiJhxy式中： 向量 中的第 个元素；()ixxi向量 中的第 个元素；()iyy模型假设；hx训练集的数量。m因此，我们只要求得使上述函数的值最小时的 的值。01和1.2 目标函数的求解我们首先从一组 开始，利用最速下降法不断改变 的值来减小 ，直到达到我01和 01和 01J（ ， ）们希望得到的最小值 3。对于最速下降法，可利用以下公式而来求解：(3)j 01j=-J： （ ， ）式中： 学习速率（步长） ；下降方向。01(,)jJ其中 为修正学习速率，即搜索的步长，当 过小，会延长搜索的时间，当 过大，可能会错过极值点，导致不收敛。还有要注意最速下降法会收敛到局部最小值的情况，这种情况即便修正学习速率也不会改变结果，并且接近局部最小值时，因为斜率的减小最速下降法会自动减少每步的补偿 4。因此，我们要现注意目标函数的极小值情况，如果有局部最小值可通过调节 参数的大小来跳过局部最小值情况，否则，就不需跳过，最终通过多次迭代可以得到最终的理想结果。- 3 -2 Matlab 求解2.1 最速下降法的循环语句在 Matlab 程序中，最速下降法是一个不断迭代循环的过程，且 确保同时更新，具体循环如下：01和001temp=-J： （ ， ）101t-： （ ， ）01=tepm：2.2 ComputeCost函数的定义J 函数即目标函数，m 代表数据的个数，具体程序如下：function J = computeCost(X, y, theta)m = length(y);J = 0;predictions=X*theta;J=1/(2*m)*(predictions-y)*(predictions-y);end2.3 gradientDescent函数的定义GradientDescent 函数就是最速下降法的迭代与循环过程 5，具体如下：function theta, J_history = gradientDescent(X, y, theta, alpha, num_iters)m = length(y);J_history = zeros(num_iters, 1);for iter = 1:num_iterstemp1 = theta(1) - (alpha / m) * sum(X * theta - y).* X(:,1);temp2 = theta(2) - (alpha / m) * sum(X * theta - y).* X(:,2);theta(1) = temp1;theta(2) = temp2;J_history(iter) = computeCost(X, y, theta);endend2.4 散点图的绘制具体的数据已经根据一定格式记录在 txt 文件中，因此只需调取其中的数据并将散点绘制到图中 6，具体程序如下：function plotData(x, y)figure; data = load(ex1data1.txt);X = data( : , 1 );Y = data( : , 2);X = ones(size(X,1),1),X; - 4 -plot(X,Y,rx,MarkerSize, 4); axis(4 24 -5 25); xlabel(x); ylabel(y); end运行程序后得到散点图如图 2 所示：图 2 散点图Fig.2 scatter diagram2.5 回归线的绘制通过机器统计学习后得到线性回归线，如图 3 所示： 图 3 回归线fig.3 The regression line运行主程序中的语句计算出 和 J 的最小值，程序如下：01和- 5 -theta,J_history = gradientDescent(X, y, theta, alpha, num_iters);Matlab 程序运行后结果报告如图 4 所示：图 4 报告窗口fig.4 The report window其中 ans 结果为 Matlab 运行的最后一步结果，下面为模型参数值。假设模型参数值： ，89571.30930.2.6 等高线及三维图的绘制将 和 J 的值绘制在三维图上，将 的值以 J 为高度绘制绘制成等高线图，并将中心点即01和 01和优化的参数用红色点标记，如下图 5 所示：图 5 三维图和等高线图Fig. 5 Graphic model and contour map等高线图中点红点代表了目标函数最小是对应的模型参数值，该参数值下的直线到每一个统计数据点的距离最短，最能代表离散统计数据的线性回归特性。2.7 结论根据已经给定的城市人口数目与饮品连锁店利润的关系，通过建立目标函数，由 Matlab 软件进行迭代求解，求解后的线性回归模型函数如下：(4)xxh19304.857.3)(3 实际问题中国人多地少，土地不可再生，土地的稀缺性导致土地价格在未来城市发展中看涨，因此，房价也自然总体趋涨 7。对于一些在沪求学的异地研究生们，在毕业后会存在是否留沪的选择。而房价是我们选择留沪的主要因素。为了了解上海的住房面积和价格的情况，我从链家网随机查找了 120 组不同区域以及不同面积的售房交易数据，将这些数据进行整理，并且希望通过统计学习的手段得出一般住房面积和价格的线性回归关系，通过设定两个期望住房面积来估计对应面积下的购房价格。- 6 -3.1 程序的设计首先，通过 EXCEL 表格将收集的数据进行整理，得出了上海地区售房面积和价格对应表格，120 组数据具体可见于附录。将 120 组数据通过 TXT 文件进行编辑后导入到 Matlab 软件中，其部分数据在软件中的表现形式如图 6 所示：图 6 售房面积和价格数据在 Matlab中表现形式Fig. 6 The data expression form of house area and price on sale in Matlab将这些数据以 txt 的格式保存，文件名为 house，以便在主程序加载数据调用。该最速下降统计学习程序主要分为调用数据、绘制离散图、梯度下降和拟合直线、可视化显示目标函数、目标函数值与迭代次数关系函数这五个部分。3.1.1 调用数据程序将 house.txt 文件中的数据加载到 data 矩阵中，将 data 矩阵中的第一列中所有值组成一个向量储存在X 向量中，将 data 矩阵中的第二列中所有值组成一个向量储存在 y 向量中，然后在 y 向量中计算训练样本数量存储到变量 m 中。其相关程序段如下：data = load(house.txt); %调用数据X = data(:, 1); y = data(:, 2); %向量存储m = length(y); %计算训练样本数量3.1.2 绘制离散图程序该过程主要包括主程序的调用以及调用函数的定义，定义函数首先要打开一个绘图窗口，然后定义数值范围，最终附上标签绘制图形或者点。通过这样的方式处理，我们可以清晰的看到上海售房面积和价格的离散分布关系。主程序部分：plotData(X, y);fprintf(Program paused.Press enter to continue.n);定义调用函数 plotData(X, y)function plotData(x, y)figure; % open a new figure windowplot(x, y, rx, MarkerSize,10);ylabel(Shanghai Commercial Housing Price in RMB one million);xlabel(Shanghai Commercial Housing Area in 10 square meters);- 7 -end3.1.3 梯度下降和拟合直线程序此处先对参数进行初始化设置，然后设置初始的迭代次数和计算步长（学习速率） ，接着进行目标函数的迭代计算，并且不断将模型参数显示在屏幕上，接着画出拟合的直线，在屏幕上显示显示“数据训练和直线拟合” ，最后根据拟合的直线估计 80和 160的商品房所对应的价格。主程序部分：fprintf(Running Gradient Descent .n)X = ones(m, 1), data(:,1); % Add a column of ones to xtheta = zeros(2, 1); % initialize fitting parameters iterations = 1000; % Some gradient descent settingsalpha =0.01; % compute and display initial costcomputeCost(X, y, theta) % run gradient descenttheta = gradientDescent(X, y, theta, alpha, iterations); % print theta to screenfprintf(Theta found by gradient descent: );fprintf(%f %f n, theta(1), theta(2); % Plot the linear fithold on; % keep previous plot visibleplot(X(:,2), X*theta, -)legend(Training data, Linear regression)hold off % dont overlay any more plots on this figurepredict1 = 1,8 *theta; % Predict values for house area of 35,000 and 70,000fprintf(For house area = 80, we predict a price of %fn,.predict1*1000000);predict2 = 1, 16 * theta;fprintf(For house area = 160, we predict a price of %fn,.predict2*1000000);fprintf(Program paused. Press enter to continue.n);pause;定义梯度下降调用函数：function theta, J_history = gradientDescent(X, y, theta, alpha, num_iters)m = length(y); % number of training examplesJ_history = zeros(num_iters, 1);for iter = 1:num_iterstemp1 = theta(1) - (alpha / m) * sum(X * theta - y).* X(:,1);%最速下降发的迭代循环 temp2 = theta(2) - (alpha / m) * sum(X * theta - y).* X(:,2);theta(1) = temp1;theta(2) = temp2;%同时更新 theta1 和 theta2J_history(iter) = computeCost(X, y, theta);endend3.1.4 可视化显示目标函数该处 Matlab 程序主要将目标函数和两模型参数的对应关系通过空间平面和等高线的形式表达出来，根据空间平面我们可以清晰看出目标函数的具体极值分布形式。根据等高线图可以看出最终迭代计算结果是否为最优解，倘若红色的点位于等高线的最中间，则该情况为最优解，否则，不是最优解。主程序部分：- 8 -fprintf(Visualizing J(theta_0, theta_1) .n)theta0_vals = linspace(-10, 10, 100);theta1_vals = linspace(-1, 4, 100);J_vals = zeros(length(theta0_vals), length(theta1_vals);for i = 1:length(theta0_vals)for j = 1:length(theta1_vals)t = theta0_vals(i); theta1_vals(j); J_vals(i,j) = computeCost(X, y, t);endendJ_vals = J_vals;figure;surf(theta0_vals, theta1_vals, J_vals)xlabel(theta_0); ylabel(theta_1);figure;contour(theta0_vals, theta1_vals, J_vals, logspace(-2, 3, 20)xlabel(theta_0); ylabel(theta_1);hold on;plot(theta(1), theta(2), rx, MarkerSize, 10, LineWidth, 10);定义目标函数：function J = computeCost(X, y, theta)m = length(y); J = 0;predictions=X*theta;J=1/(2*m)*(predictions-y)*(predictions-y);end3.2 确定计算步长我们在这里迭代计算时，并不是一开始就能迭代成功，正如前面所说步长过小，会延长搜索的时间，步长过大，可能会错过极值点，导致不收敛。故我们可通过建立目标函数值和迭代次数的关系函数，将不同的步长以一定的数量级排列，分别作出目标函数值和已运算迭代次数的关系函数曲线。为了更好的展示目标函数值和已运算迭代次数的对应关系，可以设置足够的迭代次数，设定足够的迭代次数为3000。目标函数值与迭代次数关系函数程序段：iterations=3000;alpha=1;%输入第一个步长theta,J1 = gradientDescent(X, y, zeros(2,1), alpha, iterations);alpha=2;%输入第二个步长theta,J2 = gradientDescent(X, y, zeros(2,1), alpha, iterations);alpha=3;%输入第三个步长theta,J3 = gradientDescent(X, y, zeros(2,1), alpha, iterations);figure;plot(1:numel(J1), J1, -k, LineWidth, 2);%第一条线显示连续线hold on;plot(1:numel(J2), J2, - -k, LineWidth, 2);%第二条线显示虚线- 9 -hold on;plot(1:numel(J3), J3, -.k, LineWidth, 2);%第三条线显示点画线axis(0 3000 0 15);%限定坐标值显示范围xlabel(迭代次数); ylabel(J 的值 );legend(alpha=1 , alpha=2, alpha=3)hold off现在选取步长 为 1、0.1、0.01 和 0.001、0.0001、0.00001 两组数据进行迭代运算，通过 Matlab作出目标函数值和已运算迭代次数的关系函数曲线，判断目标函数值的收敛性。 3.2.1 第一组迭代运算判断敛散性将第一组步长值 1、0.1 和 0.01 带入目标函数值与迭代次数关系函数程序段中进行运算作图，为了让图形表示清晰，其坐标值可再一次限定其范围，其运行结果如图 7 所示。根据图中结果可知：步长为 1和 0.1 时，J 的函数是发散的，而步长为 0.01 时，J 的函数是收敛的。图 7 第一组步长中目标函数值的敛散性Fig.7 The convergence curves of first group alpha 3.2.2 第二组迭代运算判断敛散性将第一组步长值 0.001、0.0001 和 0.00001 带入目标函数值与迭代次数关系函数程序段中进行运算作图，为了让图形表示清晰，其坐标值可再一次限定其范围，其运行结果如图 8 所示。根据图中结果可知：在该组所有步长下，J 的函数都是收敛的。-