0 前言

??這兩年開始畢業(yè)設(shè)計(jì)和畢業(yè)答辯的要求和難度不斷提升，傳統(tǒng)的畢設(shè)題目缺少創(chuàng)新和亮點(diǎn)，往往達(dá)不到畢業(yè)答辯的要求，這兩年不斷有學(xué)弟學(xué)妹告訴學(xué)長自己做的項(xiàng)目系統(tǒng)達(dá)不到老師的要求。

為了大家能夠順利以及最少的精力通過畢設(shè)，學(xué)長分享優(yōu)質(zhì)畢業(yè)設(shè)計(jì)項(xiàng)目，今天要分享的是

?? 大數(shù)據(jù)醫(yī)學(xué)數(shù)據(jù)分析 - 人體腸道細(xì)菌數(shù)據(jù)分析

??學(xué)長這里給一個(gè)題目綜合評(píng)分(每項(xiàng)滿分5分)

• 難度系數(shù)：3分

• 工作量：3分

• 創(chuàng)新點(diǎn)：4分

畢設(shè)幫助，選題指導(dǎo)，技術(shù)解答，歡迎打擾，見B站個(gè)人主頁

https://space.bilibili.com/33886978

1 課題背景

關(guān)于人體腸道細(xì)菌的研究表明，腸道菌群與人體的年齡、健康狀況、飲食習(xí)慣等都存在關(guān)聯(lián)，為了比較精確的探尋腸道菌群與人體年齡的關(guān)系，我們基于腸道細(xì)菌2388 個(gè)微生物特征(OTUs 和genus) 以及1638 個(gè)樣本數(shù)據(jù)(重慶、泉州樣本) 進(jìn)行數(shù)據(jù)分析。

分析目標(biāo)

• 對(duì)2000 多人的腸道細(xì)菌16S 數(shù)據(jù)進(jìn)行初步探索, 將原始數(shù)據(jù)整合成具有分析意義的數(shù)據(jù)。

• 對(duì)數(shù)據(jù)進(jìn)行特征分析, 篩選出一些能有效預(yù)測(cè)人體年齡的特征。

• 通過前期處理后的數(shù)據(jù), 建立一個(gè)合理模型, 使之具有很好的預(yù)測(cè)效果。

2.數(shù)據(jù)預(yù)處理

導(dǎo)入工具包

? ?import numpy as np
? ?import pandas as pd
? ?import matplotlib.pyplot as plt
? ?%matplotlib inline

預(yù)處理

? ?path = './data'

? ?databg = pd.read_csv(path+'/2054samples backgroud information.csv')
? ?datagenus = pd.read_csv(path+'/2054samples profiling of genus.csv')
? ?dataotu = pd.read_csv(path+'/2054samples profiling of OTUs.csv')

? ?databg = databg.iloc[:,0:18]

? ?databg = databg[(databg['Sampled Loci'] == 'chongqing')|(databg['Sampled Loci'] == 'quanzhou')]

? ?datagenus = pd.merge(databg,datagenus,on = 'SampleID',how = 'inner')

? ?dataotu = pd.merge(databg,dataotu,on = 'SampleID',how = 'inner')

缺失值情況

?databg.isnull().sum() ?#q缺失值情況

? ?datagenus.isnull().sum()[datagenus.isnull().sum()!=0]

? ?dataotu.isnull().sum()[dataotu.isnull().sum()!=0]

3 數(shù)據(jù)可視化

要研究因變量年齡(或BMI) 與腸道細(xì)菌之間的聯(lián)系以及更好地進(jìn)行預(yù)測(cè)，則有必要對(duì)年齡指標(biāo)進(jìn)行可視化分析，下圖展示了年齡以及BMI 指標(biāo)的分布情況：

從上圖中可看出，年齡數(shù)據(jù)的取值在18 與66 之間，分布主要集中在35-50 歲之間，而BMI 數(shù)據(jù)則更為集中，BMI 指數(shù)主要集中在27 左右；因此從一定程度上可說明預(yù)測(cè)年齡的均方誤差比預(yù)測(cè)BMI 更大；但是從另一個(gè)角度來說，BMI 數(shù)據(jù)與年齡數(shù)據(jù)存在一定的正相關(guān)關(guān)系，因此我們最終把年齡作為主要研究對(duì)象進(jìn)行分析。

