【蒟蒻科技】肝信息预测年龄

Author： jvruo
发布时间：August 18, 2020
5113views
4 comments
8804字数
Categories：蒟蒻科技

在上学期的网课的学习中，我参加了本校的数据分析实战课程，作为计算机学院的一名大一新生，这是我第一次接触数据分析相关内容项目。下面让我们来看看这次的预测小项目吧。

问题描述

在这次的预测中，我们需要通过134个人的年龄、性别和15938个基因的表达值来推算给出147个人性别和相同的15938个基因的表达值的年龄。

其中给我们的数据一共是三个文件，也提供给大家下载交流学习：

liver_RNAseq_training.txt - 训练集，包含134个个体的年龄、性别和15938个基因的表达值，通过RNA测序进行量化。
liver_RNAseq_validation.txt - 验证/测试集，包含147个个体的性别和相同的15938个基因的表达值，通过RNA测序进行量化。
ExampleSubmission.csv - 正确格式的示例提交文件。

由于基因数量较多，我首先想到的就是通过筛选关键基因进行线性回归预测。

项目实现

数据预处理

在进行预测之前我发现训练集和测试集的数据差异极大，为了让预测更加精准，我讲所有数据的每种基因都进行了归一化处理，让其平均值等于1，这样就能进行下面的预测了。

线性回归

在之前的BLOG里我曾经讲到过线性回归模型的实现原理，这次我首先实现了C++代码实现了基于梯度下降法的第一版程序：

void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;//维数
    num = 10000000;//下降次数 

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = ans[i];
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
                //if(m == 0) printf("QwQ\n"); 
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            ans[i] = temp[i];
            //printf("%d ", ans[i]);
        } 
    }//梯度下降
    return ;
}

由于听老师说特征量越多，其误差就会越大，我开始也没多想，就直接选取了前20个基因表达来作为我的特征值，这样就得到了Kaggle上0.23315的结果。接着我又对特征值维度和α值及下降次数进行了调整，最后在8维Alpha = 3 * 1e-7时得到这种方法的最好成绩0.33037;

筛选数据的线性回归

接着我考虑对算法进行优化。由于PCA还不熟练，所以我首先考虑的是选取较为“显著”的数据，那么什么才算“显著”呢？在前面我已经提到过，为了让数据规模相近，我对所有基因都进行了归一化处理，那么我就对处理后的所有基因都求其方差，如果方差大，代表其差异大，其对年龄的贡献应该也是更加“显著”的。毕竟如果一个基因在每个人身上表现的都一样，那就没有必要研究它和年龄是否有关系了。通过这种思路，我求出了前八大显著的基因分别是第94、53、36、21、92、124、39和97号基因并对其在读入前进行标记：

memset(can, false, sizeof(can));
can[94] = can[53] = can[36] = can[21] = can[92] = true;
can[124] = can[39] = can[87] = true;

但是结果不是很好，只到了Kaggle上0.22252，所以我又放弃了这一方法。

带PCA的线性回归

过后我在课堂上听到了老师讲解PCA，结合之前写过的BLOG中PCA的原理，也自己去学习了C++下的PCA实现方法，并参考博客完成了PCA的相关函数代码：

void featurenormalize(MatrixXd &X) {
    //计算每一维度均值
    MatrixXd meanval = X.colwise().mean();
    RowVectorXd meanvecRow = meanval;
    //样本均值化为0
    X.rowwise() -= meanvecRow;
}
void computeCov(MatrixXd &X, MatrixXd &C) {
    //计算协方差矩阵C = XTX / n-1;
    C = X.adjoint() * X;
    C = C.array() / X.rows() - 1;
}
void computeEig(MatrixXd &C, MatrixXd &vec, MatrixXd &val) {
    //计算特征值和特征向量，使用selfadjont按照对阵矩阵的算法去计算，可以让产生的vec和val按照有序排列
    SelfAdjointEigenSolver<MatrixXd> eig(C);

    vec = eig.eigenvectors();
    val = eig.eigenvalues();
}
int computeDim(MatrixXd &val) {
    int dim;
    double sum = 0;
    for (int i = val.rows()-1; i >= 0; --i)
    {
        sum += val(i, 0);
        dim = i;

        if (sum / val.sum() >= 0.95)
            break;
    }
    return val.rows() - dim;
}

