數(shù)學(xué)建模：齒輪箱故障診斷

表單

gearbox00

為齒輪箱正常工況下采集到的振動(dòng)信號(hào);表單

gearbox10

為故障狀態(tài)

1

下采集到的振動(dòng)信號(hào);表單

gearbox20

為故障狀態(tài)

2

下采集到的故障信號(hào);表單

gearbox30

為故障狀態(tài)

3

下采集到的故障信號(hào);表單

gearbox40

為故障狀態(tài)

4

下采集到的振動(dòng)信號(hào)。1、對(duì)齒輪箱各個(gè)狀態(tài)下的振動(dòng)數(shù)據(jù)進(jìn)行分析，研究正常和不同故障狀態(tài)下振動(dòng)數(shù)據(jù)的變化規(guī)律及差異，并給出刻畫(huà)這些差異的關(guān)鍵特征。這題每個(gè)狀態(tài)有四個(gè)指標(biāo)，即對(duì)每個(gè)狀態(tài)的各指標(biāo)首先做圖觀察，看看每個(gè)狀態(tài)下的數(shù)據(jù)變化趨勢(shì)，明天我簡(jiǎn)單做做給出部分圖和代碼，不過(guò)用spss或者Excel也是一樣可以的。至于差異的關(guān)鍵特征，對(duì)數(shù)據(jù)進(jìn)行特征分析即可。代碼也簡(jiǎn)單，主要就是看寫(xiě)論文吧。2、建立齒輪箱的故障檢測(cè)模型，對(duì)其是否處于故障狀態(tài)進(jìn)行檢測(cè)，并對(duì)模型的性能進(jìn)行評(píng)價(jià)。

根據(jù)第一問(wèn)提取的數(shù)據(jù)特征通過(guò)計(jì)算故障強(qiáng)度系數(shù)再進(jìn)行后續(xù)計(jì)算，見(jiàn)下文《普通公路隧道機(jī)電設(shè)施故障檢測(cè)方法研究》用小波分析法檢測(cè)故障，見(jiàn)下文《基于向量自回歸模型和小波分析法的列車(chē)充電機(jī)電流傳感器故障檢測(cè)方法》3、建立齒輪箱的故障診斷模型，對(duì)其處于何種故障狀態(tài)進(jìn)行判斷，并對(duì)模型的性能進(jìn)行評(píng)價(jià)。根據(jù)前一問(wèn)的結(jié)果先檢測(cè)哪些是故障，后用一種基于樹(shù)形卷積神經(jīng)網(wǎng)絡(luò)模型診斷是哪種故障。見(jiàn)下文《鐵路信號(hào)聯(lián)鎖故障診斷模型構(gòu)建及仿真》4、結(jié)合所建立的故障檢測(cè)和診斷模型對(duì)附件

2

中另行采集的

12

組測(cè)試數(shù)據(jù)進(jìn)行檢測(cè)和診斷分析，將分析結(jié)果填寫(xiě)到下表中(注:測(cè)試數(shù)據(jù)中可能存在除以上

4

種故障之外的故障狀態(tài)，若存在，則將對(duì)應(yīng)的診斷結(jié)果標(biāo)記為:其它故障)，并將此表格放到論文的正文中第四問(wèn)直接運(yùn)用前兩問(wèn)的模型即可。明天繼續(xù)更新，今天時(shí)間有點(diǎn)倉(cāng)促，只能大概提供這一點(diǎn)點(diǎn)。明天我會(huì)進(jìn)行修改和補(bǔ)充，大家有問(wèn)題評(píng)論或私信即可，明天我的課還是比較多的。所以應(yīng)該周末才會(huì)給出詳細(xì)的思路和部分代碼。一樣的，對(duì)于數(shù)據(jù)題我們首先要做的是數(shù)據(jù)預(yù)處理，接下來(lái)再進(jìn)行下面的操作。關(guān)于問(wèn)題一：

