补充9月数据

入模数据/1.py

@@ -0,0 +1,121 @@
+import xgboost as xgb
+import pandas as pd
+import os
+from sklearn.metrics import r2_score
+from sklearn.model_selection import train_test_split
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+mpl.rcParams['font.sans-serif']=['kaiti']
+
+pd.set_option('display.width',None)
+
+def hf_season(x):
+    list1= []
+    for i in range(1,13):
+        if x.loc[f'2021-{i}'].mean() >= x.describe()['75%']:
+            list1.append(i)
+    return list1
+
+
+
+def season(x):
+    if str(x)[5:7] in ('06','07','08','12','01','02'):
+        return 1
+    else:
+        return 0
+def month(x):
+    if str(x)[5:7] in ('08','09','10','12','01','02'):
+        return 1
+    else:
+        return 0
+def normal(nd):
+    high = nd.describe()['75%'] + 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    low = nd.describe()['25%'] - 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    return nd[(nd<high)&(nd>low)]
+
+
+data = pd.read_excel(r'C:\python-project\pytorch3\入模数据\杭州数据.xlsx',index_col='dtdate')
+data.index = pd.to_datetime(data.index,format='%Y-%m-%d')
+data = data.loc[normal(data['售电量']).index]
+# for i in range(1,13):
+# plt.plot(range(len(data['售电量'][f'2022-{i}'])),data['售电量'][f'2022-{i}'])
+# plt.show()
+print(data['售电量']['2022-9'])
+plt.plot(range(len(data['售电量']['2022-7'])),data['售电量']['2022-7'])
+plt.plot(range(len(data['售电量']['2022-7']),len(data['售电量']['2022-7'])+len(data['售电量']['2023-7'])),data['售电量']['2023-7'])
+# plt.plot(range(len(data['售电量'][['2022-9','2023-9']])),data['售电量'][['2022-9','2023-9']])
+plt.show()
+
+# print(hf_season(data.loc['2021']['售电量']))
+
+data['month'] = data.index.strftime('%Y-%m-%d').str[6]
+data['month'] = data['month'].astype('int')
+data['season'] = data.index.map(season)
+print(data.head(50))
+
+df_eval = data.loc['2023-7']
+df_train = data.loc['2021-1':'2023-6']
+# df_train = df[500:850]
+print(len(df_eval),len(df_train),len(data))
+
+print(data.drop(columns='city_name').corr(method='pearson')['售电量'])
+
+df_train = df_train[['tem_max','tem_min','24ST','rh','rh_max','prs','prs_max','prs_min','售电量','month','holiday','season']]
+
+
+# IQR = df['售电量'].describe()['75%'] - df['售电量'].describe()['25%']
+# high = df['售电量'].describe()['75%'] + 1.5*IQR
+# low = df['售电量'].describe()['25%'] - 1.5*IQR
+# print('异常值数量：',len(df[(df['售电量'] >= high) | (df['售电量'] <= low)]))
+#
+# df_train = df_train[(df['售电量'] <= high) & (df['售电量'] >= low)]
+
+
+X = df_train[['tem_max','tem_min','24ST','holiday','season']]
+X_eval = df_eval[['tem_max','tem_min','24ST','holiday','season']]
+y = df_train['售电量']
+print(y.describe())
+# best_goal = 1
+# best_i = {}
+# for i in range(400):
+
+x_train,x_test,y_train,y_test = train_test_split(X,y,test_size=0.15,random_state=42)
+model = xgb.XGBRegressor(max_depth=6, learning_rate=0.05, n_estimators=150)
+model.fit(x_train,y_train)
+
+y_pred = model.predict(x_test)
+result_test = pd.DataFrame({'test':y_test,'pred':y_pred},index=y_test.index)
+
+# 指标打印
+print(abs(y_test - y_pred).mean() / y_test.mean())
+eval_pred = model.predict(X_eval)
+
+result_eval = pd.DataFrame({'eval':df_eval['售电量'],'pred':eval_pred},index=df_eval['售电量'].index)
+
+print((result_eval['eval'].sum()-result_eval['pred'].sum())/result_eval['eval'].sum())
+
+goal = (result_eval['eval'][-3:].sum()-result_eval['pred'][-3:].sum())/result_eval['eval'].sum()
+print('goal:',goal)
+
+goal2 = (result_eval['eval'][-23:].sum()-result_eval['pred'][-23:].sum())/result_eval['eval'].sum()
+
+print('goal2:',goal2)
+print(result_eval)
+print('r2:',r2_score(y_test,y_pred))
+    # if abs(goal) < best_goal:
+    #     best_goal = abs(goal)
+    #     best_i['best_i'] = i
+    #     x = goal2
+# print(best_i,best_goal,x)
+
+
+
+# result_eval.to_csv(r'C:\Users\user\Desktop\9月各地市日电量预测结果\杭州.csv')
+# with open(r'C:\Users\user\Desktop\9月各地市日电量预测结果\偏差率.txt','a',encoding='utf-8') as f:
+#     f.write(f'杭州月末3天偏差率：{round(goal,5)},9号-月底偏差率：{round(goal2,5)}\n')
+# # 保存模型
+# model.save_model('hangzhou.bin')
+# loaded_model = xgb.XGBRegressor()
+# loaded_model.load_model('hangzhou.bin')
+# model.predict(X_eval)
+