但是经过PCA压到10维后也只跑出了Kaggle上0.28324的成绩，所以我也放弃了这种做法。

带正则化的线性回归

最后通过课堂上老师和同学们的互动，我突然想起了由于维数的骤减可能过拟合的问题，我开始反思我只取8-10维的特征，会不会出现过拟合，所以我就学习了出力过拟合的最简单的方法——增加正则化系数，通过增加正则化系数这个惩罚项，来让所有系数的和尽量小，来尽量避免过拟合的情况，我也随之对我的梯度下降进行了相关修改：

double lambd = 1; 
void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;
    num = 8000000;

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = (1 - Alpha * lambd / (double)n) * ans[i];
        //正则化！！！！ 
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) ans[i] = temp[i];
    }//梯度下降
    return ;
}

就这样我跑出了我的个人最好成绩Kaggle上0.34173，直接冲进当时前十（不过后面掉了几名）只可惜到最后才开始尝试这个，各项系数还没有调整到最佳……

所以最终代码就是上面所提到的带正则化的线性回归啦，下为代码：

#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
using namespace std;
int n, m;
double ans[11000];//各项系数
double temp[11000], h[11000];//辅助数组
double x[1100][1100], y[11000];//训练集特征值和答案
double xx[1100] = {0}, yy[1100] = {0};//训练集和测试集各基因的平均值
const double Alpha = 3 * 1e-6;//单次梯度下降的幅度
void read() {//读入及相关操作
    n = 134;
    m = 10;

    for(int i = 1; i <= n; i++) x[i][0] = 1;
    char s[100];
    FILE *fin = fopen("liver_RNAseq_training.txt", "r");
    for(int k = 1; k <= 15940; k++) fscanf(fin, "%s", s);
    for(int i = 1; i <= n; i++) { 
        fscanf(fin, "%s", s);
        if(s[0] == 'M') x[i][1] = 1;
        fscanf(fin, "%lf", &y[i]);
        //printf("%lf ", y[i]);
        for(int j = 2; j <= m; j++) {
            fscanf(fin, "%lf", &x[i][j]);
            xx[j] += x[i][j] / double(n);
        }
        for(int k = 1; k <= 15940 - m - 1; k++) fscanf(fin, "%s", s);
    }//数据读入
    xx[1] = 1;
    for(int i = 1; i <= n; i++) {
        for(int j = 1; j <= m; j++) {
            x[i][j] /= xx[j];
        }
    }//数据均一化
    fclose(fin); 
    return;
}
double lambd = 1; //正则化系数
void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;//维数
    num = 8000000;//下降次数

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = (1 - Alpha * lambd / (double)n) * ans[i];
        //正则化！！！！ 
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) ans[i] = temp[i];//ans为系数
    }//梯度下降
    return ;
}

void write() {//预测及输出
    FILE *fin = fopen("liver_RNAseq_validation.txt", "r"), *fout = fopen("ans.txt", "w");
    int n = 147, m = 10;
    //for(int i = 0; i <= m; i++) printf("%lf ", ans[i]);
    char s[100];
    double read[200][200];
    yy[1] = 1;
    for(int k = 1; k <= 15940; k++) fscanf(fin, "%s", s);
    for(int i = 1; i <= n; i++) { 
        fscanf(fin, "%s", s);
        fscanf(fin, "%s", s);
        if(s[0] == 'M') read[i][1] = 1;
        read[i][0] = 1;
        for(int j = 2; j <= m; j++) {
            fscanf(fin, "%lf", &read[i][j]);
            yy[j] += read[i][j] / double(n);

        }
        //printf("\n");
        for(int k = 1; k <= 15940 - m - 1; k++) fscanf(fin, "%s", s);
    }//读入预测数据
    //for(int i = 1; i <= m; i++) printf("%lf ", yy[i]);
    for(int i = 1; i <= n; i++) {
        double x = 0;
        x = read[i][0] * ans[0];
        for(int j = 1; j <= m; j++) {
            read[i][j] /= yy[j];//正则化
            x += read[i][j] * ans[j]; //计算答案
        }
        fprintf(fout, "%0.0lf\n", x);//输出答案
    }
    fclose(fin); 
    fclose(fout);
    return;
}
int main() {

    read();
    work();
    write();
    return 0;
}

