更新代码

LIB/MIDDLE/Anomaly Detection/V1_0_0/anomalyPCA.py

@@ -0,0 +1,119 @@
+#热失控预警：PCA异常指数
+
+import pandas as pd
+import numpy as np
+from scipy.signal import savgol_filter
+from sklearn.preprocessing import RobustScaler
+from sklearn.decomposition import PCA
+
+#筛选特征
+def makedataset1(df_data):
+    df_data=df_data.drop(['Unnamed: 0','总电流[A]','GSM信号','外电压','SOH[%]','开关状态','充电状态','故障等级','故障代码','绝缘电阻','上锁状态','加热状态','单体均衡状态','总输出状态'],axis=1,errors='ignore')
+    df_data=df_data.drop(["单体温度"+str(i) for i in range(1,5)],axis=1,errors='ignore')
+    df_data=df_data.drop(["其他温度"+str(i) for i in range(1,7)],axis=1,errors='ignore')
+    listV=[s for s in list(df_data) if '单体电压' in s]
+    for i in range(1,len(listV)+1):
+        df_data=df_data[(df_data['单体电压'+str(i)]>2200) & (df_data['单体电压'+str(i)]<4800)]
+    df_data=df_data[df_data['SOC[%]']>20]
+    df_data['时间']=[df_data.loc[i,'时间戳'][0:15] for i in df_data.index]
+    data_set=df_data.groupby('时间').mean()
+    for k in data_set.columns:
+        data_set[k]=savgol_filter(data_set[k],3,2)
+    return data_set
+
+#新建统计特征
+def makedataset2(df_data):
+    data_set=makedataset1(df_data)
+    listV=[s for s in list(df_data) if '单体电压' in s]
+    data_set["最低单体电压"]=data_set[["单体电压"+str(i) for i in range(1,len(listV)+1)]].min(axis=1)
+    data_set["最高单体电压"]=data_set[["单体电压"+str(i) for i in range(1,len(listV)+1)]].max(axis=1)
+    data_set["平均单体电压"]=data_set[["单体电压"+str(i) for i in range(1,len(listV)+1)]].mean(axis=1)
+    data_set["最大单体压差"]=[data_set.loc[k,"最高单体电压"]-data_set.loc[k,"最低单体电压"] for k in data_set.index]
+    data_set["低压差"]=[data_set.loc[k,"平均单体电压"]-data_set.loc[k,"最低单体电压"] for k in data_set.index]
+    data_set=data_set.drop(["单体电压"+str(i) for i in range(1,len(listV)+1)],axis=1)
+    return data_set
+
+#标准化
+def process(data_set):
+    features=data_set.columns
+    sX=RobustScaler(copy=True)
+    data_set2=data_set.copy()
+    data_set2.loc[:,features]=sX.fit_transform(data_set2[features])
+    return data_set2
+
+#异常指数函数
+def anomalyScores(originalDF,reducedDF):
+    loss=np.sum((np.array(originalDF)-np.array(reducedDF))**2,axis=1)
+    loss=pd.Series(data=loss,index=originalDF.index)
+    loss=(loss-np.min(loss))/(np.max(loss)-np.min(loss))
+    return loss
+
+#建立PCA模型
+def anomalyPCA(x_train_pro):
+    n_components=4
+    whiten=True
+    random_state=2
+    pca=PCA(n_components=n_components,whiten=whiten,random_state=random_state)
+    pca.fit(x_train_pro)
+    return pca
+
+#判断PCA异常指数
+def transform(df_data_pro,model,df_data):
+    #降维
+    X_train=model.transform(df_data_pro)
+    X_train=pd.DataFrame(data=X_train,index=df_data_pro.index)
+    #还原
+    X_train_inverse=model.inverse_transform(X_train)
+    X_train_inverse=pd.DataFrame(data=X_train_inverse,index=df_data_pro.index)
+    #异常指数
+    anomalyScoresModel=anomalyScores(df_data_pro,X_train_inverse)
+    anomalyScoresModel=savgol_filter(anomalyScoresModel,15,3)
+    df_data2=df_data.copy()
+    df_data2['anomalyScores_'+str(model)]=anomalyScoresModel
+    return df_data2
+
+#判断离群
+def detect_outliers(data,pred,threshold=3):
+    anomaly=data['anomalyScores_PCA(n_components=4, random_state=2, whiten=True)']
+    anomalypred=pred['anomalyScores_PCA(n_components=4, random_state=2, whiten=True)']
+    mean_d=np.mean(anomaly.values)
+    std_d=np.std(anomaly.values)
+    max_score=np.max(anomaly.values)
+    outliers2=pd.DataFrame()
+    for k in anomalypred.index:
+        z_score= (anomalypred[k]-mean_d)/std_d
+        if (np.abs(z_score) >threshold) & (anomalypred[k]>max_score):
+            outliers2=outliers2.append(pred[anomalypred.values==anomalypred[k]])
+    return outliers2
+
+#训练模型
+def train_model(data_train):
+    x_train1=makedataset1(data_train) 
+    x_train2=makedataset2(data_train)  
+    x_train_pro1=process(x_train1) 
+    x_train_pro2=process(x_train2) 
+    pca1=anomalyPCA(x_train_pro1) 
+    pca2=anomalyPCA(x_train_pro2) 
+    res1=transform(x_train_pro1,pca1,x_train1)
+    res2=transform(x_train_pro2,pca2,x_train2)
+    return pca1,pca2,res1,res2
+
+#预测
+def prediction(data_test,pca1,pca2):
+    x_test1=makedataset1(data_test) 
+    x_test2=makedataset2(data_test) 
+    x_test_pro1=process(x_test1) 
+    x_test_pro2=process(x_test2) 
+    pred1=transform(x_test_pro1,pca1,x_test1)
+    pred2=transform(x_test_pro2,pca2,x_test2)
+    return pred1,pred2
+
+#判定异常
+def check_anomaly(outliers1,outliers2):
+    if (len(outliers1)>0) & (len(outliers2)>0):
+        outliers=pd.merge(outliers1,outliers2,on='时间')
+        outliers=outliers[outliers['SOC[%]_x']>45]
+        outliers=outliers.drop(['总电压[V]_y','单体压差_y','SOC[%]_y'],axis=1)
+        return outliers
+    
+

