一、实验内容
- 根据梯度下降法、闭式解完成一元线性回归实验，并比较两种解下的实验结果。
划分训练集和测试集比例(1:2 划分）。通过 randperm 函数得到一个打乱序列的 b，该序列包含（1，n），本文中 n 为合并后矩阵的行数。由于题目要求的 1:2 划分容易出现不完美分割的情况，因此对划分的数量进行取整操作。然后 c,d 两个序列分别取打乱序列的前 33%、后 66%，根据 c,d 两个序列的数值取原矩阵对应的列数，最后将按顺序将每次划分结果对应保存于对应的训练集和测试集中，并分别绘制训练集和测试集的散点图。
遍历测试集与训练集，计算两者两两之间的欧式距离。绘图展示梯度下降线性回归，最小二乘线性回归的图像。
- 根据训练数据，学习投影矩阵。然后将 LDA 在训练样本上的低维表示结果可视化。最后使用距离最短对测试样本进行分类。
首先定义线性判别回归函数 LDA。其中，mappedX 为低维数据降维后的坐标，mapping 为映射的信息，M 为，X 为样本矩阵，label 为样本标签，no_dims: 投影空间(低维空间)的维度。
初始化。①设置特征默认维数为 2，提高函数的复用性。然后通过②bsxfun 函数对函数进行标准化，优点为当 A 和 B 的维度不相等，并且 A 和 B各自有一个维度为 1 时，函数自动调整维度，目的也是提高函数复用性。之后③去掉矩阵标签中重复的元素，将标签类别保存于 unique 中。最后初始化用到的矩阵。
通过①计算总协方差矩阵（全局散度矩阵）St = np.dot((X - mean).T, X - mean)，②St = np.dot((X - mean).T, X - mean)，③全局散度矩阵 St=Sw+Sb，所以 Sb = St-Sw 三个式子计算类内离散度矩阵 Sw。同时确保降维后的特征点维度小于 max，提高代码的安全性。
按降序排列特征值和特征向量。M 变换矩阵由 v 的最大的 nc-1 个特征值所对应的特征向量构成。对 inv(Sw)*Sb 进行特征分解[V,D]=eig(A,B) eig(A,B)返回方阵 A 和 B 的 N 个广义特征值，构成 N×N 阶对角阵 D， 其对角线上的N 个元素即为相应的广义特征值，同时将返回相应的特征向量构成 N×N 阶满秩矩阵，且满足 AV=BVD。之后判断数组的元素是否是 NaN，并提取前no_dims 个特征向量。
最后计算投影后的（二维）坐标值。
主函数部分，先通过 xlsread 或 csread 函数读取 csv 文件数据。
在 LDA 多分类后，将返回的四维降维成的二维坐标值反应到图像上，通过不同颜色和形状区分不同的标签，计算坐标两两之间的欧氏距离。
二、实验结果
将其改为 length 函数。
梯度下降得到的 theta 的值:
w:1.159369,
b:-3.369070
梯度下降法一元回归的欧式距离为 591.925823,
闭式解为 594.443011
可以得到结果，梯度下降和最小二乘线性回归的直线几乎重合，差别不大。
训练集[日本，中国]的标签分别为 1,2。
结果很合适，可以看到点就在散点图聚集处。
三、实验源码
实验六：梯度下降一元回归
% 1.根据梯度下降法完成一元线性回归实验。 % 2.根据闭式解完成一元线性回归实验。 % 3.比较两种解下的实验结果。 % 定义多元线性回归计算代价损失函数 % 参数进行线性回归拟合X和y中的数据点 function J = computeCost(X, y, theta) %计算使用theta作为的成本 m = length(y); % number of training examples J = 0; % 返回正确的变量 % J <= 10e-5; %保证误差小于10的-5次方 predictions = X * theta; J = 1/(2*m)*(predictions - y)'*(predictions - y); % J = 1/(2*m)*(predictions - y).^2; %线性函数的求解公式 end %定义一个实现梯度下降的函数 function [theta, J_history] = gradientDescent(X, y, theta, alpha, num_iters) m = length(y);%取数据的长度 J_history = zeros(num_iters, 1);%定义J_history为num_iters行1列的向量 %其中，num_iters是在调用该梯度下降函数时的参数 for iter = 1:num_iters S = (1 / m) * (X' * (X * theta - y));%相当于求导 theta = theta - alpha .* S;%theta的更新 J_history(iter) = computeCost(X, y, theta);%在这里调用computeCost(X, y, theta) %相当于在慢慢的减小代价函数 end end % 主函数 % 1.根据梯度下降法完成一元线性回归实验。 % 2.根据闭式解完成一元线性回归实验。 % 3.比较两种解下的实验结果。 clc;clear; close all; %% 数据导入 data = load('D:\Desktop\大三上\神经网络\实验6\ex1data1.txt'); %读入文件数据 x=data(:,1); %取data中第一列，人口每行一个样本 y=data(:,2); %取data中第二列，利润每行一个样本 % plotData(x,y); %显示样本数据图 %% 初始化值 m=length(x); %取x的长度，即数据的个数 x=[ones(m,1),data(:,1)]; %将原来m*1的x，多拓展了一个元素列 %全为1的一列，变为m*2的矩阵，其中的一列为全1列，主要是为了方便后面的运算 theta=zeros(2,1); %给theta先赋值为0，此theta相当于theta0和theta1 %组成的2*1的向量 iterations=1500; %设迭代次数为1500 alpha=0.01; %设梯度下降步长为0.01 %% 划分训练集和测试集比例(1:2划分） b=randperm(length(x(:,1))); %打乱列序列 e=round(length(x(:,1))*1/3); c=b(1:e); %取打乱序列的前33% d=b(e+1:end); %取打乱序列的后66% %end直到数组结尾 atrainx=x(c,:); %(:,:)为取原矩阵行和列，(:,1:50)为取原矩阵列，行取前33% atrainy=y(c,:); atestx=x(d,:); %(:,d)为取原矩阵列，行取剩余的随机66%的行 atesty=y(d,:); plot(atrainx(:,2),atrainy,'ro','MarkerSize',10) %显示了训练数据样本数据散点图 hold on; %图保留 plot(atestx(:,2),atesty,'yo','MarkerSize',10) %显示了测试数据样本数据散点图 hold on; %图保留 %% 梯度下降 theta=gradientDescent(atrainx,atrainy,theta,alpha,iterations); %调用梯度 %下降更新theta % matlab索引从1开始 fprintf('梯度下降得到的theta的值:\nw:%f,\nb:%f\n',theta(2),theta(1)) %theta(2)为b %% 闭式解 p = polyfit(atrainx(:,2),atrainy,1); %% 计算测试集与训练集两者之间的欧式距离 predicty=zeros(2,length(atesty)); %初始化两列预测y的数据集 % for i=1:length(d) for i=1:d predicty(1,i)=[1,i]*theta; %第一列预测y为梯度下降算法求得的解 predicty(2,i)=polyval(p,i); %第一列预测y为闭式解求得的解 % predicty(2,i)=i*p(1)+p(2); %polyval(p,i)就是i*p(1)+p(2)的意思 end distance1=sqrt(sum(abs(predicty(1,:)-atesty)).^2); %sum(x.^2)计算数组的平方和 distance2=sqrt(sum(abs(predicty(2,:)-atesty)).^2); %sum(x.^2)计算数组的平方和 fprintf('梯度下降法一元回归的欧式距离为 %f,\n闭式解为 %f\n',mean(distance1),mean(distance2)); %% 绘图展示结果 hold on; %保持图像不被刷新 plot(atrainx(:,2),atrainx*theta,'b-','MarkerSize',10) %显示了梯度下降拟合的直线 hold on; %保持图像不被刷新 plot(atrainx(:,2),atrainx(:,2)*p(1)+p(2),'b-','MarkerSize',10) %显示了最小二乘拟合的直线 legend('训练数据','测试数据','梯度下降线性回归','最小二乘线性回归') %图例标注
实验七：LDA一元回归
clear;clc %训练数据集数据和标签 %% 读取路径下表一数据 %方法一：xlsread % M = xlsread('D:\Desktop\实验七\EduRecivedFiles\LDA数据.csv',1,'C2:G23'); %方法二：csvread % 文件名，开始读取位置的行号（row）和列号（col）注意从0开始 % range = [R1 C1 R2 C2]，这里（R1，C1）是读取区域的左上角，（R2，C2）是读取区域的右下角 M1 = csvread('D:\Desktop\实验七\EduRecivedFiles\LDA数据.csv', 1, 2, [1 2 22 6]); x=M1(:,1:end-1); y=M1(:,end); %% LDA多分类 [mappedX,mapping,W]=FisherLDA2(x,y,2); %% 结果可视化 figure temp1=find(y==1); plot(mappedX(temp1,1),mappedX(temp1,2),'r^','MarkerSize',18) hold on temp2=find(y==2); plot(mappedX(temp2,1),mappedX(temp2,2),'bo','MarkerSize',18) hold on temp3=find(y==3); plot(mappedX(temp3,1),mappedX(temp3,2),'black*','MarkerSize',18) hold on % Z=-4:5; % plot(Z,Z*W,'-k') % axis([-4 5 -4 2]); %% 测试样本分类 % 绘制测试样本 temp0=find(y==0); plot(mappedX(21,1),mappedX(21,2),'y.','MarkerSize',18) plot(mappedX(22,1),mappedX(22,2),'green.','MarkerSize',18) legend('标签1','标签2','标签3','测试样本1日本','测试样本2中国'); hold on % 计算两者之间的欧式距离 tempdis1=inf;tempdis2=inf; %初始化最短距离 tempj=0;%初始化最短距离点的坐标 for i=1:length(mappedX)-2 distance1=sqrt(sum(abs(mappedX(i,:)-mappedX(21,:))).^2); %sum(x.^2)计算数组的平方和 distance2=sqrt(sum(abs(mappedX(i,:)-mappedX(22,:))).^2); %sum(x.^2)计算数组的平方和 if(distance1<tempdis1) tempdis1=distance1;tempj1=i; % 如果距离最小则判断为同一类型 res1=y(i,:); % 将最短距离的训练集标签赋值给测试集 end if(distance2<tempdis2) tempdis2=distance2;tempj2=i; % 如果距离最小则判断为同一类型 res2=y(i,:); % 将最短距离的训练集标签赋值给测试集 end end fprintf('训练集[日本，中国]的标签分别为%d,%d\n',res1,res2) % 绘制一维图像() % [mappedX1,mapping1,W1]=FisherLDA2(x,y,1); % Z=-4:5; % plot(Z,Z*W(1),'-k') % axis([-4 5 -4 2]); % hold on % legend('标签1','标签2','标签3','测试样本','一维直线'); function [mappedX, mapping, M] = FisherLDA2(X, labels, no_dims) % X: n*d matrix,each row is a sample; % labels: n*k matrix,each is the class label % no_dims: 投影空间(低维空间)的维度(no_dimms >= 1, default = 2，max = len(labels)-1 % mappedX： 低维数据降维后的坐标 % mapping： 映射的信息 %% 初始化 % 设置特征点默认的维数为2 if ~exist('no_dims', 'var') || isempty(no_dims) no_dims = 2; end % 标准化数据集 mapping.mean = mean(X, 1); % 当A和B的维度不相等，并且A和B各自有一个维度为1时，函数自动调整维度 X = bsxfun(@minus, X, mapping.mean); % 减法操作 % 去掉矩阵标签中重复的元素 [classes, bar, labels] = unique(labels); nc = length(classes); % 初始化类内离散度矩阵 Sw Sw = zeros(size(X, 2), size(X, 2)); % 计算总协方差矩阵（全局散度矩阵）St = np.dot((X - mean).T, X - mean) St = cov(X); %% 全局散度矩阵St=Sw+Sb，so Sb = St-Sw % St = np.dot((X - mean).T, X - mean) % 计算类内离散度矩阵Sw for i=1:nc cur_X = X(labels == i,:); % 更新类内离散度矩阵Sw C = cov(cur_X); p = size(cur_X, 1) / (length(labels) - 1); Sw = Sw + (p * C); end % 计算类间离散度矩阵Sb = St-Sw Sb = St - Sw; Sb(isnan(Sb)) = 0; Sw(isnan(Sw)) = 0; Sb(isinf(Sb)) = 0; Sw(isinf(Sw)) = 0; % 确保嵌入特征空间中特征点的维数小于max if nc <= no_dims no_dims = nc - 1; warning(['目标维度降为： ' num2str(no_dims) ]); end %% M变换矩阵由v的最大的nc-1个特征值所对应的特征向量构成 % 对inv(Sw)*Sb进行特征分解[V,D]=eig(A,B) % eig(A,B)返回方阵A和B的N个广义特征值，构成N×N阶对角阵D， % 其对角线上的N个元素即为相应的广义特征值，同时将返回相应的特征向量构成N×N阶满秩矩阵，且满足AV=BVD。 [M, lambda] = eig(Sb, Sw); % 按降序排列特征值和特征向量 lambda(isnan(lambda)) = 0; % 判断数组的元素是否是NaN [lambda, ind] = sort(diag(lambda), 'descend'); %diag(lambda)提取lambda矩阵的对角线 M = M(:,ind(1:min([no_dims size(M, 2)]))); %提取前no_dims个特征向量，ind——最小值下标 %% 计算投影后的中心值 mappedX = X * M; % 存储测试集投影中心值 mapping.M = M; mapping.val = lambda;