一些小问题

在我们的代码中，最后输出的是ans.txt，但是答案要求输出的是.csc结构的，所以我们还需要把txt格式的文件转换成csv，由于我输出的是一列答案，第i行表示的是第i个人年龄的预测值。所以我还需要一个程序来进行转换：

#include <fstream>
#include <cstdio>
using namespace std;
int ans[200];
int main()
{        
    FILE *fin = fopen("ans.txt", "r");
    for(int i = 1; i <= 147; i++) fscanf(fin, "%d", &ans[i]);
    fclose(fin);
    ofstream oFile;

        //打开要输出的文件
    oFile.open("ans.csv", ios::out | ios::trunc);
    oFile << "ID" << "," << "Age"  << endl;
    for(int i = 1; i <= 147; i++) {
        oFile << "test" << i << "," << ans[i]  << endl;
    }

    oFile.close();
    return 0;
}

个人总结

这是我第一次参与数据分析相关项目，在工程的最开始自己瞎尝试太多天（比如多层神经网络）才开始导致到后来尝试机会有点不够用……虽然最后还是没有到达我算法的最优情况，也没有使用更加高级的算法，但是这一次从零开始学习统计学习算法完成编程的过程确实让我收获很多。最后希望您能喜欢我这个预测项目！喜欢这篇BLOG！

Last modification：September 5th, 2021 at 09:15 pm

If you think my article is useful to you, please feel free to appreciate

4 comments

ls
November 19th, 2020 at 05:47 am

大一？？？tql，xd中学搞竞赛的吗

Reply
1. jvruoAuthor
  February 14th, 2021 at 05:20 pm
  
  @ls
  嗯嗯，中学有接触OI这一块的
  
  Reply
木薯
October 16th, 2020 at 03:04 am

tql

Reply
1. jvruoAuthor
  February 14th, 2021 at 05:21 pm
  
  @木薯
  木薯tql
  
  Reply

【蒟蒻科技】肝信息预测年龄

jvruo • 2020 年 08 月 18 日

在上学期的网课的学习中，我参加了本校的数据分析实战课程，作为计算机学院的一名大一新生，这是我第一次接触数据分析相关内容项目。下面让我们来看看这次的预测小项目吧。

问题描述

在这次的预测中，我们需要通过134个人的年龄、性别和15938个基因的表达值来推算给出147个人性别和相同的15938个基因的表达值的年龄。

其中给我们的数据一共是三个文件，也提供给大家下载交流学习：

liver_RNAseq_training.txt - 训练集，包含134个个体的年龄、性别和15938个基因的表达值，通过RNA测序进行量化。
liver_RNAseq_validation.txt - 验证/测试集，包含147个个体的性别和相同的15938个基因的表达值，通过RNA测序进行量化。
ExampleSubmission.csv - 正确格式的示例提交文件。

由于基因数量较多，我首先想到的就是通过筛选关键基因进行线性回归预测。

项目实现

数据预处理

线性回归

在之前的BLOG里我曾经讲到过线性回归模型的实现原理，这次我首先实现了C++代码实现了基于梯度下降法的第一版程序：

