【Python项目实战】基于时间卷积网络(Temporal Convolution Network ,TCN)的发动机剩余寿命预测_tcn预测剩余寿命

摘要航空发动机结构复杂,状态变量多且相互之间存在着严重非线性特征,传统的基于物理失效模型的方法难以精确地预测发动机的剩余寿命(RUL)。针对此问题,采用时间卷积网络(Temporal Convolution Network ,TCN)作为一种最新出现的序列神经网络,被证明在序列数据预测上有良好的效果。采用TCN实现对发动机剩余寿命进行预测,预测过程通过建立退化模型,给每个训练样本添加RUL标签;将特征输入构建的卷积神经网络得到剩余寿命的预测值。为了验证方法的有效性,在NASA提供的涡轮风扇发动机仿真数据集(C-MAPSS)上进行了测试,结果表明采用TCN算法拥有更高的精度。

航空发动机是航空器最重要的设备之一,为航空器提供动力来源,也是日常航空维护作业中的主要关注对象。有研究表明，飞机机械故障导致的飞行事故数量仅低于飞行员操作失误引起的飞行事故数,而发动机故障是机械故障里占比最高的一项,这种现象在军用飞行器里尤为明显。预测发动机剩余寿命对于及时开展视情维修、避免飞行安全事故具有重要的意义。从大体上来分,目前常用的寿命预测方法主要为两类，即基于发动机的运行规律采用物理模型的方法进行预测;从数据挖掘的角度出发采用数据驱动的方法。由于航空发动机包含许多部件,各部件功能差别较大,结构组成复杂,并且发动机失效故障种类繁多,常常是多种故障的复合出现,且各个变量之间的耦合程度很深,从故障表征现象中难以直接和某个变量对应起来,因此依靠原始数据确定特征指标建立精确物理失效模型较为困难。传统的数据驱动预测方法如依据统计性能退化监测数据进行寿命预测、采用退化特征相似性的寿命预测[以及采用非线性维纳过程进行寿命预测的建模研究等,是从统计的角度出发,依据监测数据,根据统计性能指标,结合历史数据得到预测结果。这些方法更多的是从整体上寻找规律,再从整体上进行预测,而事实上,个体与个体之间的差别不可忽略,以整体规律衡量个体属性的准确度不高。
随着深度学习的兴起,数据驱动的方法迎来了新的选择,目前深度学习在图像识别和文本分析等领域都显现了较高的应用价值。近年来,深度学习的方法在其他领域也得到了广泛的应用,长短记忆神经网络、深度置信网络,卷积神经网络等在故障诊断、寿命预测方面均有较好的应用。通过这种基于数据驱动的方法通过挖掘数据和故障状态之间的隐含关系进行预测,对物理机理依赖程度低，可以直接将原始采样数据作为输入,通过提取出重要的特征用于预测,其通用性更强。

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
from sklearn.metrics import r2_score ,mean_squared_error , mean_absolute_error

import warnings
from tcn import TCN, tcn_full_summary
from tensorflow.keras.layers import Dense
from tensorflow.keras.models import Sequential

warnings.filterwarnings("ignore")

np.set_printoptions(suppress=True)

import numpy as np

def mape(actual, pred):
    actual, pred = np.array(actual), np.array(pred)
    return np.mean(np.abs((actual - pred) / actual))*100


# %%加入外生变量进行预测
def data_main(df_train, df_test, chosen_scale):
    data_train = df_train.copy()
    data_test = df_test.copy()
    # 首先确定测试集的长度
    x_train_con = data_train.iloc[:, :-1].values
    x_test_con = data_test.iloc[:, :-1].values
    y_train = data_train.iloc[:, -1].values
    y_test = data_test.iloc[:, -1].values
    # 先用训练集的数据计算标准化的指标，然后对整个数据集进行转换
    if chosen_scale == "NormalizedScaler":
        scale_con = NormalizedScaler(min=np.min(x_train_con, axis=0), max=np.max(x_train_con, axis=0))
        scale_label = NormalizedScaler(min=np.min(y_train), max=np.max(y_train))
    elif chosen_scale == "StandardScaler":
        scale_con = StandardScaler(mean=np.mean(x_train_con, axis=0), std=np.std(x_train_con, axis=0))
        scale_label = StandardScaler(mean=np.mean(y_train), std=np.std(y_train))
    else:
        print("please chosen a scale between NormalizedScaler and StandardScaler")
        return None
    x_train_con_nor = scale_con.transform(x_train_con)
    x_test_con_nor = scale_con.transform(x_test_con)
    y_train_nor = scale_label.transform(y_train)
    y_test_nor = scale_label.transform(y_test)
    return scale_con, scale_label, x_train_con_nor, x_test_con_nor, y_train_nor, y_test_nor


class StandardScaler:
    def __init__(self, mean, std):
        self.mean = mean
        self.std = std

    def transform(self, data):
        return (data - self.mean) / self.std

    # 主要用于后面进行单变量的损失函数操作
    def inverse_transform(self, data):
        return data * self.std + self.mean


# 进行归一化
class NormalizedScaler:
    def __init__(self, min, max):
        self.min = min
        self.max = max

    def transform(self, data):
        return (data - self.min) / (self.max - self.min)

    def inverse_transform(self, data):
        return data * (self.max - self.min) + self.min



if __name__ == "__main__":

    # %%首先读取数据
    train_data = pd.read_excel("data.xlsx", sheet_name="train")
    test_data = pd.read_excel("data.xlsx", sheet_name="test")
    test_label = pd.read_csv('RUL_FD001.txt', sep='\s+', header=None, names=['RUL'])
    test_data.drop_duplicates(subset=["引擎"], keep="last", inplace=True)
    test_data.reset_index(inplace=True, drop=True)
    test_data["y"] = test_label["RUL"].values
    train_data_group = train_data.groupby("引擎")
    train_data_list = []
    for yinqing, tmp in train_data_group:
        tmp["y"] = tmp["周期"].max() - tmp["周期"]
        train_data_list.append(tmp)
    train_data = pd.concat(train_data_list)
    use_columns = [2, 3, 4, 7, 8, 9, 11, 12, 13, 14, 15, 17, 20, 21, "y"]
    train_data = train_data[use_columns]
    test_data = test_data[use_columns]

    scale_con, scale_label, x_train_con_nor, x_test_con_nor, y_train_nor, y_test_nor = data_main(train_data, test_data,
                                                                                                 "StandardScaler")

# TCN
    batch_size, time_steps, input_dim = None, 14, 1
    tcn_layer = TCN(input_shape=(time_steps, input_dim))
    m = Sequential([
        tcn_layer,
        Dense(1)
    ])

    m.compile(optimizer='adam', loss='mse')

    tcn_full_summary(m, expand_residual_blocks=False)

    x_train_con_nor=np.reshape(x_train_con_nor,[np.size(x_train_con_nor,0),np.size(x_train_con_nor,1),1])
    m.fit(x_train_con_nor, np.ravel(y_train_nor), epochs=1000, validation_split=0.2)

    x_test_con_nor_constant=np.reshape(x_test_con_nor,[np.size(x_test_con_nor,0),np.size(x_test_con_nor,1),1])
    tcn_test_pred = m.predict(x_test_con_nor_constant)
    # 对标准化之后预测的结果进行还原
    y_pre_test_inv = scale_label.inverse_transform(tcn_test_pred).reshape(-1, 1)
    y_true_test_inv = scale_label.inverse_transform(y_test_nor).reshape(-1, 1)
    error=np.abs(tcn_test_pred-y_test_nor)
    np.savetxt('arraypred.txt',y_pre_test_inv,delimiter=',')
    np.savetxt('arraytest.txt', y_true_test_inv, delimiter=',')
    np.savetxt('error.txt', error, delimiter=',')
    print("****TCN****")

    print('R2:',r2_score(tcn_test_pred ,y_test_nor))
    print('MSE:', mean_squared_error(tcn_test_pred ,y_test_nor))
    print('MAE:', mean_absolute_error(tcn_test_pred ,y_test_nor))
    print('MAPE:',mape(y_test_nor ,tcn_test_pred))

    fig1=plt.figure(figsize=(8,6))
    plt.plot(y_pre_test_inv)
    plt.plot(y_true_test_inv)
    plt.legend('TCN-Predict data','Real data')


1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133

**指标RMSE 达到0.001 **

安装TCN库

具体说明：https://github.com/philipperemy/keras-tcn

代码下载：https://download.csdn.net/download/qq_45047246/86059805