- 拟合用到的数据是老师给的鲍鱼数据<abalone.csv>,想做一个“整体重量”关于“长度”、“直径”、“高度”的多元回归分析。下面是多元回归的基础代码:
clc;clear
A1=importdata('abalone.csv');
% A1.data(101:end,:)=[]; %% 弄十几个数据看一眼
X1=A1.data(:,1:3); %% 三个自变量:length,diametar,heigth
y=A1.data(:,4); %% 一个因变量:whole weigth
X=[ones(length(y),1),X1(:,1).^2,X1(:,2).^2.5,X1(:,3)]; %%点乘构造高次项
%%研究过图像,比较符合"Y=a*X1^2+b*X2^2.5+X3"的形状
[b, bint, r, rint, stats]=regress(y,X);
stats; %%stats,最关键的度量
rcoplot(r,rint) %%刻画目前data的拟合残差
- 下面贴原来的残差图
- 下面打算删除残差不合格的点,思考了半天构造了一个删除函数(untitled.m),功能是将因变量y,自变量矩阵X放进去,以残差为“canshu”为阈值,删除不合格点,返回删除后的y,X矩阵。
function [laterX,latery]=untitled(y,X,canshu)
[~, ~, r, ~, ~]=regress(y,X);
[chang,~]=size(X);
for i=1:chang %%数据的指针删除code
if (abs(r(i))>canshu) %%绝对值的设定也是无奈之举,看上面的残差图,不合格点在0点上下都有
X(i,:)=[];
y(i,:)=[]; %%这两句删除
laterX=X(:,:);
latery=y(:,:); %%这两句提取结果
chang=chang-1;
end
if(i>=chang)
break %%指针操作防止 i 溢出
end
end
end
- 构造了空矩阵box存放观察y函数参数b的所有变化
number=1;
box=zeros(700,7);
- 这个是主函数里面调用function删除残差的
for i=0.01:0.001:0.7 %%0.01-0.7,由上面残差图设置的,超过0.9就不合理了,步长0.001还是有点误差的
[boxX1,boxy]=untitled(y,X,i); %%boxX1和boxy存放去噪后的数据
[b1,~, ~, ~,~]=regress(boxy,boxX1);
lost=((b1(1)-b(1))^2+(b1(2)-b(2))^2*(boxX1(2).^2)...
+(b1(3)-b(3))^2*(boxX1(3).^2.5)+(b1(4)-b(4))^2*boxX1(4)).^2;
%%构造损失函数lost,经验风险最小化,未考虑结构最小化,也就是没考虑过拟合的问题。
%%因为第一轮box观察了去噪后的stats值,从数理分析的角度都是合适的,mse与R2都在允许范围内
%%所以第二轮只观察参数b的变化,来使折线图更加逼近真实,具体思辨过程请看我前一篇文章
ZZH=log(10000*lost); %%log正则化,因为lost有10^8阶的变化
box(number,1)=i; %%box放参数b开始规划
box(number,2:5)=b1;
box(number,6)=ZZH;
number=number+1;
end
- 下面就是为了求lost的min
figure(50) %%看ZZH图像
plot(1:length(box(:,5)),box(:,6),'r');
- 贴ZZH图,可以看出有一个最小值
[h,hh]=min(box(:,6)) %%用min函数求ZZH的min,因为ZZH的min也就是lost的min,及想求的lost的min
%%h,hh存放这个min值时的阈值,求最优阈值下的y,X矩阵,再次拟合运算及为所求的最优参数b11
[boxX11,boxy1]=untitled(y,X,box(hh,1));
[b11,~, r11, rint11,~]=regress(boxy1,boxX11);
figure(60)
Z0=b(1)+b(2)*boxX11(:,1).^2+b(3)*boxX11(:,2).^2.5+b(4)*boxX11(:,3);
%%未去噪时的“拟合—原数据对比图”
Z=b11(1)+b11(2)*boxX11(:,1).^2+b11(3)*boxX11(:,2).^2.5+b11(4)*boxX11(:,3);
%%去噪后,求最优参数后的 “拟合—原数据对比图”
plot(boxX11,Z,'b',boxX11,Z0,'y')
figure(70)
rcoplot(r11,rint11)
%%去噪后,求最优参数后的残差图
- 下面贴图对比
(两个函数我设置了两种颜色的,看不出来是因为拟合的太好而盖住了)
感觉残差不合格点少了一点哈,想完全删除是不可能的,效果图的大致覆盖说明了这一次拟合是非常好的,做到了拟合现在,预测未来。
(希望我的代码能帮到大家,由于我还是个小白,code有错,设计思路,数理分析有问题也很正常,欢迎大家提出来。有错误才有进步)
最后凑不要脸的求一个赞hhhh
- (2020/4/30)更新:
上午老师看了我的博客给我提了几点意见
1.模型的级数是怎么看出来的,2.数据是不能随便删除的
我想了一下确实。首先讲第二点:
data不能随便删除,以为每个数据都有存在的意义。我的拟合优度已经达到0.94了已经非常好了,所以data的额外删除是没必要的,可以容忍。你把部分data删除而留下自己需要的数据,强求模型的优度而这样做,不是鸵鸟行为么?
而如果拟合优度R2过小时,可以考虑删除误差比较大的data来达到比较好的拟合结果。我写的delete函数就有用武之地了。
然后讲第一点:
X1适合平方项,X2我给了2.5次方,X3当成了直线来拟合。
而用直线拟合Y=a*X+b来拟合时优度只有0.81,还是比较合适的。