void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;//维数
    num = 10000000;//下降次数 

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = ans[i];
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
                //if(m == 0) printf("QwQ\n"); 
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            ans[i] = temp[i];
            //printf("%d ", ans[i]);
        } 
    }//梯度下降
    return ;
}

筛选数据的线性回归

memset(can, false, sizeof(can));
can[94] = can[53] = can[36] = can[21] = can[92] = true;
can[124] = can[39] = can[87] = true;

但是结果不是很好，只到了Kaggle上0.22252，所以我又放弃了这一方法。

带PCA的线性回归

过后我在课堂上听到了老师讲解PCA，结合之前写过的BLOG中PCA的原理，也自己去学习了C++下的PCA实现方法，并参考博客完成了PCA的相关函数代码：

void featurenormalize(MatrixXd &X) {
    //计算每一维度均值
    MatrixXd meanval = X.colwise().mean();
    RowVectorXd meanvecRow = meanval;
    //样本均值化为0
    X.rowwise() -= meanvecRow;
}
void computeCov(MatrixXd &X, MatrixXd &C) {
    //计算协方差矩阵C = XTX / n-1;
    C = X.adjoint() * X;
    C = C.array() / X.rows() - 1;
}
void computeEig(MatrixXd &C, MatrixXd &vec, MatrixXd &val) {
    //计算特征值和特征向量，使用selfadjont按照对阵矩阵的算法去计算，可以让产生的vec和val按照有序排列
    SelfAdjointEigenSolver<MatrixXd> eig(C);

    vec = eig.eigenvectors();
    val = eig.eigenvalues();
}
int computeDim(MatrixXd &val) {
    int dim;
    double sum = 0;
    for (int i = val.rows()-1; i >= 0; --i)
    {
        sum += val(i, 0);
        dim = i;

        if (sum / val.sum() >= 0.95)
            break;
    }
    return val.rows() - dim;
}

但是经过PCA压到10维后也只跑出了Kaggle上0.28324的成绩，所以我也放弃了这种做法。

带正则化的线性回归

double lambd = 1; 
void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;
    num = 8000000;

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = (1 - Alpha * lambd / (double)n) * ans[i];
        //正则化！！！！ 
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) ans[i] = temp[i];
    }//梯度下降
    return ;
}

所以最终代码就是上面所提到的带正则化的线性回归啦，下为代码：

#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
using namespace std;
int n, m;
double ans[11000];//各项系数
double temp[11000], h[11000];//辅助数组
double x[1100][1100], y[11000];//训练集特征值和答案
double xx[1100] = {0}, yy[1100] = {0};//训练集和测试集各基因的平均值
const double Alpha = 3 * 1e-6;//单次梯度下降的幅度
void read() {//读入及相关操作
    n = 134;
    m = 10;

    for(int i = 1; i <= n; i++) x[i][0] = 1;
    char s[100];
    FILE *fin = fopen("liver_RNAseq_training.txt", "r");
    for(int k = 1; k <= 15940; k++) fscanf(fin, "%s", s);
    for(int i = 1; i <= n; i++) { 
        fscanf(fin, "%s", s);
        if(s[0] == 'M') x[i][1] = 1;
        fscanf(fin, "%lf", &y[i]);
        //printf("%lf ", y[i]);
        for(int j = 2; j <= m; j++) {
            fscanf(fin, "%lf", &x[i][j]);
            xx[j] += x[i][j] / double(n);
        }
        for(int k = 1; k <= 15940 - m - 1; k++) fscanf(fin, "%s", s);
    }//数据读入
    xx[1] = 1;
    for(int i = 1; i <= n; i++) {
        for(int j = 1; j <= m; j++) {
            x[i][j] /= xx[j];
        }
    }//数据均一化
    fclose(fin); 
    return;
}
double lambd = 1; //正则化系数
void work() {
    memset(ans, 0, sizeof(ans));
    double delta = 1;
    int num, m = 10;//维数
    num = 8000000;//下降次数

    for(int o = 1; o <= num; o++) {
        for(int i = 0; i <= m; i++) temp[i] = (1 - Alpha * lambd / (double)n) * ans[i];
        //正则化！！！！ 
        //==========================================
        for(int i = 1; i <= n; i++) {
            h[i] = 0;
            for(int j = 0; j <= m; j++)
            {
                h[i] += x[i][j] * ans[j];
            }
            h[i] -= y[i];
        }
        //==========================================
        for(int i = 0; i <= m; i++) {
            for(int j = 1; j <= n; j++) {
                temp[i] = temp[i] - (Alpha / (double)m) * h[j] * x[j][i];
            }
        }
        //==========================================
        for(int i = 0; i <= m; i++) ans[i] = temp[i];//ans为系数
    }//梯度下降
    return ;
}

