多项式拟合——正规方程法

c++ 专栏收录该内容
1 篇文章 0 订阅

前言

这几天师哥给我布置多项式拟合的作业,可我只是一个大一菜鸡,只会一点点c++(T_T)胖虎好帅
在b站看了几天吴老大(吴恩达)的机器学习视频,好不容易打出了梯度下降结果
在这里插入图片描述

连我的电脑都觉得南

在这里插入图片描述

面对精度问题手足无措的我选择投向正规方程法的怀抱

正文开始 如果看不懂那就去看吴恩达机器学习的视频吧

一、一个方程

在这里插入图片描述
其中θ为参数矩阵,X为训练样本自变量的矩阵,Y为训练样本的答案矩阵
有了这个方程,我们可以直接求解出多项式拟合的参数
为了实现这个方程,我们得写出以下几个函数

  1. 矩阵乘法
  2. 矩阵转置
  3. 矩阵求逆
    我自定义了一个结构体
struct ma
{
 double matrix[nn][mm];
 int line,row;
};

line与row代表矩阵的行和列

二、算法实现

1、矩阵乘法

一个大小为N * U 的矩阵A与一个U * M的矩阵B相乘,会得到一个N*M的矩阵C
C内元素与A、B内元素的关系为:
在这里插入图片描述
看不懂的小同学补补线性代数就可以了
因此,可以得到

for(int i=0;i<ll;i++)
  for(int j=0;j<rr;j++)
   for(int k=0;k<uu;k++)
    C.matrix[i][j]+=A.matrix[i][k]*B.matrix[k][j];

为了方便,使用重载运算符

ma operator * (ma A,ma B)//矩阵乘法 
{
 int ll=A.line,rr=B.row,uu=A.row;
 ma C;memset(C.matrix,0,sizeof(C.matrix));//初始化矩阵内元素
 C.line=ll;C.row=rr;//初始化C的大小
 for(int i=0;i<ll;i++)
  for(int j=0;j<rr;j++)
   for(int k=0;k<uu;k++)
    C.matrix[i][j]+=A.matrix[i][k]*B.matrix[k][j];
 return C;
}

要算出矩阵C,只要写C=A*B真的巨方便(自从我知道了这个东西就再也离不开它了)(●’◡’●)

2、矩阵转置

啥是矩阵转置?
在这里插入图片描述
就这个样子啦,代码也非常简单啦

ma transposition(ma A)//转置矩阵 
{
 int ll=A.line,rr=A.row;
 ma C;C.line=rr;C.row=ll;
 for(int i=0;i<ll;i++)
  for(int j=0;j<rr;j++)
   C.matrix[j][i]=A.matrix[i][j];
 return C;
}

3、矩阵求逆

因为我当时刚上线性代数没多久,没教矩阵求逆,然后就百度,发现伴随矩阵法是能用c++实现的
而且XT*X得到的肯定是方阵,满足伴随矩阵法的使用条件
里面有讲伴随矩阵法求逆矩阵自己翻翻看
求伴随矩阵,得求代数余子式。做到这个部分,我忽然想起,之前闲着无聊打了个程序来求解n阶行列式(会输出运算过程)来做线性代数作业(顺便装装B嘻嘻),所以赶紧翻出之前的代码复制粘贴哈哈!!代码如下

double determinant (ma A)//计算行列式 
{
 int ll=A.line,num;
 if(ll==2)return (A.matrix[0][0]*A.matrix[1][1]-A.matrix[1][0]*A.matrix[0][1]);
 double ans=0;
 for(int i=0;i<ll;i++)
 {
  ma B;
  for(int j=1;j<ll;j++)
  {
   num=0;
   for(int k=0;k<ll;k++)//第几列 
   {
    if(k==i)continue;
     B.matrix[j-1][num++]=A.matrix[j][k];
   }
  }
  B.line=ll-1;
  if(i&1)
   ans-=A.matrix[0][i]*determinant(B);
  else 
   ans+=A.matrix[0][i]*determinant(B);
 }
 return ans;
}

这个别问我怎么搞的,递归实现,你手算怎么算的,就让电脑怎么算(才不会告诉你当初为了装这波B打了半个小时才解了作业里的两个四阶行列式)
知道怎么求解n阶行列式后,要求逆矩阵就比较容易了,直接模拟就行了,只是这个复杂度爆棚,四重循环,最内层循环还调用了inserve
(神威太湖之光表示这个程序它跑不动┭┮﹏┭┮)
吴老大说逆矩阵只要O(n^3)就可以解决,所以大家有更好的解法就不要用我的了哈
(说什么傻话,大家都跑去用matlab、octave、python了つ﹏⊂)
代码如下

ma inverse(ma A)//求逆矩阵 
{
 int ll=A.line;
 double und;und=determinant(A);
 ma B,C;B.line=ll-1;B.row=ll-1;C.line=ll;C.row=ll;
 int hang,lie;
 for(int i=0;i<ll;i++)
 {
  for(int j=0;j<ll;j++)
  {
   hang=0;lie=0;
   for(int k=0;k<ll;k++)
   {
    if(k==i)continue;
    for(int u=0;u<ll;u++)
    {
     if(u==j)continue;
     B.matrix[hang][lie++]=A.matrix[k][u];
    }
    lie=0;hang++;
   }
   C.matrix[i][j]=determinant(B);
   if((i+j)%2)C.matrix[i][j]*=-1;
   C.matrix[i][j]/=und;
  }
   
 }
 C=transposition(C);
 return C;
  
}

大概就讲到这里了吧,第一次写博客,以上是我为了完成作业,了解了概念就开始自己瞎几把搞的结果(?),有什么错误的欢迎大佬们指正我,希望能帮到大家哈哈!!!

  • 1
    点赞
  • 0
    评论
  • 3
    收藏
  • 打赏
    打赏
  • 扫一扫,分享海报

©️2022 CSDN 皮肤主题:数字20 设计师:CSDN官方博客 返回首页

打赏作者

InSan1r

你的鼓励将是我创作的最大动力

¥2 ¥4 ¥6 ¥10 ¥20
输入1-500的整数
余额支付 (余额:-- )
扫码支付
扫码支付:¥2
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、C币套餐、付费专栏及课程。

余额充值