区域电量19年至今数据.py

@@ -0,0 +1,177 @@
+import xgboost as xgb
+import pandas as pd
+import os
+from sklearn.metrics import r2_score
+from sklearn.model_selection import train_test_split
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+import datetime
+import math
+from sklearn.preprocessing import LabelEncoder
+mpl.rcParams['font.sans-serif']=['kaiti']
+pd.set_option('display.width',None)
+
+
+
+
+
+
+def season(x):
+    if str(x)[5:7] in ['07', '08']:
+        return 2
+    elif str(x)[5:7] in ['01', '02', '03', '06', '09', '11', '12']:
+        return 1
+    elif str(x)[5:7] in['04', '05', '10']:
+        return 0
+def normal(nd):
+    high = nd.describe()['75%'] + 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    low = nd.describe()['25%'] - 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    return nd[(nd<high)&(nd>low)].index
+
+# df = pd.read_excel(r'C:\Users\鸽子\Desktop\杭州19年至今日电量及气象数据.xlsx',sheet_name=0)
+# df_elec = pd.read_excel(r'C:\Users\鸽子\Desktop\杭州19年至今日电量及气象数据.xlsx',sheet_name=1)
+# df_elec.columns = df_elec.columns.map(lambda x:x.strip())
+# df_elec['售电量'] = df_elec['售电量']/10000
+# df.columns = df.columns.map(lambda x:x.strip())
+# df = df[['dtdate','tem_max','tem_min']]
+# # print(df.head())
+# # print(df_elec.head())
+#
+# merge_df = pd.merge(df_elec,df,left_on='pt_date',right_on='dtdate')[['pt_date','tem_max','tem_min','售电量']]
+# merge_df.set_index('pt_date',inplace=True)
+# merge_df.index = pd.to_datetime(merge_df.index,format='%Y%m%d')
+#
+#
+# merge_df['month'] = merge_df.index.strftime('%Y-%m-%d').str[5:7]
+# merge_df['month'] = merge_df['month'].astype('int')
+# merge_df.to_csv('杭州入模数据.csv',encoding='gbk')
+data = pd.read_csv(r'杭州入模数据.csv',encoding='gbk')
+
+data.drop_duplicates(subset='pt_date',inplace=True)
+
+data.set_index('pt_date',inplace=True)
+data.index = pd.to_datetime(data.index)
+
+print(data.loc['2023-07'])
+def jq(y,x):
+    a=365.242 * (y - 1900) + 6.2 + 15.22 * x - 1.9 * math.sin(0.262 * x)
+    return datetime.date(1899,12,31)+datetime.timedelta(days=int(a))
+# print(jq(2020,0))
+jq_list=['小寒', '大寒', '立春', '雨水', '惊蛰', '春分', '清明', '谷雨', '立夏', '小满', '芒种', '夏至', '小暑', '大暑', '立秋', '处暑', '白露', '秋分', '寒露', '霜降', '立冬', '小雪', '大雪','冬至']
+jq_dict={}
+for j in range(2019,2024):
+    for i in range(24):
+        jq_dict[jq(j,i).strftime('%Y-%m-%d')]=jq_list[i]
+# print(jq_dict)
+
+data['24ST']=data.index
+data['24ST']=data['24ST'].astype('string').map(jq_dict)
+data['24ST'].fillna(method='ffill',inplace=True)
+data['24ST'].fillna('冬至',inplace=True)
+
+
+# data为数据集  product_tags为需要编码的特征列(假设为第一列)
+le = LabelEncoder()
+data['24ST'] = le.fit_transform(data['24ST'])
+data = data.loc[normal(data['售电量'])]
+
+data['season'] = data.index.map(season)
+print(data['售电量'].describe())
+print(data)
+# list2 = []
+# list0 = []
+# list1 = []
+# for i in ('01','02','03','04','05','06','07','08','09','10','11','12'):
+#     month_index = df.index.strftime('%Y-%m-%d').str[5:7] == f'{i}'
+#     if df.loc[month_index]['售电量'].mean() >= df['售电量'].describe()['75%']:
+#         list2.append(i)
+#     elif df.loc[month_index]['售电量'].mean() <= df['售电量'].describe()['25%']:
+#         list0.append(i)
+#     else:
+#         list1.append(i)
+# print(list0,list1,list2)
+
+
+
+# data = pd.read_excel(r'C:\python-project\pytorch3\入模数据\杭州数据.xlsx',index_col='dtdate')
+# data.index = pd.to_datetime(data.index,format='%Y-%m-%d')
+# data = data.loc[normal(data['售电量']).index]
+# plt.plot(range(len(data['售电量']['2021':'2022'])),data['售电量']['2021':'2022'])
+# plt.show()
+
+# # print(hf_season(data.loc['2021']['售电量']))
+
+# data['month'] = data.index.strftime('%Y-%m-%d').str[6]
+# data['month'] = data['month'].astype('int')
+# data['season'] = data.index.map(season)
+# print(data.head(50))
+#
+df_eval = data.loc['2023-9']
+df_train = data.loc['2019-1':'2023-8']
+print(len(df_train),len(df_eval))
+plt.plot(range(len(data.loc['2019-1':'2023-9'])),data.loc['2019-1':'2023-9'])
+plt.show()
+# df_train = df[500:850]
+print(len(df_eval),len(df_train),len(data))
+
+print(data.corr(method='pearson')['售电量'])
+
+df_train = df_train[['tem_max','tem_min','24ST','售电量','season']]
+
+
+# IQR = df['售电量'].describe()['75%'] - df['售电量'].describe()['25%']
+# high = df['售电量'].describe()['75%'] + 1.5*IQR
+# low = df['售电量'].describe()['25%'] - 1.5*IQR
+# print('异常值数量：',len(df[(df['售电量'] >= high) | (df['售电量'] <= low)]))
+#
+# df_train = df_train[(df['售电量'] <= high) & (df['售电量'] >= low)]
+
+
+X = df_train[['tem_max','tem_min','season','24ST']]
+X_eval = df_eval[['tem_max','tem_min','season','24ST']]
+y = df_train['售电量']
+print(y.describe())
+# best_goal = 1
+# best_i = {}
+# for i in range(400):
+
+x_train,x_test,y_train,y_test = train_test_split(X,y,test_size=0.2,random_state=42)
+model = xgb.XGBRegressor(max_depth=6, learning_rate=0.05, n_estimators=150)
+model.fit(x_train,y_train)
+
+y_pred = model.predict(x_test)
+result_test = pd.DataFrame({'test':y_test,'pred':y_pred},index=y_test.index)
+
+# 指标打印
+print(abs(y_test - y_pred).mean() / y_test.mean())
+eval_pred = model.predict(X_eval)
+
+result_eval = pd.DataFrame({'eval':df_eval['售电量'],'pred':eval_pred},index=df_eval['售电量'].index)
+
+print((result_eval['eval'].sum()-result_eval['pred'].sum())/result_eval['eval'].sum())
+
+goal = (result_eval['eval'][-3:].sum()-result_eval['pred'][-3:].sum())/result_eval['eval'].sum()
+print('goal:',goal)
+
+goal2 = (result_eval['eval'][-23:].sum()-result_eval['pred'][-23:].sum())/result_eval['eval'].sum()
+
+print('goal2:',goal2)
+print(result_eval)
+print('r2:',r2_score(y_test,y_pred))
+#     if abs(goal) < best_goal:
+#         best_goal = abs(goal)
+#         best_i['best_i'] = i
+#         x = goal2
+# print(best_i,best_goal,x)
+
+
+
+# result_eval.to_csv(r'C:\Users\user\Desktop\9月各地市日电量预测结果\杭州.csv')
+# with open(r'C:\Users\user\Desktop\9月各地市日电量预测结果\偏差率.txt','a',encoding='utf-8') as f:
+#     f.write(f'杭州月末3天偏差率：{round(goal,5)},9号-月底偏差率：{round(goal2,5)}\n')
+# # 保存模型
+# model.save_model('hangzhou.bin')
+# loaded_model = xgb.XGBRegressor()
+# loaded_model.load_model('hangzhou.bin')
+# model.predict(X_eval)
+