void write() {//预测及输出
    FILE *fin = fopen("liver_RNAseq_validation.txt", "r"), *fout = fopen("ans.txt", "w");
    int n = 147, m = 10;
    //for(int i = 0; i <= m; i++) printf("%lf ", ans[i]);
    char s[100];
    double read[200][200];
    yy[1] = 1;
    for(int k = 1; k <= 15940; k++) fscanf(fin, "%s", s);
    for(int i = 1; i <= n; i++) { 
        fscanf(fin, "%s", s);
        fscanf(fin, "%s", s);
        if(s[0] == 'M') read[i][1] = 1;
        read[i][0] = 1;
        for(int j = 2; j <= m; j++) {
            fscanf(fin, "%lf", &read[i][j]);
            yy[j] += read[i][j] / double(n);

        }
        //printf("\n");
        for(int k = 1; k <= 15940 - m - 1; k++) fscanf(fin, "%s", s);
    }//读入预测数据
    //for(int i = 1; i <= m; i++) printf("%lf ", yy[i]);
    for(int i = 1; i <= n; i++) {
        double x = 0;
        x = read[i][0] * ans[0];
        for(int j = 1; j <= m; j++) {
            read[i][j] /= yy[j];//正则化
            x += read[i][j] * ans[j]; //计算答案
        }
        fprintf(fout, "%0.0lf\n", x);//输出答案
    }
    fclose(fin); 
    fclose(fout);
    return;
}
int main() {

    read();
    work();
    write();
    return 0;
}

一些小问题

#include <fstream>
#include <cstdio>
using namespace std;
int ans[200];
int main()
{        
    FILE *fin = fopen("ans.txt", "r");
    for(int i = 1; i <= 147; i++) fscanf(fin, "%d", &ans[i]);
    fclose(fin);
    ofstream oFile;

        //打开要输出的文件
    oFile.open("ans.csv", ios::out | ios::trunc);
    oFile << "ID" << "," << "Age"  << endl;
    for(int i = 1; i <= 147; i++) {
        oFile << "test" << i << "," << ans[i]  << endl;
    }

    oFile.close();
    return 0;
}

【蒟蒻科技】肝信息预测年龄

问题描述

项目实现

数据预处理

线性回归

筛选数据的线性回归

带PCA的线性回归

带正则化的线性回归

一些小问题

个人总结

4 comments

Leave a Comment Cancel reply

【蒟蒻数据结构】主席树

【蒟蒻科技】实现光的波长和RGB值的转换

【蒟蒻计算几何】旋转卡壳算法

【蒟蒻数学】DFT和FFT

【蒟蒻计算几何】二维凸包

【蒟蒻科技】使用vbs脚本控制键盘和鼠标

【蒟蒻算法】反悔贪心

【蒟蒻机器学习基础】Octave中的基础操作

【蒟蒻机器学习基础】PCA的应用

【蒟蒻数学】巧妙运用二分图解决博弈论

【蒟蒻科技】肝信息预测年龄

问题描述

项目实现

数据预处理

线性回归

筛选数据的线性回归

带PCA的线性回归

带正则化的线性回归

一些小问题

个人总结