可以做每組數(shù)據(jù)的散點(diǎn)圖，我這個(gè)就給出一張例圖。Excel和SPSS都可以做的，SPSS做的估計(jì)會(huì)更好看一點(diǎn)。Excel好好調(diào)參也可以很好看。做這種圖的，第一自己觀察看看有沒(méi)有什么顯著差異，第二寫(xiě)論文。這個(gè)也算是震動(dòng)信號(hào)數(shù)據(jù)的一種分析，論文寫(xiě)寫(xiě)就行。后面的計(jì)算有比較多的方法1、通過(guò)計(jì)算方差，篩選特征。計(jì)算方法參考下文\o"??????數(shù)據(jù)篩選特征方法-方差法_gao_vip的博客-CSDN博客_方差選擇法特征篩選"??????數(shù)據(jù)篩選特征方法-方差法_gao_vip的博客-CSDN博客_方差選擇法特征篩選2、計(jì)算特征重要性，取重要性高的那類(lèi)作為特征。3、SVD（奇異值分解），這個(gè)方法比較好的。大家可以去看看下面這個(gè)論文，這一問(wèn)能搞定?！禨VD曲率譜降噪和快速譜峭度的滾動(dòng)軸承微弱故障特征提取》\o"常用數(shù)據(jù)特征提取，時(shí)域特征、頻域特征、小波特征提取匯總；特征提??；有效matlab代碼_入間同學(xué)的博客-CSDN博客_小波熵特征提取"常用數(shù)據(jù)特征提取，時(shí)域特征、頻域特征、小波特征提取匯總；特征提??；有效matlab代碼_入間同學(xué)的博客-CSDN博客_小波熵特征提取4、小波分析特征提取%裝入變換放大器輸入輸出數(shù)據(jù)%bf_150ms.dat為正常系統(tǒng)輸出信號(hào)%bf_160ms.dat為故障系統(tǒng)輸出信號(hào)loadbf_150ms.dat;loadbf_160ms.dat;s1＝bf_150ms(1:1000);%s1為正常信號(hào)s2＝bf_160ms(1:1000);%s2為故障信號(hào)%畫(huà)出正常信號(hào)與故障信號(hào)的原始波形tittle(“原始信號(hào)’);Ylabel('s1');subplot(922);plot(s2);title('故障信號(hào)');Ylabel('s2');%============================================%用dbl小波包對(duì)正常信號(hào)s1進(jìn)行三層分解[t,d]＝wpdec(sl,3,'db','shannon');%plontree(t)%畫(huà)小波包樹(shù)結(jié)構(gòu)的圖形%下面對(duì)正常信號(hào)第三層各系數(shù)進(jìn)行重構(gòu)%s130是指信號(hào)sl的[3,0]結(jié)點(diǎn)的重構(gòu)系數(shù);其他依次類(lèi)推sl30＝wprcoef(t,d,[3,0]);s13l＝wprcoef(t,d,[3,1]);s132＝wprcoef(t,d,[3,2]);sl33＝wprcoef(t,d,[3,3]);sl34＝wprcoef(t,d,[3,4]);s135＝wprcoef(t,d,[3,5]);s136＝wprcoef(t,d,[3,6]);s137＝wprcoef(t,d,[3,7]);%畫(huà)出至構(gòu)系數(shù)的波形subplot(9,2,3);plot(s130);Ylabel('S130');subpolt(9,2,5);plot(s131);Ylabel('S13l');subplot(9,2,7);plot(s132);Ylabel('S132');subplot(9,2,9);plot(s133);Ylabel('S133');subplot(9,2,11);plot(s134);Ylabel('S134');subplot(9,2,13);plot(s135);Ylabel('S135');subplot(9,2,15);plot(s136);Ylabel('S136');subplot(9,2,17);plot(s137);Ylabel('S137');%--------------------------------------%計(jì)算正常信號(hào)各重構(gòu)系數(shù)的方差%s10是指s130的方差，其他依此類(lèi)推s10＝norm(sl30);sll＝norm(s131);s12＝norm(sl32);s13＝norm(sl33);sl4＝norm(s134);s15＝norm(s135);s16＝norm(sl36);s17＝norm(sl37);%向量ssl是針對(duì)信號(hào)s1構(gòu)造的向量disp=('正常信號(hào)的輸出向量')ssl＝[sl0,s11,sl2,sl3,s14,s15,sl6,s17]%＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝%用db1小波包對(duì)故障信號(hào)s2進(jìn)行三層分解[t,d]=wpdec(s2,3,'db1','shannon');%plottree(t)%畫(huà)小波包樹(shù)結(jié)構(gòu)的圖形%s230是指信號(hào)S2的[3,0]結(jié)點(diǎn)的重構(gòu)系數(shù),其他以此類(lèi)推s230=wprcoef(t,d,[3,0]);s231=wprcoef(t,d,[3,1]);s232=wprcoef(t,d,[3,2]);s233=wprcoef(t,d,[3,3]);s234=wprcoef(t,d,[3,4]);s235=wprcoef(t,d,[3,5]);s236=wprcoef(t,d,[3,6]);s237=wprcoef(t,d,[3,7]);%畫(huà)出重構(gòu)系數(shù)的波形subplot(9,2,4);plot(s230);Ylabel('S230');subplot(9,2,6);plot(s231);Ylabel('S231');subplot(9,2,8);plot(s232);Ylabel('S232');subplot(9,2,10);plot(s233);Ylabel('S233');subplot(9,2,12);plot(s234);Ylabel('S234');subplot(9,2,14);plot(s235);Ylabel('S235');subplot(9,2,16);plot(s236);Ylabel('S236');subplot(9,2,18);plot(s237);Ylabel('S237');%----------------------------------------------------------%計(jì)算故障信號(hào)各重構(gòu)系數(shù)的方差%s20是指s230的方差，其他依次類(lèi)推s20＝norm(s230);s21＝norm(s231);s22＝norm(s232);s23＝norm(s233);s24＝norm(s234);s25＝norm(s235);s26＝norm(s236);s27＝norm(s237);%向量ss2是針對(duì)信號(hào)S1構(gòu)造的向量disp('故障信號(hào)')來(lái)源于：\o"軸承故障診斷小波分析.7z_小波軸承故障-專(zhuān)業(yè)指導(dǎo)其他資源-CSDN下載"軸承故障診斷小波分析.7z_小波軸承故障-專(zhuān)業(yè)指導(dǎo)其他資源-CSDN下載這問(wèn)的方法實(shí)在太多了，代碼也有很全的，重點(diǎn)就是寫(xiě)論文上。關(guān)于問(wèn)題二：1、PCA故障檢測(cè)clc;clear;%%1.導(dǎo)入數(shù)據(jù)%產(chǎn)生訓(xùn)練數(shù)據(jù)num_sample=100;a=10*randn(num_sample,1);x1=a+randn(num_sample,1);x2=1*sin(a)+randn(num_sample,1);x3=5*cos(5*a)+randn(num_sample,1);x4=0.8*x2+0.1*x3+randn(num_sample,1);xx_train=[x1,x2,x3,x4];%產(chǎn)生測(cè)試數(shù)據(jù)a=10*randn(num_sample,1);x1=a+randn(num_sample,1);x2=1*sin(a)+randn(num_sample,1);x3=5*cos(5*a)+randn(num_sample,1);x4=0.8*x2+0.1*x3+randn(num_sample,1);xx_test=[x1,x2,x3,x4];xx_test(51:100,2)=xx_test(51:100,2)+15*ones(50,1);%%2.數(shù)據(jù)處理Xtrain=xx_train;Xtest=xx_test;X_mean=mean(Xtrain);X_std=std(Xtrain);[X_row,X_col]=size(Xtrain);Xtrain=(Xtrain-repmat(X_mean,X_row,1))./repmat(X_std,X_row,1);%標(biāo)準(zhǔn)化處理%%3.PCA降維SXtrain=cov(Xtrain);%求協(xié)方差矩陣[T,lm]=eig(SXtrain);%求特征值及特征向量,特征值排列順序?yàn)閺男〉酱驞=flipud(diag(lm));%將特征值從大到小排列%確定降維后的數(shù)量num=1;whilesum(D(1:num))/sum(D)<0.85num=num+1;endP=T(:,X_col-num+1:X_col);%取對(duì)應(yīng)的向量P_=fliplr(P);%特征向量由大到小排列%%4.計(jì)算T2和Q的限值%求置信度為99%時(shí)的T2統(tǒng)計(jì)控制限,T=k*(n^2-1)/n(n-k)*F(k,n-k)%其中k對(duì)應(yīng)num,n對(duì)應(yīng)X_rowT2UCL1=num*(X_row-1)*(X_row+1)*finv(0.99,num,X_row-num)/(X_row*(X_row-num));%求置信度為99%時(shí)的T2統(tǒng)計(jì)控制限%求置信度為99%的Q統(tǒng)計(jì)控制限fori=1:3th(i)=sum((D(num+1:X_col)).^i);endh0=1-2*th(1)*th(3)/(3*th(2)^2);ca=norminv(0.99,0,1);QU=th(1)*(h0*ca*sqrt(2*th(2))/th(1)+1+th(2)*h0*(h0-1)/th(1)^2)^(1/h0);%置信度為99%的Q統(tǒng)計(jì)控制限%%5.模型測(cè)試n=size(Xtest,1);Xtest=(Xtest-repmat(X_mean,n,1))./repmat(X_std,n,1);%標(biāo)準(zhǔn)化處理%求T2統(tǒng)計(jì)量，Q統(tǒng)計(jì)量[r,y]=size(P*P');I=eye(r,y);T2=zeros(n,1);Q=zeros(n,1);lm_=fliplr(flipud(lm));%T2的計(jì)算公式Xtest.T*P_*inv(S)*P_*Xtestfori=1:nT2(i)=Xtest(i,:)*P_*inv(lm_(1:num,1:num))*P_'*Xtest(i,:)';Q(i)=Xtest(i,:)*(I-P*P')*Xtest(i,:)';end%%6.繪制T2和SPE圖figure('Name','PCA');subplot(2,1,1);plot(1:i,T2(1:i),'k');holdon;plot(i:n,T2(i:n),'k');title('統(tǒng)計(jì)量變化圖');xlabel('采樣數(shù)');ylabel('T2');holdon;line([0,n],[T2UCL1,T2UCL1],'LineStyle','--','Color','r');subplot(2,1,2);plot(1:i,Q(1:i),'k');holdon;plot(i:n,Q(i:n),'k');title('統(tǒng)計(jì)量變化圖');xlabel('采樣數(shù)');ylabel('SPE');holdon;line([0,n],[QU,QU],'LineStyle','--','Color','r');%%7.繪制貢獻(xiàn)圖%7.1.確定造成失控狀態(tài)的得分S=Xtest(51,:)*P(:,1:num);r=[];fori=1:numifS(i)^2/lm_(i)>T2UCL1/numr=cat(2,r,i);endend%7.2.計(jì)算每個(gè)變量相對(duì)于上述失控得分的貢獻(xiàn)cont=zeros(length(r),X_col);fori=length(r)forj=1:X_colcont(i,j)=abs(S(i)/D(i)*P(j,i)*Xtest(51,j));endend%7.3.計(jì)算每個(gè)變量的總貢獻(xiàn)CONTJ=zeros(X_col,1);forj=1:X_colCONTJ(j)=sum(cont(:,j));end%7.4.計(jì)算每個(gè)變量對(duì)Q的貢獻(xiàn)e=Xtest(51,:)*(I-P*P');%選取第60個(gè)樣本來(lái)檢測(cè)哪個(gè)變量出現(xiàn)問(wèn)題。contq=e.^2;%5.繪制貢獻(xiàn)圖figuresubplot(2,1,1);bar(contq,'g');xlabel('變量號(hào)');ylabel('SPE貢獻(xiàn)率%');holdon;subplot(2,1,2);bar(CONTJ,'r');xlabel('變量號(hào)');ylabel('T^2貢獻(xiàn)率%');來(lái)源于：\o"基于PCA的故障診斷方法（matlab）_湯憲宇的博客-CSDN博客_pca故障診斷"基于PCA的故障診斷方法（matlab）_湯憲宇的博客-CSDN博客_pca故障診斷不懂的轉(zhuǎn)原文看，花個(gè)時(shí)間學(xué)學(xué)。2、《普通公路隧道機(jī)電設(shè)施故障檢測(cè)方法研究》中的方法這個(gè)方法就是計(jì)算，看懂這篇論文就行。3、小波分析這個(gè)方法也看論文學(xué)習(xí)，找找代碼。代碼需要下載\o"使用小波進(jìn)行信號(hào)中的特征檢測(cè)：使用小波進(jìn)行特征檢測(cè)的MATLAB代碼和相關(guān)文件-matlab開(kāi)發(fā)_-互聯(lián)網(wǎng)文檔類(lèi)資源-CSDN下載"使用小波進(jìn)行信號(hào)中的特征檢測(cè)：使用小波進(jìn)行特征檢測(cè)的MATLAB代碼和相關(guān)文件-matlab開(kāi)發(fā)_-互聯(lián)網(wǎng)文檔類(lèi)資源-CSDN下載\o"基于小波分析的故障檢測(cè)與診斷碩士論文-基于小波分析的故障檢測(cè)與診斷.rar_-互聯(lián)網(wǎng)文檔類(lèi)資源-CSDN下載"基于小波分析的故障檢測(cè)與診斷碩士論文-基于小波分析的故障檢測(cè)與診斷.rar_-互聯(lián)網(wǎng)文檔類(lèi)資源-CSDN下載兩篇都有，我也沒(méi)下過(guò)，但看標(biāo)題來(lái)說(shuō)。第一個(gè)是代碼，第二個(gè)是論文。至于性能評(píng)價(jià)無(wú)非就是那幾個(gè)指標(biāo)，準(zhǔn)確率、召回率等等。%樣本標(biāo)記為0和1，num為選取前n個(gè)特征的數(shù)據(jù)用于分類(lèi)%需要安裝好SVMfunction[sens,spec,F1,pre,rec,acc]=SEERES(train,trainclass,test,testclass,num)acc=zeros(num,1);sens=zeros(num,1);spec=zeros(num,1);F1=zeros(num,1);pre=zeros(num,1);rec=zeros(num,1);FeatureNumber=zeros(num,1);[len,b]=size(testclass);forn=1:numlabel=trainclass;data=train(:,1:n);testlabel=testclass;testdata=test(:,1:n);model=svmtrain(label,data,'-s0-t0-b1');%默認(rèn)C-SVC類(lèi)型，0線性2RBF,-b會(huì)輸出概率[predictlabel,accuracy,Scores]=svmpredict(testlabel,testdata,model,'-b1');acc(n,1)=accuracy(1,1);FeatureNumber(n,1)=n;tp=0;fn=0;fp=0;tn=0;fory=1:lenifpredictlabel(y,1)==1&&testclass(y,1)==1tp=tp+1;elseifpredictlabel(y,1)==1&&testclass(y,1)==0fp=fp+1;elseifpredictlabel(y,1)==0&&testclass(y,1)==1fn=fn+1;elseifpredictlabel(y,1)==0&&testclass(y,1)==0tn=tn+1;endendsens(n,1)=tp/(tp+fn);spec(n,1)=tn/(tn+fp);pre(n,1)=tp/(tp+fp);rec(n,1)=sens(n,1);F1(n,1)=2*(pre(n,1)*rec(n,1))/(pre(n,1)+rec(n,1));end來(lái)源于：\o"查準(zhǔn)率、召回率、敏感性、特異性和F1-score的計(jì)算及Matlab實(shí)現(xiàn)_黑山白雪m的博客-CSDN博客_召回率和敏感性"查準(zhǔn)率、召回率、敏感性、特異性和F1-score的計(jì)算及Matlab實(shí)現(xiàn)_黑山白雪m的博客-CSDN博客_召回率和敏感性關(guān)于問(wèn)題三：1、LSTM故障診斷\o"基于LSTM的故障診斷_舊日之歌的博客-CSDN博客_故障檢測(cè)lstm"基于LSTM的故障診斷_舊日之歌的博客-CSDN博客_故障檢測(cè)lstm看看上面的這個(gè)博文。2、小波分析與神經(jīng)網(wǎng)絡(luò)這個(gè)方法的代碼在CSDN里需要下載，但是這個(gè)很成熟了，而且已經(jīng)有人用這個(gè)方法做過(guò)診斷齒輪箱的故障。大家可以去看看下面這篇論文?！痘谛〔ń翟牒虰P神經(jīng)網(wǎng)絡(luò)的風(fēng)力發(fā)電機(jī)組齒輪箱故障診斷研究》就是代碼比較難找，但方法應(yīng)該是比較好的。3、深度學(xué)習(xí)故障診斷算法這個(gè)方法python厲害，擅長(zhǎng)編程的可以考慮，用的是殘差收縮網(wǎng)絡(luò)。#!/usr/bin/envpython3#-*-coding:utf-8-*-"""CreatedonSatDec2823:24:052019ImplementedusingTensorFlow1.0.1andKeras2.2.1M.Zhao,S.Zhong,X.Fu,etal.,DeepResidualShrinkageNetworksforFaultDiagnosis,IEEETransactionsonIndustrialInformatics,2019,DOI:10.1109/TII.2019.2943898@author:super_9527"""from__future__importprint_functionimportkerasimportnumpyasnpfromkeras.datasetsimportmnistfromkeras.layersimportDense,Conv2D,BatchNormalization,Activationfromkeras.layersimportAveragePooling2D,Input,GlobalAveragePooling2Dfromkeras.optimizersimportAdamfromkeras.regularizersimportl2fromkerasimportbackendasKfromkeras.modelsimportModelfromkeras.layers.coreimportLambdaK.set_learning_phase(1)#Inputimagedimensionsimg_rows,img_cols=28,28#Thedata,splitbetweentrainandtestsets(x_train,y_train),(x_test,y_test)=mnist.load_data()ifK.image_data_format()=='channels_first':x_train=x_train.reshape(x_train.shape[0],1,img_rows,img_cols)x_test=x_test.reshape(x_test.shape[0],1,img_rows,img_cols)input_shape=(1,img_rows,img_cols)else:x_train=x_train.reshape(x_train.shape[0],img_rows,img_cols,1)x_test=x_test.reshape(x_test.shape[0],img_rows,img_cols,1)input_shape=(img_rows,img_cols,1)#Noiseddatax_train=x_train.astype('float32')/255.+0.5*np.random.random([x_train.shape[0],img_rows,img_cols,1])x_test=x_test.astype('float32')/255.+0.5*np.random.random([x_test.shape[0],img_rows,img_cols,1])print('x_trainshape:',x_train.shape)print(x_train.shape[0],'trainsamples')print(x_test.shape[0],'testsamples')#convertclassvectorstobinaryclassmatricesy_train=keras.utils.to_categorical(y_train,10)y_test=keras.utils.to_categorical(y_test,10)defabs_backend(inputs):returnK.abs(inputs)defexpand_dim_backend(inputs):returnK.expand_dims(K.expand_dims(inputs,1),1)defsign_backend(inputs):returnK.sign(inputs)defpad_backend(inputs,in_channels,out_channels):pad_dim=(out_channels-in_channels)//2inputs=K.expand_dims(inputs,-1)inputs=K.spatial_3d_padding(inputs,((0,0),(0,0),(pad_dim,pad_dim)),'channels_last')returnK.squeeze(inputs,-1)#ResidualShrinakgeBlockdefresidual_shrinkage_block(incoming,nb_blocks,out_channels,downsample=False,downsample_strides=2):residual=incomingin_channels=incoming.get_shape().as_list()[-1]foriinrange(nb_blocks):identity=residualifnotdownsample:downsample_strides=1residual=BatchNormalization()(residual)residual=Activation('relu')(residual)residual=Conv2D(out_channels,3,strides=(downsample_strides,downsample_strides),padding='same',kernel_initializer='he_normal',kernel_regularizer=l2(1e-4))(residual)residual=BatchNormalization()(residual)residual=Activation('relu')(residual)residual=Conv2D(out_channels,3,padding='same',kernel_initializer='he_normal',kernel_regularizer=l2(1e-4))(residual)#Calculateglobalmeansresidual_abs=Lambda(abs_backend)(residual)abs_mean=GlobalAveragePooling2D()(residual_abs)#Calculatescalingcoefficientsscales=Dense(out_channels,activation=None,kernel_initializer='he_normal',kernel_regularizer=l2(1e-4))(abs_mean)scales=BatchNormalization()(scales)scales=Activation('relu')(scales)scales=Dense(out_channels,activation='sigmoid',kernel_regularizer=l2(1e-4))(scales)scales=Lambda(expand_dim_backend)(scales)#Calculatethresholdsthres=keras.layers.multiply([abs_mean,scales])#Softthresholdingsub=keras.layers.subtract([residual_abs,thres])zeros=keras.layers.subtract([sub,sub])n_sub=keras.layers.maximum([sub,zeros])residual=keras.layers.multiply([Lambda(sign_backend)(residual),n_sub])#DownsamplingusingthepooL-sizeof(1,1)ifdownsample_strides>1:identity=AveragePooling2D(pool_size=(1,1),strides=(2,2))(identity)#Zero_paddingtomatchchannelsifin_channels!=out_channels:identity=Lambda(pad_backend,arguments={'in_channels':in_channels,'out_channels':out_channels})(identity)residual=keras.layers.add([residual,identity])returnresidual#defineandtrainamodelinputs=Input(shape=input_shape)net=Conv2D(8,3,padding='same',kernel_initializer='he_normal',kernel_regularizer=l2(1e-4))(inputs)net=residual_shrinkage_block(net,1,8,downsample=True)net=BatchNormalization()(net)net=Activation('relu')(net)net=GlobalAveragePooling2D()(net)outputs=Dense(10,activation='softmax',kernel_initializer='he_normal',kernel_regularizer=l2(1e-4))(net)model=Model(inputs=inputs,outputs=outputs)pile(loss='categorical_crossentropy',optimizer=Adam(),metrics=['accuracy'])model.fit(x_train,y_train,batch_size=100,epochs=5,verbose=1,validation_data=(x_test,y_test))#getresultsK.set_learning_phase(0)DRSN_train_score=model.evaluate(x_train,y_train,batch_size=100,verbose=0)print('Trainloss:',DRSN_train_score[0])print('Trainaccuracy:',DRSN_train_score[1])DRSN_test_score=model.evaluate(x_test,y_test,batch_size=100,verbose=0)print('Testloss:',DRSN_test_score[0])print('Testaccuracy:',DRSN_test_score[1])來(lái)源于：\o"深度學(xué)習(xí)故障診斷算法：殘差收縮網(wǎng)絡(luò)_weixin_47174159的博客-CSDN博客_故障診斷算法"深度學(xué)習(xí)故障診斷算法：殘差收縮網(wǎng)絡(luò)_weixin_47174159的博客-CSDN博客_故障診斷算法性能評(píng)價(jià)借鑒前一問(wèn)的。關(guān)于第四問(wèn)：直接運(yùn)用前面的模型就行，注意論文寫(xiě)好點(diǎn)。第三問(wèn)代碼的更新：LSTM故障診斷matlab代碼：其中第一個(gè)用于讀取并劃分原始數(shù)據(jù)