各地级市日电量模型/杭州.py

@@ -37,7 +37,7 @@ def normal(nd):
 data = pd.read_excel(r'C:\python-project\pytorch3\入模数据\杭州数据.xlsx',index_col='dtdate')
 data.index = pd.to_datetime(data.index,format='%Y-%m-%d')
 data = data.loc[normal(data['售电量']).index]
-plt.plot(range(len(data)),data['售电量'])
+plt.plot(range(len(data['售电量']['2021':'2022'])),data['售电量']['2021':'2022'])
 plt.show()
 
 # print(hf_season(data.loc['2021']['售电量']))
@@ -47,12 +47,12 @@ data['month'] = data['month'].astype('int')
 data['season'] = data.index.map(season)
 print(data.head(50))
 
-df_eval = data.loc['2023-9']
-df_train = data.loc['2021-1':'2023-8']
+df_eval = data.loc['2023-7']
+df_train = data.loc['2021-1':'2023-6']
 # df_train = df[500:850]
 print(len(df_eval),len(df_train),len(data))
 
-print(data.drop(columns='city_name').corr(method='pearson')['season'])
+print(data.drop(columns='city_name').corr(method='pearson')['售电量'])
 
 df_train = df_train[['tem_max','tem_min','24ST','rh','rh_max','prs','prs_max','prs_min','售电量','month','holiday','season']]
 