LIB/MIDDLE/Anomaly Detection/V1_0_0/main_anomalyPCA.py

@@ -0,0 +1,42 @@
+#热失控预警：PCA异常指数
+#训练模型
+
+from LIB.BACKEND import DBManager
+dbManager = DBManager.DBManager()
+from LIB.MIDDLE.CellStateEstimation.Common import log
+import pandas as pd
+from anomalyPCA import *
+import joblib
+import datetime
+
+dataSOH = pd.read_excel('sn-20210903.xlsx',sheet_name='sn-20210903')
+fileNames = dataSOH['sn']
+fileNames = list(fileNames)
+l = len(fileNames)
+
+now_time=datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S')   #type: str
+now_time=datetime.datetime.strptime(now_time,'%Y-%m-%d %H:%M:%S')     #type: datetime
+start_time=now_time-datetime.timedelta(days=365)
+end_time=str(now_time)
+start_time=str(start_time)
+
+mylog=log.Mylog('log.txt','error')
+mylog.logcfg()
+
+for k in range(l): 
+    try: 
+        sn = fileNames[k]
+        df_data = dbManager.get_data(sn=sn, start_time=start_time, end_time=end_time, data_groups=['bms'])
+        data_train = df_data['bms']
+        
+        if len(data_train)>0:
+            pca1,pca2,res1,res2=train_model(data_train) 
+            joblib.dump(pca1,'pca1_'+sn+'.m')  
+            joblib.dump(pca2,'pca2_'+sn+'.m')  
+            res1.to_csv('res1_'+sn+'.csv',encoding='gbk')
+            res2.to_csv('res2_'+sn+'.csv',encoding='gbk')
+    
+    except Exception as e:
+        print(repr(e))
+        mylog.logopt(sn,e)
+        pass