下圖展示了背景信息中兩個(gè)離散變量的分布情況：

從上圖發(fā)現(xiàn)來自重慶地區(qū)的樣本和泉州地區(qū)的樣本數(shù)量相當(dāng)。而性別變量，性別1 的數(shù)量則是性別2 數(shù)量的一倍之多。 ? 細(xì)菌數(shù)據(jù)總共包含了2388 個(gè)特征，但大多數(shù)特征對(duì)預(yù)測(cè)年齡是幾乎沒有幫助的，并且特征數(shù)量過多對(duì)后續(xù)模型的精度均有較大影響，因此我們必須進(jìn)行有效的特征篩選。

下面是實(shí)現(xiàn)數(shù)據(jù)可視化的R代碼：

? multiplot <- function(..., plotlist=NULL, file, cols=1, layout=NULL) {
? ? ?library(grid)
? ? ?

? ? ?# Make a list from the ... arguments and plotlist
? ? ?plots <- c(list(...), plotlist)
? ? ?
? ? ?numPlots = length(plots)
? ? ?
? ? ?# If layout is NULL, then use 'cols' to determine layout
? ? ?if (is.null(layout)) {
? ? ? ?# Make the panel
? ? ? ?# ncol: Number of columns of plots
? ? ? ?# nrow: Number of rows needed, calculated from # of cols
? ? ? ?layout <- matrix(seq(1, cols * ceiling(numPlots/cols)),
? ? ? ? ? ? ? ? ? ? ? ? ncol = cols, nrow = ceiling(numPlots/cols))
? ? ?}
? ? ?
? ? ?if (numPlots==1) {
? ? ? ?print(plots[[1]])
? ? ? ?
? ? ?} else {
? ? ? ?# Set up the page
? ? ? ?grid.newpage()
? ? ? ?pushViewport(viewport(layout = grid.layout(nrow(layout), ncol(layout))))
? ? ? ?
? ? ? ?# Make each plot, in the correct location
? ? ? ?for (i in 1:numPlots) {
? ? ? ? ?# Get the i,j matrix positions of the regions that contain this subplot
? ? ? ? ?matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))
? ? ? ? ?
? ? ? ? ?print(plots[[i]], vp = viewport(layout.pos.row = matchidx$row,
? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?layout.pos.col = matchidx$col))
? ? ? ?}
? ? ?}
? ?}
? ?#使用示例multiplot(p1,p2,p3,p4,cols=xxx)
? ?library(ggplot2)
? ?path = 'F:/研習(xí)部/第二次案例分析/腸道16S測(cè)序數(shù)據(jù)/'
? ?databg = read.csv(paste0(path,'2054samples backgroud information.csv'))
? ?databg = databg[,-((ncol(databg)-3):ncol(databg))]
? ?databg = databg[databg$Sampled.Loci %in% c('chongqing','quanzhou'),]
? ?databg$Gender<-as.factor(databg$Gender)
? ?databg$Sampled.Loci<-as.factor(databg$Sampled.Loci)
? ?#年齡分布圖
? ?p1<-ggplot(data = databg)+geom_histogram(aes(x = Age,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = Age),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#BMI
? ?p2<-ggplot(data = databg)+geom_histogram(aes(x = BMI,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = BMI),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#性別
? ?p3<-ggplot(data = databg)+geom_bar(aes(x = Gender),fill = 'lightblue')
? ?#腰圍
? ?p4<-ggplot(data = databg)+geom_histogram(aes(x = Waistline,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = Waistline),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#sp
? ?ggplot(data = databg)+geom_histogram(aes(x = SP,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = SP),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#dp
? ?ggplot(data = databg)+geom_histogram(aes(x = DP,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = DP),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#fbg
? ?ggplot(data = databg)+geom_histogram(aes(x = FBG,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = FBG),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#膽固醇
? ?p5<-ggplot(data = databg)+geom_histogram(aes(x = Cholesterol,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = Cholesterol),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#甘油三酯
? ?p6<-ggplot(data = databg)+geom_histogram(aes(x = Triglyceride,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = Triglyceride),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#ldl
? ?ggplot(data = databg)+geom_histogram(aes(x = LDL,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = LDL),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#hdl
? ?ggplot(data = databg)+geom_histogram(aes(x = HDL,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = HDL),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#ua
? ?ggplot(data = databg)+geom_histogram(aes(x = UA,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = UA),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#肌酐
? ?ggplot(data = databg)+geom_histogram(aes(x = Creatinine,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = Creatinine),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#egfr
? ?ggplot(data = databg)+geom_histogram(aes(x = eGFR,y = ..density..),fill = 'lightblue')+
? ? ?stat_density(aes(x = eGFR),geom = 'line',position = 'identity',color = 'red',alpha = 0.8)
? ?#位置
? ?p6<-ggplot(data = databg)+geom_bar(aes(x = Sampled.Loci),fill = 'lightblue')
? ?multiplot(p1,p2,cols = 2)
? ?multiplot(p3,p6,cols = 2)

4 特征工程

首先利用最大信息系數(shù)和距離相關(guān)系數(shù)這兩個(gè)統(tǒng)計(jì)量對(duì)特征進(jìn)行初步的篩選，再對(duì)得到的特征集合利用遞歸的XGBOOST 模型，得到最優(yōu)的特征子集。

人工預(yù)先篩選

? ?ndep = databg.shape[1]

? ?dataotu.head()

最大信息系數(shù)MIC

最大信息系數(shù)對(duì)于一般的變量之間關(guān)系具有普適效應(yīng)，不僅可以發(fā)現(xiàn)變量之間的線性關(guān)系，也可以發(fā)現(xiàn)變量之間的非線性關(guān)系。 ? 最大信息系數(shù)的計(jì)算方式完全基于互信息的計(jì)算方式。對(duì)于兩個(gè)連續(xù)型隨機(jī)變量，首先將其所在的二維空間使用m 乘以n 的網(wǎng)格劃分，則可以將落在第 ( x , y ) (x, y) (x,y) 格子中的數(shù)據(jù)點(diǎn)的頻率作為 p ( x , y ) p(x, y) p(x,y) 的估計(jì)：

? ?\begin{equation}
? ?p(x,y) = \frac{n_{x,y}}{N}
? ?\end{equation}

其中 n x , y n_{x,y} nx,y為格子 ( x , y ) (x, y) (x,y)中的數(shù)據(jù)個(gè)數(shù)， N N N為總樣本數(shù)。

則根據(jù)互信息量的計(jì)算公式，可以得到隨機(jī)變量 X X X， Y Y Y 的互信息 I I I 為：

將Age作為因變量，對(duì)genus數(shù)據(jù)集和otu數(shù)據(jù)集中的所有變量計(jì)算與年齡的MIC，并繪制成散點(diǎn)圖：

index1 = list(range(ndep,datagenus.shape[1]))
? ?index2 = list(range(ndep,dataotu.shape[1]))


? ?from minepy import MINE
? ?def calmic(a,b):
? ? ? ?m = MINE()
? ? ? ?m.compute_score(a,b)
? ? ? ?return m.mic()
? ?MIC1 = []
? ?MIC2 = []
? ?for i in range(ndep,datagenus.shape[1]):
? ? ? ?MIC1.append(calmic(datagenus['Age'],datagenus.iloc[:,i]))
? ?for i in range(ndep,dataotu.shape[1]):
? ? ? ?MIC2.append(calmic(dataotu['Age'],dataotu.iloc[:,i]))


? ?fig1 = plt.figure()
? ?fig2 = plt.figure()
? ?ax1 = fig1.add_subplot(1,1,1)
? ?ax1.scatter(range(len(MIC1)),MIC1)
? ?fig1.savefig('F:/研習(xí)部/第二次案例分析/ElegantNote 2.00/image/c1.png')
? ?ax2 = fig2.add_subplot(1,1,1)
? ?ax2.scatter(range(len(MIC2)),MIC2)
? ?fig2.savefig('F:/研習(xí)部/第二次案例分析/ElegantNote 2.00/image/c2.png')

從上面兩張圖（第一張為數(shù)據(jù)集genus的 M I C MIC MIC，第二張為數(shù)據(jù)集otu的 M I C MIC MIC）可以看到，大部分特征的MIC值都十分集中，接近于0(即與Age沒有明顯的聯(lián)系)。若將MIC值等于0.04作為閾值，可以發(fā)現(xiàn)大于該值的特征數(shù)量較少，MIC值均較大，與MIC值小于0.04的特征有明顯的分化趨勢(shì)。因此將 M I C = 0.04 MIC=0.04 MIC=0.04作為篩選特征的一個(gè)衡量標(biāo)準(zhǔn)。

距離相關(guān)系數(shù)dcorr

與最大信息系數(shù)類似，距離相關(guān)系數(shù)也是用于衡量特征之間相關(guān)性的一個(gè)統(tǒng)計(jì)量，且其具有較好的魯棒性。

距離相關(guān)系數(shù)也是介于0與1之間的統(tǒng)計(jì)量。距離相關(guān)系數(shù)為0，說明兩個(gè)特征之間相互獨(dú)立，相反其值越大，說明兩個(gè)特征之間的相關(guān)程度越大。

依舊是將Age作為因變量，對(duì)genus數(shù)據(jù)集和otu數(shù)據(jù)集中的所有變量計(jì)算與年齡之間的最大相關(guān)系數(shù)，并將最大相關(guān)系數(shù)繪制成散點(diǎn)圖如下圖所示。

?from scipy.spatial.distance import pdist, squareform
? ?import numpy as np
? ?def distcorr(X, Y):
? ? ? ?X = np.atleast_1d(X)
? ? ? ?Y = np.atleast_1d(Y)
? ? ? ?if np.prod(X.shape) == len(X):
? ? ? ? ? ?X = X[:, None]
? ? ? ?if np.prod(Y.shape) == len(Y):
? ? ? ? ? ?Y = Y[:, None]
? ? ? ?X = np.atleast_2d(X)
? ? ? ?Y = np.atleast_2d(Y)
? ? ? ?n = X.shape[0]
? ? ? ?if Y.shape[0] != X.shape[0]:
? ? ? ? ? ?raise ValueError('Number of samples must match')
? ? ? ?a = squareform(pdist(X))
? ? ? ?b = squareform(pdist(Y))
? ? ? ?A = a - a.mean(axis=0)[None, :] - a.mean(axis=1)[:, None] + a.mean()
? ? ? ?B = b - b.mean(axis=0)[None, :] - b.mean(axis=1)[:, None] + b.mean()
? ? ? ?

? ? ? ?dcov2_xy = (A * B).sum()/float(n * n)
? ? ? ?dcov2_xx = (A * A).sum()/float(n * n)
? ? ? ?dcov2_yy = (B * B).sum()/float(n * n)
? ? ? ?dcor = np.sqrt(dcov2_xy)/np.sqrt(np.sqrt(dcov2_xx) * np.sqrt(dcov2_yy))
? ? ? ?return dcor

? ?dcorr1 = []
? ?dcorr2 = []
? ?for i in range(ndep,datagenus.shape[1]):
? ? ? ?dcorr1.append(distcorr(datagenus['Age'],datagenus.iloc[:,i]))
? ?for i in range(ndep,dataotu.shape[1]):
? ? ? ?dcorr2.append(distcorr(dataotu['Age'],dataotu.iloc[:,i]))

? ?fig1 = plt.figure()
? ?fig2 = plt.figure()
? ?ax1 = fig1.add_subplot(1,1,1)
? ?ax1.scatter(range(len(dcorr1)),dcorr1)
? ?ax2 = fig2.add_subplot(1,1,1)
? ?ax2.scatter(range(len(dcorr2)),dcorr2)
? ?fig1.savefig('F:/研習(xí)部/第二次案例分析/ElegantNote 2.00/image/c3.png')
? ?fig2.savefig('F:/研習(xí)部/第二次案例分析/ElegantNote 2.00/image/c4.png')

從上圖（第一張為數(shù)據(jù)集genus的 d c o r dcor dcor，第二張為數(shù)據(jù)集otu的 d c o r dcor dcor）可以看到，大部分特征的最大相關(guān)系數(shù)值都十分集中，接近于0(即與Age沒有明顯的聯(lián)系)。若將 d c o r dcor dcor值等于0.06作為閾值，可以發(fā)現(xiàn)大于該值的特征數(shù)量較少， d c o r dcor dcor值均較大，與 d c o r dcor dcor值小于0.06的特征有明顯的分化趨勢(shì)。因此將 d c o r = 0.06 dcor=0.06 dcor=0.06也作為篩選特征的一個(gè)衡量標(biāo)準(zhǔn)。

篩選 MIC是0.04，dcorr是0.06 邏輯關(guān)系是且(或 也可以嘗試)

? ?select_index1 = [index1[i] for i in range(len(index1)) if (MIC1[i]>0.04 and dcorr1[i]>0.06)] #or
? select_index2 = [index2[i] for i in range(len(index2)) if (MIC2[i]>0.04 and dcorr2[i]>0.06)] #or

? ?print(len(select_index1),len(select_index2))

遞歸xgboost

XGBOOST

前面兩個(gè)小節(jié)將最大信息系數(shù)值大于0.04以及最大相關(guān)系數(shù)值大于0.06的特征初步篩選出來，作為與因變量年齡有一定相關(guān)關(guān)系的特征來進(jìn)一步研究。在數(shù)據(jù)集genus中總共選出了18個(gè)特征，在數(shù)據(jù)集otu中總共選出了70個(gè)特征。

為了進(jìn)一步研究腸道細(xì)菌特征和年齡之間的關(guān)系，找到具有更好預(yù)測(cè)效果的變量，本文將XGBOOST模型 [ 2 ] ^{[2]} [2]作為基模型，采用遞歸的方式對(duì)特征進(jìn)行進(jìn)一步篩選。

XGBOOST是一種集成學(xué)習(xí)算法，每一個(gè)子樹是一棵CART樹。每一棵子樹利用基尼純度來衡量每一個(gè)特征的重要性，并且以此來選擇分裂特征，而某個(gè)特征在所有子樹中作為分裂結(jié)點(diǎn)的出現(xiàn)次數(shù)就是其對(duì)應(yīng)的特征重要性。

本文將年齡作為因變量，將之前篩選出的88個(gè)特征作為輸入，得到的特征重要性如下圖所示(為了清楚的展示得分概況，下圖只隨機(jī)選取了20個(gè)特征來作為樣例)：

? ?##考慮AGE是響應(yīng)變量
? ?datagenus.head()
? ?dataxgb = pd.merge(datagenus.iloc[:,[0,2]+select_index1],dataotu.iloc[:,[0]+select_index2],on = 'SampleID',how = 'inner')

? ?dataxgb.head()

? ?import xgboost as xgb
? ?from sklearn.feature_selection import RFE
? ?from sklearn.model_selection import train_test_split
? ?from sklearn.model_selection import cross_val_score
? ?Y = dataxgb.iloc[:,1]
? ?X = dataxgb.iloc[:,2:]

? ?model = xgb.XGBRegressor()
? ?model.fit(X.iloc[:,list(set(np.random.randint(0,X.shape[1]+1,20)))],Y)
? ?from xgboost import plot_importance ?#畫出重要性圖片
? ?plot_importance(model)
? ?fig=plt.gcf()
? ?fig.savefig('F:/研習(xí)部/第二次案例分析/ElegantNote 2.00/image/c5.png')

遞歸篩選

在選擇XGBOOST為基模型后，利用遞歸篩選的方式來選擇最優(yōu)的特征子集。即初始子集是原始特征集合(包含88個(gè)特征)，每次刪去當(dāng)前集合中一個(gè)最不重要的特征，利用3折交叉驗(yàn)證計(jì)算模型的預(yù)測(cè)得分。最終選取得分最高的子集作為最終的特征集合，具體算法流程算法一所示：

算法一：遞歸篩選特征 ?

step1：初始特征集合為F ? step2: 選出其中最不重要的特征 F i F{i} Fi，并且記錄此時(shí)的模型得分 s s s ? step3: 將特征集合更新為 F ? F i F-F{i} F?Fi，跳轉(zhuǎn)回step2直至 F = ? F=\varnothing F=? ? step4: 根據(jù)記錄的模型得分選出最優(yōu)得分對(duì)應(yīng)的集合作為最終特征集合 S S S ?

將88個(gè)特征輸入算法，最終篩選得到39個(gè)特征，genus數(shù)據(jù)集包含7個(gè)特征，otu數(shù)據(jù)集包含32個(gè)特征。這些特征如下表所示(列舉15個(gè))：

重要特征（列舉15個(gè)） ?

Alistipes Anaerostipes Blautia Otu67 Otu68 ? Erysipelatoclostridium Granulicatella Otu10 Otu8 Otu80 ? Salmonella Turicibacter Otu119 Otu102 Otu99 ?

由于genus數(shù)據(jù)集與otu數(shù)據(jù)集是按照兩種方式來分類細(xì)菌的，因此兩個(gè)數(shù)據(jù)集中的特征會(huì)有所重合，而重合的特征會(huì)使得數(shù)據(jù)噪聲增加，使模型產(chǎn)生過擬合現(xiàn)象。因此本文對(duì)篩選得到的39個(gè)特征利用特征注釋表進(jìn)行人工查重，發(fā)現(xiàn)此39個(gè)特征未包含重復(fù)特征。

def ?rec_select(a,b): #定義遞歸選取變量的函數(shù)，每次刪去最不重要的變量，直至刪完，選取精度最好的子集
? ? ? ?global ans
? ? ? ?ans = []
? ? ? ?for i in range(a,b+1): #a
? ? ? ? ? ?model = xgb.XGBRegressor()
? ? ? ? ? ?selector = RFE(estimator=model, n_features_to_select=i) #特征選取數(shù)量從小到大
? ? ? ? ? ?X_choice = selector.fit_transform(X, Y) #選出新的數(shù)據(jù)集
? ? ? ? ? ?X_train, X_test, Y_train, Y_test = train_test_split(X_choice,Y,test_size=0.3, random_state=0)
? ? ? ? ? ?ans.append(np.mean(cross_val_score(model,X_choice,Y,cv = 3))) #交叉驗(yàn)證作為得分
? ? ? ? ? ?if i%10==0:
? ? ? ? ? ? ? ?print(i)
? ? ? ?index = ans.index(max(ans))
? ? ? ?return list(range(a,b+1))[index]
? ?rec_select(20,88)

? ?selector = RFE(estimator=model, n_features_to_select=39)
? ?selector.fit(X, Y)
? ?selector.support_

? ?X.iloc[:,selector.support_].head()

5 神經(jīng)網(wǎng)絡(luò)模型

網(wǎng)絡(luò)的架構(gòu)

本文選擇只包含全連接層的神經(jīng)網(wǎng)絡(luò)并且結(jié)合dropout層來進(jìn)行年齡預(yù)測(cè)。網(wǎng)絡(luò)示例如下圖所示：

為了使樣本各個(gè)特征的比重相同，且使得神經(jīng)網(wǎng)絡(luò)的反向傳播算法得以很好的收斂，本文首先對(duì)數(shù)據(jù)進(jìn)行如下處理： ?

其中 X X X表示特征數(shù)據(jù)集，隨后再對(duì)數(shù)據(jù)進(jìn)行歸一化處理。

超參數(shù)的確定

在確定完網(wǎng)絡(luò)架構(gòu)后，模型依舊有幾個(gè)重要的超參數(shù)需要確定。

學(xué)習(xí)率的確定，學(xué)習(xí)率需要根據(jù)模型的輸出結(jié)果手動(dòng)調(diào)節(jié)，若模型損失函數(shù)變化很小，則需要適當(dāng)調(diào)大學(xué)習(xí)率，若模型損失函數(shù)反復(fù)變化，則需要降低學(xué)習(xí)率，并采用自適應(yīng)學(xué)習(xí)率，本文設(shè)置為0.001。

學(xué)習(xí)批次的確定，由于樣本量較大，若完整訓(xùn)練的數(shù)據(jù)集再更新權(quán)重會(huì)耗費(fèi)大量運(yùn)算時(shí)間，選擇訓(xùn)練批次是合理的解決方法，一般批次大小取2的指數(shù)次方倍，本文的批次設(shè)置為32,64,128或256。從中選取損失最小的值作為最終批次，經(jīng)過參數(shù)調(diào)試選擇256為訓(xùn)練批次。

迭代次數(shù)的確定，一般選取若繼續(xù)增大迭代次數(shù)時(shí)，損失函數(shù)繼續(xù)下降而測(cè)試集上誤差卻上升的最小的迭代次數(shù)，本文選取為100。

模型層數(shù)的確定，層數(shù)設(shè)定為3層，對(duì)于卷積神經(jīng)網(wǎng)絡(luò)更多的隱藏層意味著更強(qiáng)的解釋能力，但對(duì)于全連接層網(wǎng)絡(luò)，更深的層對(duì)網(wǎng)絡(luò)解釋能力幫助并不大，但是3層網(wǎng)絡(luò)普遍優(yōu)于2層網(wǎng)絡(luò)。

隱藏層結(jié)點(diǎn)數(shù)確定，一般將隱藏層結(jié)點(diǎn)數(shù)設(shè)置為一個(gè)大于輸入值的數(shù)，本文將三個(gè)隱藏層結(jié)點(diǎn)數(shù)均設(shè)置為64。

激活函數(shù)的確定，本文的輸出層利用Relu函數(shù)來激活，三個(gè)隱藏層的激活函數(shù)均使用PRelu函數(shù)來激活。

Dropout層丟棄率的確定，為了防止網(wǎng)絡(luò)過擬合，為每一個(gè)隱藏層再加上Dropout層，即每一次訓(xùn)練時(shí)隨機(jī)丟棄該層中一定比例的神經(jīng)元，丟棄率本文設(shè)置為0.2

網(wǎng)絡(luò)權(quán)重的初始值均采0與1之間的均勻分布。

6 結(jié)果與討論

模型結(jié)果

將特征工程處理后的數(shù)據(jù)帶入神經(jīng)網(wǎng)絡(luò)模型訓(xùn)練，并進(jìn)行10折交叉驗(yàn)證，將平均絕對(duì)誤差(MAE)作為損失衡量標(biāo)準(zhǔn)。最終得到的誤差為6.81歲。為了更客觀的考量預(yù)測(cè)結(jié)果的合理性，本文分別計(jì)算了訓(xùn)練數(shù)據(jù)集年齡與BMI的相關(guān)系數(shù)以及預(yù)測(cè)年齡與對(duì)應(yīng)BMI的相關(guān)系數(shù)，分別為0.224與0.152?？梢钥吹絻烧叨汲尸F(xiàn)一定的正相關(guān)關(guān)系，預(yù)測(cè)結(jié)果很好的保持了數(shù)據(jù)原有的性質(zhì)。

取交叉驗(yàn)證中的某一測(cè)試數(shù)據(jù)集來進(jìn)行可視化，將真實(shí)年齡作為橫坐標(biāo)，預(yù)測(cè)年齡作為縱坐標(biāo)，并將真實(shí)年齡作為基準(zhǔn)線，再考慮預(yù)測(cè)樣本之間的性別差異，畫出預(yù)測(cè)散點(diǎn)圖:

從上圖中可以看到大多數(shù)樣本點(diǎn)都落在平均誤差的范圍內(nèi)，有一部分點(diǎn)的預(yù)測(cè)誤差較大，誤差大于平均意義上的6.29歲。且這些樣本點(diǎn)主要集中在30歲以下和55歲以上的樣本，這與訓(xùn)練數(shù)據(jù)集兩端年齡的樣本數(shù)量不足有關(guān)。同時(shí)也可以發(fā)現(xiàn)不同的性別其年齡預(yù)測(cè)的結(jié)果大致相同，性別這一變量并不會(huì)對(duì)年齡的預(yù)測(cè)產(chǎn)生重大影響。

若考慮預(yù)測(cè)樣本地區(qū)之間的差異，預(yù)測(cè)散點(diǎn)圖如下圖所示：

從上圖可以發(fā)現(xiàn)，對(duì)于50歲以上的樣本，對(duì)于重慶地區(qū)的年齡預(yù)測(cè)誤差要比泉州地區(qū)的預(yù)測(cè)誤差更大。若要繼續(xù)提高預(yù)測(cè)精度，可以考慮對(duì)50歲以上的樣本進(jìn)行分地區(qū)預(yù)測(cè)。

由于BMI變量也具有重要的實(shí)際意義，因此本文也將BMI變量作為因變量，采用同樣的模型進(jìn)行特征篩選和數(shù)據(jù)訓(xùn)練，最終的預(yù)測(cè)誤差(MAE)為1.94，結(jié)果也較為理想，其預(yù)測(cè)結(jié)果如下圖所示(BMI的高預(yù)測(cè)精度主要源于其數(shù)據(jù)方差較小，年齡數(shù)據(jù)分布的標(biāo)準(zhǔn)差為10.1，而BMI數(shù)據(jù)分布的標(biāo)準(zhǔn)差3.01，僅為年齡標(biāo)準(zhǔn)差的三分之一！)

連續(xù)SMOTE上采樣

將利用特征工程小節(jié)中的方法篩選得到的數(shù)據(jù)劃分為訓(xùn)練集與測(cè)試集，并且畫出訓(xùn)練數(shù)據(jù)的年齡數(shù)據(jù)分布，如下圖右圖所示：

通過上圖右圖可以發(fā)現(xiàn)，樣本數(shù)據(jù)主要集中在30-55這個(gè)年齡段，可以預(yù)見到如此的訓(xùn)練樣本分布會(huì)導(dǎo)致最終的預(yù)測(cè)結(jié)果也會(huì)集中分布在這個(gè)年齡段，從而無法很好預(yù)測(cè)30歲以下，55歲以上的樣本。

為了解決這一問題，對(duì)訓(xùn)練集合中29歲以下，58歲以上的樣本進(jìn)行了SMOTE過采樣^{[1]}的方式(采樣時(shí)將年齡也歸入自變量一同進(jìn)行采樣)，近鄰樣本數(shù)選取為5個(gè)。最后得到的數(shù)據(jù)分布如上圖左圖所示。

經(jīng)過上采樣，各年齡段的樣本分布明顯均勻了很多，這也有利于模型做出更為合理的判斷。

模型比較

將利用上采樣數(shù)據(jù)集訓(xùn)練的模型(所有人工生成數(shù)據(jù)均未落入測(cè)試集中)與利用原始數(shù)據(jù)訓(xùn)練的模型進(jìn)行比較，其結(jié)果如表：

數(shù)據(jù)集MAE預(yù)測(cè)值域SMOTE上采樣數(shù)據(jù)集6.81[ 28.9 , 60.1 ] [28.9, 60.1] [28.9,60.1]原始數(shù)據(jù)集6.8130.1 , 56.6 {30.1, 56.6} 30.1,56.6

可以看到盡管SMOTE上采樣未對(duì)模型的預(yù)測(cè)精度有所改進(jìn)，卻使得模型的預(yù)測(cè)值域更大。而在精度相同的情況下，我們更希望獲得一個(gè)有更多響應(yīng)情況的模型，這也證明了進(jìn)行SMOTE上采樣的合理性。

本文的訓(xùn)練特征既包含了genus數(shù)據(jù)集中的特征，又包含了OTU數(shù)據(jù)集中的特征，將訓(xùn)練結(jié)果與使用單獨(dú)數(shù)據(jù)集的訓(xùn)練結(jié)果進(jìn)行對(duì)比，如下表所示：

數(shù)據(jù)集MAEGenus數(shù)據(jù)集8.03Otu數(shù)據(jù)集7.14合并數(shù)據(jù)集6.81

從上表中可以看到，模型在合并訓(xùn)練集上得到的預(yù)測(cè)精度最高。因此與年齡相關(guān)的細(xì)菌類別并不一定是基于一種分類方式得到的類別的，可能與多種分類方式的多個(gè)細(xì)菌類別息息相關(guān)。

將模型與其他模型進(jìn)行橫向?qū)Ρ?，結(jié)果如下表所示：

模型MAE(原始數(shù)據(jù)集)MAE(采樣數(shù)據(jù)集)神經(jīng)網(wǎng)絡(luò)6.816.81隨機(jī)森林8.258.61XGBOOST7.838.83線性回歸7.809.22

可以看到與一些常用模型比較，本文的預(yù)測(cè)精度擁有較為明顯的優(yōu)勢(shì)。

畢設(shè)幫助，選題指導(dǎo)，技術(shù)解答，歡迎打擾，見B站個(gè)人主頁

https://space.bilibili.com/33886978