@@ -73,7 +73,7 @@ print(y.describe())
 # best_i = {}
 # for i in range(400):
 
-x_train,x_test,y_train,y_test = train_test_split(X,y,test_size=0.15,random_state=209)
+x_train,x_test,y_train,y_test = train_test_split(X,y,test_size=0.15,random_state=42)
 model = xgb.XGBRegressor(max_depth=6, learning_rate=0.05, n_estimators=150)
 model.fit(x_train,y_train)
 
@@ -85,14 +85,16 @@ print(abs(y_test - y_pred).mean() / y_test.mean())
 eval_pred = model.predict(X_eval)
 
 result_eval = pd.DataFrame({'eval':df_eval['售电量'],'pred':eval_pred},index=df_eval['售电量'].index)
-# print(result_eval)
+
 print((result_eval['eval'].sum()-result_eval['pred'].sum())/result_eval['eval'].sum())
 
 goal = (result_eval['eval'][-3:].sum()-result_eval['pred'][-3:].sum())/result_eval['eval'].sum()
 print('goal:',goal)
 
 goal2 = (result_eval['eval'][-23:].sum()-result_eval['pred'][-23:].sum())/result_eval['eval'].sum()
+
 print('goal2:',goal2)
+print(result_eval)
 print('r2:',r2_score(y_test,y_pred))
     # if abs(goal) < best_goal:
     #     best_goal = abs(goal)