第二個(gè)用于完成劃分訓(xùn)練集測(cè)試集，特征提取+分類(lèi)等工作//本函數(shù)Input為//interval-數(shù)據(jù)劃分長(zhǎng)度，默認(rèn)為6400，即每6400個(gè)數(shù)據(jù)點(diǎn)劃為一個(gè)樣本//ind_column-通道，默認(rèn)為2，即選取第二個(gè)通道//Output為//label1label2,...,label8，分別為劃分好的8種故障類(lèi)型的樣本%function[label1label2label3label4label5label6label7label8]=read_data_1800_High(interval,ind_column)ifnargin<2ind_column=2;//如果傳遞的實(shí)參小于2個(gè)，默認(rèn)ind_column為2endifnargin<1interval=6400;//默認(rèn)interval=6400endfile_rul='E:\Datasets\PHMdatachallenge\2009PHMSocietyConferenceDataChallenge-gearbox\spur_30hz_High\';//以下為獲取file_rul路徑下.mat格式的所有文件file_folder=fullfile(file_rul);dir_output=dir(fullfile(file_folder,'*.mat'));file_name={dir_}';num_file=max(size(file_name));//num_file為文件數(shù)，本例中num_file=8，8個(gè)文件，分別存儲(chǔ)齒輪箱的8種故障數(shù)據(jù)fori=1:num_filefile=[file_rul,file_name{i}];load(file);[filepath,name,ext]=fileparts(file);raw=eval(name); //每6400個(gè)點(diǎn)劃分為一個(gè)樣本n=1;left_index=1+(n-1)*interval;right_index=n*interval;whileright_index<=size(raw,1)temp=raw(left_index:right_index,ind_column);//eval函數(shù)構(gòu)造label1,label2,...等變量名eval(['label'num2str(i)'(:,n)=temp;']);n=n+1;left_index=1+(n-1)*interval;right_index=n*interval;endendend//讀取數(shù)據(jù)，label1為一個(gè)6400*83的數(shù)組，83為每種故障類(lèi)型所得到的樣本數(shù)[label1,label2,label3,label4,label5,label6,label7,label8]=read_data_1800_High();num_categories=8;//由于matlab中LSTM建模需要，用num2cell函數(shù)將label1轉(zhuǎn)為cell型，label_x_cell為一個(gè)1×83的cell型數(shù)組，每個(gè)cell存儲(chǔ)6400個(gè)數(shù)據(jù)點(diǎn)label1_x_cell=num2cell(label1,1);label2_x_cell=num2cell(label2,1);label3_x_cell=num2cell(label3,1);label4_x_cell=num2cell(label4,1);label5_x_cell=num2cell(label5,1);label6_x_cell=num2cell(label6,1);label7_x_cell=num2cell(label7,1);label8_x_cell=num2cell(label8,1);num_1=length(label1_x_cell);num_2=length(label2_x_cell);num_3=length(label3_x_cell);num_4=length(label4_x_cell);num_5=length(label5_x_cell);num_6=length(label6_x_cell);num_7=length(label7_x_cell);num_8=length(label8_x_cell);//創(chuàng)建用于存儲(chǔ)每種故障類(lèi)型的標(biāo)簽的數(shù)據(jù)結(jié)構(gòu)，由于matlab中l(wèi)stm建模需要，也需要cell型數(shù)據(jù)。例如，label1_y為一個(gè)83×1的cell型數(shù)組，目前其值為空l(shuí)abel1_y=cell(num_1,1);label2_y=cell(num_2,1);label3_y=cell(num_3,1);label4_y=cell(num_4,1);label5_y=cell(num_5,1);label6_y=cell(num_6,1);label7_y=cell(num_7,1);label8_y=cell(num_8,1);//創(chuàng)建故障類(lèi)型的標(biāo)簽，用1,2,3,...,8表示8種故障標(biāo)簽，給對(duì)應(yīng)標(biāo)簽賦值。fori=1:num_1;label1_y{i}='1';endfori=1:num_2;label2_y{i}='2';endfori=1:num_3;label3_y{i}='3';endfori=1:num_4;label4_y{i}='4';endfori=1:num_5;label5_y{i}='5';endfori=1:num_6;label6_y{i}='6';endfori=1:num_7;label7_y{i}='7';endfori=1:num_8;label8_y{i}='8';end//用dividerand函數(shù)將每種故障類(lèi)型的數(shù)據(jù)隨機(jī)劃分為4:1的比例，分別用作訓(xùn)練和測(cè)試[trainInd_label1,~,testInd_label1]=dividerand(num_1,0.8,0,0.2);[trainInd_label2,~,testInd_label2]=dividerand(num_2,0.8,0,0.2);[trainInd_label3,~,testInd_label3]=dividerand(num_3,0.8,0,0.2);[trainInd_label4,~,testInd_label4]=dividerand(num_4,0.8,0,0.2);[trainInd_label5,~,testInd_label5]=dividerand(num_5,0.8,0,0.2);[trainInd_label6,~,testInd_label6]=dividerand(num_6,0.8,0,0.2);[trainInd_label7,~,testInd_label7]=dividerand(num_7,0.8,0,0.2);[trainInd_label8,~,testInd_label8]=dividerand(num_8,0.8,0,0.2);//構(gòu)建每種故障類(lèi)型的訓(xùn)練數(shù)據(jù)xTrain_label1=label1_x_cell(trainInd_label1);yTrain_label1=label1_y(trainInd_label1);xTrain_label2=label2_x_cell(trainInd_label2);yTrain_label2=label2_y(trainInd_label2);xTrain_label3=label3_x_cell(trainInd_label3);yTrain_label3=label3_y(trainInd_label3);xTrain_label4=label4_x_cell(trainInd_label4);yTrain_label4=label4_y(trainInd_label4);xTrain_label5=label5_x_cell(trainInd_label5);yTrain_label5=label5_y(trainInd_label5);xTrain_label6=label6_x_cell(trainInd_label6);yTrain_label6=label6_y(trainInd_label6);xTrain_label7=label7_x_cell(trainInd_label7);yTrain_label7=label7_y(trainInd_label7);xTrain_label8=label8_x_cell(trainInd_label8);yTrain_label8=label8_y(trainInd_label8);//構(gòu)建每種故障類(lèi)型的測(cè)試數(shù)據(jù)xTest_label1=label1_x_cell(testInd_label1);yTest_label1=label1_y(testInd_label1);xTest_label2=label2_x_cell(testInd_label2);yTest_label2=label2_y(testInd_label2);xTest_label3=label3_x_cell(testInd_label3);yTest_label3=label3_y(testInd_label3);xTest_label4=label4_x_cell(testInd_label4);yTest_label4=label4_y(testInd_label4);xTest_label5=label5_x_cell(testInd_label5);yTest_label5=label5_y(testInd_label5);xTest_label6=label6_x_cell(testInd_label6);yTest_label6=label6_y(testInd_label6);xTest_label7=label7_x_cell(testInd_label7);yTest_label7=label7_y(testInd_label7);xTest_label8=label8_x_cell(testInd_label8);yTest_label8=label8_y(testInd_label8);//將每種故障類(lèi)型的數(shù)據(jù)整合，構(gòu)建完整的訓(xùn)練集和測(cè)試集xTrain=[xTrain_label1xTrain_label2xTrain_label3xTrain_label4xTrain_label5xTrain_label6xTrain_label7xTrain_label8];yTrain=[yTrain_label1;yTrain_label2;yTrain_label3;yTrain_label4;yTrain_label5;yTrain_label6;yTrain_label7;yTrain_label8];num_train=size(xTrain,2);xTest=[xTest_label1xTest_label2xTest_label3xTest_label4xTest_label5xTest_label6xTest_label7xTest_label8];yTest=[yTest_label1;yTest_label2;yTest_label3;yTest_label4;yTest_label5;yTest_label6;yTest_label7;yTest_label8];num_test=size(xTest,2);//================================================================================//以下分別對(duì)每個(gè)樣本，提取三種特征：1.瞬時(shí)頻率，2.瞬時(shí)譜熵，3.小波包能量，//上述三種特征后面會(huì)被送入分類(lèi)器中進(jìn)行分類(lèi)，實(shí)驗(yàn)結(jié)果表明，將小波包能量作為特征，//能夠取得最高的分類(lèi)精度//提取瞬時(shí)頻率：用matlab的pspectrum對(duì)每個(gè)樣本進(jìn)行譜分解，再用instfreq函數(shù)計(jì)算瞬時(shí)頻率FreqResolu=25;TimeResolu=0.12;//theoutputofpspectrum'p'containsanestimateoftheshort-term,time-localizedpowerspectrumofx.//Inthiscase,pisofsizeNf×Nt,whereNfisthelengthoffandNtisthelengthoft.[p,f,t]=cellfun(@(x)pspectrum(x,fs,'TimeResolution',TimeResolu,'spectrogram'),xTrain,'UniformOutput',false);instfreqTrain=cellfun(@(x,y,z)instfreq(x,y,z)',p,f,t,'UniformOutput',false);[p,f,t]=cellfun(@(x)pspectrum(x,fs,'TimeResolution',TimeResolu,'spectrogram'),xTest,'UniformOutput',false);instfreqTest=cellfun(@(x,y,z)instfreq(x,y,z)',p,f,t,'UniformOutput',false);//提取瞬時(shí)譜熵：用matlab的pspectrum對(duì)每個(gè)樣本進(jìn)行譜分解，再用pentropy函數(shù)計(jì)算瞬時(shí)頻率[p,f,t]=cellfun(@(x)pspectrum(x,fs,'TimeResolution',TimeResolu,'spectrogram'),xTrain,'UniformOutput',false);pentropyTrain=cellfun(@(x,y,z)pentropy(x,y,z)',p,f,t,'UniformOutput',false);[p,f,t]=cellfun(@(x)pspectrum(x,fs,'TimeResolution',TimeResolu,'spectrogram'),xTest,'UniformOutput',false);pentropyTest=cellfun(@(x,y,z)pentropy(x,y,z)',p,f,t,'UniformOutput',false);//提取小波包能量//num_level=5表示進(jìn)行小波包五層分解，共獲得2^5=32個(gè)值組成的特征向量。num_level=5;index=0:1:2^num_level-1;//wpdec為小波包分解函數(shù)treeTrain=cellfun(@(x)wpdec(x,num_level,'dmey'),xTrain,'UniformOutput',false);treeTest=cellfun(@(x)wpdec(x,num_level,'dmey'),xTest,'UniformOutput',false);fori=1:num_trainforj=1:length(index) //wprcoef為小波系數(shù)重構(gòu)函數(shù)reconstr_coef=wprcoef(treeTrain{i},[num_level,index(j)]);//計(jì)算能量energy(j)=sum(reconstr_coef.^2);endenergyTrain_doule(i,:)=energy;endenergyTrain=num2cell(energyTrain_doule,2);energyTrain=energyTrain';fori=1:num_testforj=1:length