浙江电压等级电量/电压等级_输出为3.py

@@ -0,0 +1,165 @@
+import numpy as np
+import pandas as pd
+import torch
+from torch import nn
+from multiprocessing import Pool
+import matplotlib.pyplot as plt
+import os
+os.environ["KMP_DUPLICATE_LIB_OK"]="TRUE"
+DAYS_FOR_TRAIN = 10
+torch.manual_seed(42)
+class LSTM_Regression(nn.Module):
+
+    def __init__(self, input_size, hidden_size, output_size=1, num_layers=2):
+        super().__init__()
+
+        self.lstm = nn.LSTM(input_size, hidden_size, num_layers)
+        self.fc = nn.Linear(hidden_size, output_size)
+
+    def forward(self, _x):
+        x, _ = self.lstm(_x)  # _x is input, size (seq_len, batch, input_size)
+        s, b, h = x.shape  # x is output, size (seq_len, batch, hidden_size)
+        x = x.view(s * b, h)
+        x = self.fc(x)
+        x = x.view(s, b, -1)  # 把形状改回来
+        return x
+
+
+def create_dataset(data, days_for_train=5) -> (np.array, np.array):
+    dataset_x, dataset_y = [], []
+    for i in range(len(data) - days_for_train-3):
+        _x = data[i:(i + days_for_train)]
+        dataset_x.append(_x)
+        dataset_y.append(data[i + days_for_train:i + days_for_train+3])
+    return (np.array(dataset_x), np.array(dataset_y))
+
+def normal(nd):
+    high = nd.describe()['75%'] + 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    low = nd.describe()['25%'] - 1.5*(nd.describe()['75%']-nd.describe()['25%'])
+    return nd[(nd<high)&(nd>low)]
+
+def run(file_dir,excel):
+    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')
+
+    data = pd.read_excel(os.path.join(file_dir,excel), sheet_name=0,index_col=' stat_date ')
+
+    data.columns = data.columns.map(lambda x: x.strip())
+
+    data.index = pd.to_datetime(data.index,format='%Y%m%d')
+    data.sort_index(inplace=True)
+    print(data.head())
+    data = data.loc['2021-01':'2023-09']
+
+    data.drop(columns=[i for i in data.columns if (data[i] == 0).sum() / len(data) >= 0.5], inplace=True)  # 去除0值列
+    print('len(data):', len(data))
+    list_app = []
+    for industry in data.columns:
+        df = data[industry]
+        df = df[df.values != 0]  # 去除0值行
+        df = normal(df)
+        df = df.astype('float32').values  # 转换数据类型
+
+
+        # 标准化到0~1
+        max_value = np.max(df)
+        min_value = np.min(df)
+        df = (df - min_value) / (max_value - min_value)
+
+        dataset_x, dataset_y = create_dataset(df, DAYS_FOR_TRAIN)
+        print('len(dataset_x:)', len(dataset_x))
+
+        # 划分训练集和测试集
+        train_size = len(dataset_x) - 3
+        train_x = dataset_x[:train_size]
+        train_y = dataset_y[:train_size]
+
+        # 将数据改变形状，RNN 读入的数据维度是 (seq_size, batch_size, feature_size)
+        train_x = train_x.reshape(-1, 1, DAYS_FOR_TRAIN)
+        train_y = train_y.reshape(-1, 1, 3)
+
+        # 转为pytorch的tensor对象
+        train_x = torch.from_numpy(train_x).to(device)
+        train_y = torch.from_numpy(train_y).to(device)
+
+        model = LSTM_Regression(DAYS_FOR_TRAIN, 32, output_size=3, num_layers=2).to(device)  # 导入模型并设置模型的参数输入输出层、隐藏层等
+
+
+        train_loss = []
+        loss_function = nn.MSELoss()
+        optimizer = torch.optim.Adam(model.parameters(), lr=0.005, betas=(0.9, 0.999), eps=1e-08, weight_decay=0)
+        for i in range(1500):
+            out = model(train_x)
+            loss = loss_function(out, train_y)
+            loss.backward()
+            optimizer.step()
+            optimizer.zero_grad()
+            train_loss.append(loss.item())
+            # print(loss)
+        # 保存模型
+        # torch.save(model.state_dict(),save_filename)
+        # torch.save(model.state_dict(),os.path.join(model_save_dir,model_file))
+
+        # for test
+        model = model.eval()  # 转换成测试模式
+        # model.load_state_dict(torch.load(os.path.join(model_save_dir,model_file)))  # 读取参数
+        dataset_x = dataset_x.reshape(-1, 1, DAYS_FOR_TRAIN)  # (seq_size, batch_size, feature_size)
+        dataset_x = torch.from_numpy(dataset_x).to(device)
+
+        pred_test = model(dataset_x)  # 全量训练集
+        # 模型输出 (seq_size, batch_size, output_size)
+        pred_test = pred_test.view(-1)
+        pred_test = np.concatenate((np.zeros(DAYS_FOR_TRAIN), pred_test.cpu().detach().numpy()))
+
+        # plt.plot(pred_test, 'r', label='prediction')
+        # plt.plot(df, 'b', label='real')
+        # plt.plot((train_size, train_size), (0, 1), 'g--')  # 分割线 左边是训练数据 右边是测试数据的输出
+        # plt.legend(loc='best')
+        # plt.show()
+
+
+        # 创建测试集
+        # result_list = []
+        # 以x为基础实际数据，滚动预测未来3天
+        x = torch.from_numpy(df[-14:-4]).to(device)
+        pred = model(x.reshape(-1,1,DAYS_FOR_TRAIN)).view(-1).detach().numpy()
+
+
+        # for i in range(3):
+        #     next_1_8 = x[1:]
+        #     next_9 = model(x.reshape(-1,1,DAYS_FOR_TRAIN))
+        #     # print(next_9,next_1_8)
+        #     x = torch.concatenate((next_1_8, next_9.view(-1)))
+        #     result_list.append(next_9.view(-1).item())
+
+
+        # 反归一化
+        pred = pred * (max_value - min_value) + min_value
+        df = df * (max_value - min_value) + min_value
+
+        print(pred)
+        # 打印指标
+        print(abs(pred - df[-3:]).mean() / df[-3:].mean())
+        result_eight = pd.DataFrame({'pred': np.round(pred,1),'real': df[-3:]})
+        target = (result_eight['pred'].sum() - result_eight['real'].sum()) / df[-31:].sum()
+        result_eight['loss_rate'] = round(target, 5)
+        result_eight['industry'] = industry
+        list_app.append(result_eight)
+        print(target)
+        print(result_eight)
+    final_df = pd.concat(list_app,ignore_index=True)
+    final_df.to_csv('市行业电量.csv',encoding='gbk')
+    print(final_df)
+
+        # result_eight.to_csv(f'./月底预测结果/9月{excel[:2]}.txt', sep='\t', mode='a')
+        # with open(fr'./偏差/9月底偏差率.txt', 'a', encoding='utf-8') as f:
+        #     f.write(f'{excel[:2]}{industry}:{round(target, 5)}\n')
+
+if __name__ == '__main__':
+    file_dir = r'C:\Users\user\Desktop\浙江各地市分电压日电量数据'
+
+    run(file_dir,'杭州.xlsx')
+    # p = Pool(4)
+    # for excel in os.listdir(file_dir):
+    #     p.apply_async(func=run,args=(file_dir,excel))
+    # p.close()
+    # p.join()