• 2023年 国赛 数学建模C 基于遗传算法和神经网络的销量定价模型

    一、写在开头

            阅读者可能需要先阅读2023年国赛C题才能读懂下面的内容。

            文章着重于解题方向指引和经历分享，只解释部分核心代码。        

            10月11日更新下，国家奖二等奖。

    二、内容概述

            刚刚做完比赛，对这段经历和对问题的处理方法进行下记录。

    三、个人经历

            今年大三，第一次参加数学建模，队友是大二9月份确定好的，徐同学和周同学（两个人性格都超级好）。大二上的时候学院有一门课就是数模，没怎么认真听，不过最后的大作业比较认真的拿sklearn的BP神经网络和k-means写了一个关于中草药的国赛题。后来学院有组织数模培训，但是听了听就是机器学习，都学过就没有去。赛前没有准备太多，为了防止手生也为了CSP在leecode上刷了二三十道dp。

            周同学不会用python，沟通能力很好，负责论文（是计科的，不会写但是能看懂代码，word很会用也很认真），徐同学和我一个专业，代码不如我但是自学能力、抗压能力沟通能力很强，也很认真，他论文和模型都做了一些，是队员1所以也充当队长，我负责代码和模型（最后提交的代码都是我敲的）。

            我们第一天出题很快就确定了选C题，因为另外两个都读不懂，关于物理的问题也不是我们擅长的（全员网计学院）。开题晚上进行了数据处理，没有做题，11点回去休息。第二天因为读取数据太慢（附件2 87万条数据）花了一上午没做题把处理后数据存到了文件里，当时觉得特别焦虑一上午啥也没写，但是后来做题因为提前存了文件节约了相当多的时间。然后第一天做完了第一题，第二题开了个头。第二天白天做完了第二题，第三题没碰，晚上才开始写第三题，我的队友们在我写第二题的时候在看第三和第四题，当时感觉进度相当慢，有点焦虑，感觉C题比感觉到的难。然后通宵写第三题，我现学python解最优化，但是题目是非凸的，相当不好解。队友用SPSS预测了销量，然后我再继续处理，第二天上午做完了，解最优化用的是现学的遗传算法（之前听说过也了解他的原理，但是是第一次实用）。第四题是开放性问题，我没太重视，写完第三题感觉任务不太多了，摆烂状态帮忙改论文。那天下午继续修改论文，但是下午四点我突然发现我的第三问有一个点写的相当有问题，徐同学觉得别改了，没时间了。但是我感觉时间还够，和他说了我可以后他没再说啥，继续改论文，我就开始高度集中精神修改第四问，然后改的时候也挺焦虑的，特别怕队友怪我（队友什么都没说），不过最后我用1个半小时成功把代码修改好了，结果也很不错。然后就是和队友修修改改论文，八点交MD5回去休息，任务完成。

            开题晚上和第一天晚上都是正常作息。第二天晚上熬夜干题，我没有纯熬，12点去找了个桌子椅子睡了3个小时，3点又爬起来写到6点开始恶心又去休息了30分钟起来继续写，徐同学晚上过来我那屋躺地上睡了3小时，周同学断断续续休息了3多小时左右。

            当时我们在怎么解第二题和第三题是有分歧的。开题晚上时候第二题我认为定价策略应该算完全成本加成的系数，徐同学认为应该直接用定价，我俩谁也说服不了谁所以决定先去做题，我去做数据处理和第一题，徐同学去做第二题，但是第二天上午他发现他处理后的数据并不和他想的那样，定价和销量负相关，然后俩个人又想了想他觉得是应该用成本加成定价系数，所以我去重做第二问定价策略的部分，他去学ARIMA来做销量预测。我在第二天下午埋头磕第二题的时候徐同学尝试解第三题，我看他写了一个最优化的式子，但是我俩最优化成绩都稀烂（绩点4/5，最优化考70）所以当时讨论怎么解的时候俩人都否定了这么解，然后我继续闷头干题2，他去想新解法。晚上徐同学认为可以用背包加贪心策略去做，但是我感觉这么做一是方法简单难拿高分，而是没法用到前面的结论，所以否定这种做法，我的想法是把之前解第二问的思想加到第三问中，然后现学，磕最优化解法。但是徐同学觉得无法用python实现（我们三个都不擅长matlab，python之前没做过），然后我尝试说服他，然后他觉得我的想法并不合理，也难实现，剩下的时间来不及，提了好几个不合理的点（里面确实有一个地方相当不合理，关于如何预测第三问菜品定价的，后面我会细说）一直不支持我这么做。之后又尝试说服他了好一会儿，他说你就按你的想法做，我就开始干。

            下面写下关于那个相当不合理的点的处理。第二天（最后一天）下午四点半我突然发现，第三问的菜品定价系数预测是直接使用第二问用大类的销量和定价系数数据进行训练的神经网络。也就是说我的目标是通过对神经网络输入小类的销量，对具体小类的定价系数进行预测，但是我这里用来预测的神经网络是用大类的数据训练的，而大类的两个值：销量和定价系数，前者是大类总销量，后者是大类定价系数加权平均值。所以我此时建立的实际映射是：通过大类中的大类总销量值约为目前单个蔬菜的销量值对当前单个蔬菜的定价系数进行预测（也就是拿大类全天买的销量极小的某些位置对单品类定价系数进行预测）。这个映射就相当有问题，所以我还是坚持要改，挑选数据，重新训练小品类模型把大品类模型换掉。

    四、具体解法

    4.1数据处理

            首先引入处理excel文件的pandas包和处理矩阵用的numpy包。

            读入附件1，这个文件就是小品id、小品类名称、所在的大品类id、大品类名称的对应，相当适合做成hash表，我用python dict()实现。并顺序存储了两个名称、id在list()中，以后遍历就按照这个遍历和确定对立关系（小品类251种，大品类6种）。

            读入附件2，文件太大肯定不能每次都来读取他，需要压缩转存。浏览四个问题我们可以知道，题目中提问的时间最小单位为“天”，而附件二每一笔交易是不同时刻的，即“X年X月X日X时X分”，交易的时和分我们是不需要考虑的，所以我们可以读入附件2后，首先遍历确定其天数有多少（函数名：function_time_count(self)，result:1085天）。然后初始化两个空矩阵，一个251*1085用于存放每个小品类的每天销量，一个6*1085用于存放每个大品类的每天销量。

    1. def read_datas(path_data1, path_data2, path_data3, path_data2_sort):
    2.     data2 = pd.read_excel(path_data2)
    3.     data1 = pd.read_excel(path_data1)
    4.     data3 = pd.read_excel(path_data3)
    5.     data22 = pd.read_excel(path_data2_sort)
    6.     return data1, data2, data3, data22
    7.  
    8.  
    9. class Data():
    10.     def __init__(self, path_data1=r"D:\win10浏览器\CUMCM2023Problems\C题\附件1.xlsx",
    11.                  path_data2=r"D:\win10浏览器\CUMCM2023Problems\C题\附件2_做第一题简化版.xlsx",
    12.                  path_data3=r"D:\win10浏览器\CUMCM2023Problems\C题\附件3_去除改版.xlsx",
    13.                  path_data2_sort=r"D:\win10浏览器\CUMCM2023Problems\C题\附件2_打折与否排序.xlsx",
    14.                  read_ori_plus=False
    15.                  ):
    16.         self.data1, self.data2, self.data3, self.data22 = read_datas(path_data1, path_data2, path_data3,
    17.                                                                      path_data2_sort)
    18.         self.types = dict()
    19.         self.data1_values = self.data1.values
    20.         self.data2_values = self.data2.values
    21.         self.data3_values = self.data3.values
    22.  
    23.         self.data22_values = self.data22.values
    24.         self.data2_values_fullSort = self.data22_values
    25.         self.read_ori_plus = read_ori_plus
    26.  
    27.     def function_integration(self):
    28.         self.function1()
    29.         self.function2()
    30.         # self.record_one_food_function()
    31.  
    32.     # dict[单位编号]=序列号
    33.     def function1(self):
    34.         rows_type = dict()  # dict[单位编号]=序列号
    35.         i = -1
    36.         type_names = list()
    37.         for row in self.data1_values:
    38.             i += 1
    39.             rows_type[int(row[0])] = i
    40.             type_names.append(row[1])
    41.  
    42.         self.type_names = type_names
    43.         self.rows_type = rows_type
    44.  
    45.         # 保存文件
    46.         tf = open(get_path("names_type.json"), "w")
    47.         json.dump(type_names, tf)
    48.         tf.close()
    49.  
    50.         tf = open(get_path("findIndexFromID.json"), "w")
    51.         json.dump(rows_type, tf)
    52.         tf.close()
    53.  
    54.     # 共6种大类 dict[分类]=list(编码)
    55.     # dict[小id]=爸爸id
    56.     def function2(self):
    57.         rows_type_big = dict()  # dict[分类]=list(编码)
    58.         typeBig_names = list()
    59.         find_bigtype_fromsmallID = dict()  # dict[小id]=爸爸id
    60.         typeBig_ids = list()  # 大类的id
    61.  
    62.         before_typeBig = 123
    63.         for row in self.data1_values:
    64.             big_id = int(row[2])
    65.             big_name = row[3]
    66.             id = int(row[0])
    67.             if before_typeBig != big_id:
    68.                 before_typeBig = big_id
    69.                 typeBig_names.append(big_name)
    70.                 typeBig_ids.append(big_id)
    71.                 rows_type_big[big_id] = list()
    72.             rows_type_big[big_id].append(id)
    73.             find_bigtype_fromsmallID[id] = big_id
    74.  
    75.         self.typeBig_ids = typeBig_ids
    76.         self.rows_type_big = rows_type_big  # dict[分类]=list(编码)
    77.         self.typeBig_names = typeBig_names
    78.         self.find_bigtype_fromsmallID = find_bigtype_fromsmallID
    79.  
    80.         # 保存文件
    81.         tf = open(get_path("findIndexFromBigID.json"), "w")
    82.         json.dump(rows_type_big, tf)
    83.         tf.close()
    84.         tf = open(get_path("names_bigType.json"), "w")
    85.         json.dump(typeBig_names, tf)
    86.         tf.close()
    87.         tf = open(get_path("listBig_ids.json"), "w")
    88.         json.dump(typeBig_ids, tf)
    89.         tf.close()
    90.         tf = open(get_path("find_bigtype_fromsmallID.json"), "w")
    91.         json.dump(find_bigtype_fromsmallID, tf)
    92.         tf.close()
    93.  
    94.     # 用于求得单个品种食物每天的销量
    95.     # 用于求得单个大品类每天的销量
    96.     def record_one_food_function(self):
    97.         time = '0'
    98.         record_one_type = np.zeros(shape=(251, 1085))
    99.         record_one_bigtype = np.zeros(shape=(6, 1085))
    100.         i = -1
    101.         record_time = list()
    102.         for current in self.data2_values:
    103.             current_time = current[0]
    104.             current_id = int(current[1])
    105.             current_weight = float(current[2])
    106.             if time != current_time:
    107.                 time = current_time
    108.                 i += 1  # index
    109.                 record_time.append(time)
    110.  
    111.             food_number = self.rows_type[current_id]
    112.             # print(food_number,self.type_names[food_number],current_id)
    113.             record_one_type[food_number][i] += current_weight
    114.  
    115.             big_id = self.find_bigtype_fromsmallID[current_id]
    116.             index_bigId = self.typeBig_ids.index(big_id)
    117.             # print(current_id,index_bigId,big_id,self.typeBig_names[index_bigId])
    118.             record_one_bigtype[index_bigId][i] += current_weight
    119.  
    120.         self.record_one_type = record_one_type
    121.         self.record_one_bigtype = record_one_bigtype
    122.         self.record_time = record_time
    123.  
    124.     # 用于求得单个品种的每天销量
    125.     # result=1085
    126.     def function_time_count(self):
    127.         # 看有几天
    128.         sum_time = 0
    129.         current = '0'
    130.         for row in self.data2_values:
    131.             if current != row[0]:
    132.                 current = row[0]
    133.                 sum_time += 1
    134.         print(sum_time)
    135.         # 按行遍历
    136.         current_time = '0'
    137.         for row in self.data2_values:
    138.             if current_time != row[0]:
    139.                 current_time = row[0]
    140.  
    141.     # 计算每日物品的进货价
    142.     def function_deal_data3(self):
    143.         values3 = self.data3_values
    144.         times3_list = list()
    145.         record_cost = np.zeros(shape=(251, 1085))
    146.  
    147.         time = '0'
    148.         i = -1
    149.         for current in values3:
    150.             current_time = current[0]
    151.             current_id = int(current[1])
    152.             current_cost = float(current[2])
    153.             if time != current_time:
    154.                 time = current_time
    155.                 i += 1
    156.                 times3_list.append(time)
    157.  
    158.             food_number = self.rows_type[current_id]
    159.             record_cost[food_number][i] = current_cost
    160.         # save
    161.         np.save(r'D:\win10浏览器\CUMCM2023Problems\C题\record_cost.npy', record_cost)
    162.         np.save(r"D:\win10浏览器\CUMCM2023Problems\C题\time3_list.npy", np.array(times3_list))
    163.  
    164.     # 时间Int化
    165.     def deal_with_time_save(self):
    166.         self.time_Int = time_Int(self.time3_list)
    167.         np.save(r'D:\win10浏览器\CUMCM2023Problems\C题\time_Int.npy', self.time_Int)
    168.  
    169.     # 将每天的每个品类和每个单品的销量进行记录
    170.     def save_array(self):
    171.         np.save(r'D:\win10浏览器\CUMCM2023Problems\C题\record_one_type.npy', self.record_one_type)
    172.         np.save(r'D:\win10浏览器\CUMCM2023Problems\C题\record_one_bigtype.npy', self.record_one_bigtype)
    173.  
    174.     def load_array(self):
    175.         # record_one_food_function
    176.         self.record_one_type = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\record_one_type.npy')
    177.         self.record_one_bigtype = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\record_one_bigtype.npy')
    178.  
    179.         self.record_cost = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\record_cost.npy')
    180.         self.time3_list = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\time3_list.npy', allow_pickle=True)
    181.  
    182.         if self.read_ori_plus == True:
    183.             self.plus_type = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\plus_type.npy')
    184.             self.weight_type_noDiscount = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\weight_type_noDiscount.npy')
    185.             self.weight_type_Discount = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\weight_type_Discount.npy')
    186.             self.plus_lower_type = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\plus_lower_type.npy')
    187.  
    188.         self.plus_Bigtype = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\plus_Bigtype.npy')
    189.         self.weight_Bigtype_noDiscount = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\weight_Bigtype_noDiscount.npy')
    190.         self.weight_Bigtype_Discount = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\weight_Bigtype_Discount.npy')
    191.         self.plus_lower_Bigtype = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\plus_lower_Bigtype.npy')
    192.  
    193.         self.time_Int = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\time_Int.npy')
    194.         self.time_Int_list = list(self.time_Int)
    195.         # function1
    196.         tf = open(get_path("names_type.json"), "r")
    197.         new_dict = json.load(tf)  # list
    198.         self.names_type = new_dict
    199.         tf.close()
    200.         self.type_names = new_dict
    201.  
    202.         tf = open(get_path("findIndexFromID.json"), "r")
    203.         new_dict = json.load(tf)
    204.         self.findIndexFromID = new_dict  # dict
    205.         tf.close()
    206.         self.rows_type = new_dict
    207.  
    208.         # function2
    209.         tf = open(get_path("findIndexFromBigID.json"), "r")
    210.         new_dict = json.load(tf)  #
    211.         self.findIndexFromBigID = new_dict
    212.         tf.close()
    213.         self.rows_type_big = new_dict  # dict[分类]=list(编码)
    214.  
    215.         tf = open(get_path("names_bigType.json"), "r")
    216.         new_dict = json.load(tf)
    217.         self.names_bigType = new_dict  # dict
    218.         tf.close()
    219.         self.typeBig_names = new_dict
    220.  
    221.         tf = open(get_path("listBig_ids.json"), "r")
    222.         new_dict = json.load(tf)  # list
    223.         self.listBig_ids = new_dict
    224.         tf.close()
    225.         self.typeBig_ids = new_dict
    226.  
    227.         tf = open(get_path("find_bigtype_fromsmallID.json"), "r")
    228.         new_dict = json.load(tf)
    229.         self.find_bigtype_fromsmallID = new_dict  # dict
    230.         tf.close()
    231.         self.find_bigtype_fromsmallID = new_dict  # 没变
    232.  
    233.     # 计算加成值
    234.     def addition(self):
    235.         pass
    236.  
    237.     def save_array2(self):
    238.         pass

            其中的核心代码是（获得251种小类和6大类1085天每日销量）：

    1.     # 用于求得单个品种食物每天的销量
    2.     # 用于求得单个大品类每天的销量
    3.     def record_one_food_function(self):
    4.         time = '0'
    5.         record_one_type = np.zeros(shape=(251, 1085))
    6.         record_one_bigtype = np.zeros(shape=(6, 1085))
    7.         i = -1
    8.         record_time = list()
    9.         for current in self.data2_values:
    10.             current_time = current[0]
    11.             current_id = int(current[1])
    12.             current_weight = float(current[2])
    13.             if time != current_time:
    14.                 time = current_time
    15.                 i += 1  # index
    16.                 record_time.append(time)
    17.  
    18.             food_number = self.rows_type[current_id]
    19.             # print(food_number,self.type_names[food_number],current_id)
    20.             record_one_type[food_number][i] += current_weight
    21.  
    22.             big_id = self.find_bigtype_fromsmallID[current_id]
    23.             index_bigId = self.typeBig_ids.index(big_id)
    24.             # print(current_id,index_bigId,big_id,self.typeBig_names[index_bigId])
    25.             record_one_bigtype[index_bigId][i] += current_weight
    26.  
    27.         self.record_one_type = record_one_type
    28.         self.record_one_bigtype = record_one_bigtype
    29.         self.record_time = record_time

    4.2第一问求解

            引入做二维图的matplotlib.pyplot。

            对于大类，一共就6个，直接使用np.corrcoef()计算他们之间的相关系数（皮尔逊相关系数），然后用plt把相关系数们当作热力图输出出来就成，效果如下：

    0c1fd5d73aef4da4b30cff63b14bcafc.png

            核心代码（热力图）： 

    1. # 热力图
    2. def heat_map_type(d):
    3.     np.set_printoptions(threshold=np.inf)
    4.  
    5.     # print(d.record_one_type[:,0])
    6.     # print(d.record_one_type[:,-1])
    7.     # print(d.record_one_bigtype[0,:])
    8.     corr = pd.DataFrame(d.record_one_type)
    9.     corr = corr.corr()
    10.     # print(corr)
    11.  
    12.     fig = plt.figure()
    13.     ax = fig.add_subplot(1, 1, 1)
    14.     p1 = sns.heatmap(corr, xticklabels=d.type_names, yticklabels=d.type_names, cmap='viridis')
    15.     ax.set_title('Heat Map2')
    16.  
    17.     s1 = p1.get_figure()
    18.     plt.show()
    19.     s1.savefig(r'D:\win10浏览器\CUMCM2023Problems\C题\HeatMap_small.jpg', dpi=300, bbox_inches='tight')
    20.     print('ok')

            对于小类，251种啊，画出个热力图一个小格快成像素点了不能这么弄，所以用聚类，我用的k-means，为了获取相关度尽量大的把k设为80，将与其他蔬菜相关度低的蔬菜单独聚成一簇，那些簇内元素不为1的就是我实际需要的相关度较高的小类蔬菜。聚类的结果是有相当多一部分的蔬菜是聚到了一类（超过20个菜品聚到了一个簇），剩余的都是一个簇内不超过15个，以四个以下居多，举两个例子：

    44a001c771f4458e8915d3ce2953a02c.png

            对于那个相当多元素的簇，我的解释是那些是蔬菜可能可以种大棚里，或者居民在各个季节需求都比较平均，没啥淡旺季，所以进一簇了（就比如小米辣椒）。

            核心代码（小品类寻找有联系的蔬菜）：

    1. def KMS_type_show(d, k=80):
    2.     kms = KMeans(n_clusters=k, random_state=0).fit(d.record_one_type)
    3.     labels = kms.labels_
    4.     # record_cluster = labels.sorted()
    5.     # 首先剔除常用少量蔬菜
    6.     count = np.bincount(labels)
    7.     biggest_count = count.max()
    8.     biggest_count_number = np.where(biggest_count == count)[0][0]
    9.     # 其次剔除单个元素
    10.     cluster_need_index = list()
    11.     cluster_need_count = list()
    12.     index = -1
    13.     for i in count:
    14.         index += 1
    15.         if i <= 1 or i == biggest_count:
    16.             pass
    17.         else:
    18.             cluster_need_index.append(index)  # 簇名
    19.             cluster_need_count.append(i)  # 出现次数
    20.     print('最大簇的簇的簇标记为{},其内商品有:'.format(biggest_count_number))
    21.     # goods序号
    22.     cluster_Biggest = np.where(labels == biggest_count_number)[0]
    23.     cluster_set = np.zeros(shape=(biggest_count, 1085))
    24.     names_cs = list()
    25.     for i in range(len(cluster_Biggest)):
    26.         print(d.type_names[cluster_Biggest[i]])
    27.         cluster_set[i, :] = d.record_one_type[cluster_Biggest[i], :]
    28.         names_cs.append(d.type_names[cluster_Biggest[i]])
    29.     # show
    30.     x = np.arange(1085)
    31.     x = np.expand_dims(x, 0).repeat(biggest_count, axis=0)
    32.     plt.plot(x.T, cluster_set.T)
    33.     plt.legend(names_cs, loc="best")
    34.     plt.show()
    35.  
    36.     # 处理有效簇
    37.     print('以下为有效簇的信息：')
    38.     for i in range(len(cluster_need_index)):
    39.         print('簇标记为{},簇内有商品{}个,分别为：'.format(cluster_need_index[i], cluster_need_count[i]))
    40.         cluster_current = np.where(labels == cluster_need_index[i])[0]
    41.         names_cs = list()
    42.         cluster_set = np.zeros(shape=(cluster_need_count[i], 1085))
    43.  
    44.         index = -1
    45.         for j in cluster_current:
    46.             index += 1
    47.             print('{}'.format(d.type_names[j], end=' '))
    48.             names_cs.append(d.type_names[j])
    49.             cluster_set[index, :] = d.record_one_type[j, :]
    50.             # show
    51.         x = np.arange(1085)
    52.         x = np.expand_dims(x, 0).repeat(cluster_need_count[i], axis=0)
    53.         plt.plot(x.T, cluster_set.T)
    54.         plt.legend(names_cs, loc="best")
    55.         plt.show()

    4.3第二问求解

            首先要做出大品类销量和成本加成定价的关系，我的做法是找单个大品类每日的销量和成本加成系数的关系，而题目中说明了"从需求侧来看，蔬菜类 商品的销售量与时间往往存在一定的关联关系",所以认为应该考虑大品类每日的销量是受时间的影响的，因此我对每一种大类都使用聚类聚成两类，将其分为淡旺季，然后再进一步对每一季处理。

            核心代码（分两季）：

    1. def show_bigtype(d):
    2.     x = np.arange(1085)
    3.     x = np.expand_dims(x, 0).repeat(6, axis=0)
    4.     plt.plot(x.T, d.record_one_bigtype.T)
    5.     plt.legend(d.typeBig_names, loc="best")
    6.     plt.show()
    7. def mounth_analysis():
    8.     mounth_record = np.zeros(shape=(6, 12))
    9.     two_season = np.zeros(shape=(6, 12))
    10.  
    11.     record = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\record_one_bigtype.npy')
    12.     day_record = np.load(r'D:\win10浏览器\CUMCM2023Problems\C题\time3_list.npy', allow_pickle=True)
    13.  
    14.     for i in range(6):
    15.         current_month = 7
    16.         index = 0
    17.         for j in range(len(day_record)):
    18.             if current_month != day_record[j].month:
    19.                 current_month = day_record[j].month
    20.             mounth_record[i][current_month - 1] += record[i][j]
    21.     for i in range(6):
    22.         kms = KMeans(n_clusters=2, random_state=0).fit(mounth_record[i, :].reshape(-1, 1))
    23.         two_season[i, :] = kms.labels_
    24.     np.save(r"D:\win10浏览器\CUMCM2023Problems\C题\two_season.npy", two_season)
    25.     print(two_season)

            分为两季后，对于每一季的销量和定价系数分别处理、显示、分析，以下举两个图片：

    82ef4bdf49f244feb20d0c63f9a580ab.png

            可以看出，销量受加成系数影响较低，不呈现负相关，商超在某种菜品销量高时更倾向于提高定价以获取更多利润。

            然后要给出各蔬菜品类未来一周(2023 年 7 月 1-7 日)的日补货总量和定价策略，我们首先用季节性ARIMA时间序列分析对这一周的每日销量进行预测，然后将销量作为这几天的日补货总量。再用季节销量做属性，加成系数做标签，用神经网络对他俩做拟合，以后通过销量预测定价用。部分拟合结果：

    400bcbc0c4ca459f88bceb3bb2f401f2.png49919db24a0f4d179564652802e36066.png

            之后将ARIMA预测到的销量输入刚才训练好的bp网络中，确定加成定价系数。

            核心代码（训练神经网路并展示模型）：

    1. sklearn.neural_network import MLPRegressor
    2. from sklearn.model_selection import train_test_split
    3. import joblib
    4. def Regression_show2(save_sn):
    5.     # 防止刷新模型
    6.     if save_sn != None:
    7.         return
    8.     model = list(list(list() for j in range(2)) for i in range(6))
    9.     for foodBig in range(6):
    10.         for reason in range(2):
    11.             # x是加成定价，y是销量
    12.             X = save_sn[foodBig][reason][0].reshape(-1)
    13.             y = save_sn[foodBig][reason][1].reshape(-1, 1)
    14.             X, y = y, X
    15.  
    16.             X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.38, random_state=42)
    17.             mlp_regressor = MLPRegressor(hidden_layer_sizes=(400), max_iter=5000, random_state=42)
    18.             mlp_regressor.fit(X_train, y_train)
    19.             y_pred = mlp_regressor.predict(X_test)
    20.             plt.scatter(X_test, y_test, c='b', label='Data')
    21.             plt.scatter(X_test, y_pred, c='r', label='Predictions')
    22.             plt.xlabel('X')
    23.             plt.ylabel('y')
    24.             plt.legend()
    25.             plt.show()
    26.             # save model
    27.             joblib.dump(mlp_regressor,
    28.                         r'D:\win10浏览器\CUMCM2023Problems\C题\model2\model{}_{}.pkl'.format(foodBig, reason))
    29. #提选出时间

    4.4第三问求解

            首先写上约束公式，出优化公式：

    b50d588be50540eab362732c9142c32a.png

            公式解释：

    1fa6b31c661b4da0b950533b56c11795.png   

            用scikit-opt(国人写的包)里的GA(遗传算法)和sklearn.neural_network里的MLPRegressor。先筛选出6月24日到6月30日的可以进的蔬菜，然后用第二题的相同的方法每一个小类训练训练60个输入销量输出加成系数的神经网络，之后引入sko的GA包把优化公式写进去迭代求解就行，附上官方API:

    文档 (scikit-opt.github.io)

            迭代四百次的结果（论文里好像忘记加这个图了呜呜）：

    641fc44f404c4281a4a30e1b1007ef7d.png

            可以看出最后最优值（纯利润）在1150RMB收敛了。

            核心代码（训练小类网络和用GA算法解优化问题）：

    1. def Regression_show():
    2.     plus = np.load(get_path('plus_2_624707.npy'))
    3.     wei  = np.load(get_path('wei_2_624707.npy'))
    4.     #food_sale可用
    5.  
    6.     for index in range(max_choice):
    7.         #id = food_sale[index]
    8.         #type_index = findIndexFromID[str(id)]
    9.         # x是加成定价，y是销量
    10.         X = plus[index][:].reshape(-1)
    11.         y = wei[index][:].reshape(-1,1)
    12.         X, y = y, X
    13.  
    14.         X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.38, random_state=42)
    15.         mlp_regressor = MLPRegressor(hidden_layer_sizes=(400), max_iter=5000, random_state=42)
    16.         mlp_regressor.fit(X_train, y_train)
    17.         if index<10:
    18.             y_pred = mlp_regressor.predict(X_test)
    19.             plt.scatter(X_test, y_test, c='b', label='Data')
    20.             plt.scatter(X_test, y_pred, c='r', label='Predictions')
    21.             plt.xlabel('X')
    22.             plt.ylabel('y')
    23.             plt.legend()
    24.             plt.show()
    25.         # save model
    26.         joblib.dump(mlp_regressor,r'D:\win10浏览器\CUMCM2023Problems\C题\model2\model{}.pkl'.format(index))
    27. def read_modelAnd_predict2(predict_list, bigFoodIndex=0, reason=0):
    28.     # load model
    29.     # print('predict_list',predict_list)
    30.     rfc2 = joblib.load(r"D:\win10浏览器\CUMCM2023Problems\C题\model\model{}_{}.pkl".format(bigFoodIndex, reason))
    31.     result = rfc2.predict(np.array(predict_list).reshape(-1, 1))
    32.     # print('result',result)
    33.     return result
    34. def read_modelAnd_predict(predict_list,Index=0):
    35.     # load model
    36.     # print('predict_list',predict_list)
    37.     rfc2 = joblib.load(r"D:\win10浏览器\CUMCM2023Problems\C题\model2\model{}.pkl".format(Index))
    38.     result = rfc2.predict(np.array(predict_list).reshape(-1, 1))
    39.     # print('result',result)
    40.     return result
    1. #优化目标
    2. def fun_question(x):
    3.     #loss_npList = list(loss_npList)
    4.     es = x[:max_choice:]#选取61个值     进货量
    5.     qs = x[max_choice::]#选取61个值     是否选择
    6.     w=0.0
    7.     for i in range(max_choice):#遍历每一个可取值
    8.         #确定id,选择父类
    9.         id = food_sale[i]
    10.         #print(data4_index_list,'\n',id)
    11.  
    12.         loss_index = data4_index_list.index(id)
    13.         loss = loss_npList[loss_index]*0.01*0.75#增加
    14.  
    15.         e_loss = es[i]*(1-loss)
    16.         #print('es[i]={},loss={}'.format(es[i],loss))
    17.         #寻找物品当前的批发价
    18.         cost = cost_71[i]
    19.         if qs[i]>0.5:
    20.             one = read_modelAnd_predict(e_loss,i)*e_loss*qs[i]*cost-cost*es[i]
    21.             w += one
    22.             #print('one',one)
    23.  
    24.     return -1*w[0]
    25. from sko.GA import GA
    26. #GA计算
    27. def optimize_function():
    28.     #进货量
    29.     #是否选择
    30.     constraint_ueq = (
    31.         lambda x: sum(x[max_choice::]) - 33,
    32.         lambda x: 27 - sum(x[max_choice::]),
    33.     )
    34.     arima_one = read_excel_71_predict()
    35.     std_one   = np.load(get_path('std_canSale.npy'))
    36.     std_two   = list( 2*i for i in std_one)
    37.     lb_0 = [0]*(max_choice)
    38.     ub_1 = [1]*(max_choice)
    39.     #只看销量
    40.     # 最小值
    41.     min_put = [2.5]*max_choice
    42.     #arima +- 2*std
    43.     upper = np.add(arima_one,std_one)
    44.     lower = np.subtract(arima_one,std_one)
    45.  
    46.     lb = find_lb(min_put,lower)
    47.     ub = finb_up(min_put,upper)
    48.     lb = list(lb)+list(lb_0)
    49.     ub = list(ub)+list(ub_1)
    50.  
    51.     precision = [1e-7 for i in range(max_choice)]+[1 for i in range(max_choice)]
    52.     ga = GA(func=fun_question,n_dim=max_choice*2,
    53.             size_pop=50,max_iter=300,prob_mut=0.001,
    54.             lb = lb,ub=ub,
    55.             constraint_ueq=constraint_ueq,
    56.             precision=precision
    57.             )
    58.     best_x, best_y = ga.run()
    59.     print('best_x:', best_x, '\n', 'best_y:', best_y)
    60.     np.save(get_path(r'answer3\best_x.npy'),best_x)
    61.     np.save(get_path(r'answer3\best_y.npy'), best_y)
    62.     Y_history = pd.DataFrame(ga.all_history_Y)
    63.     fig, ax = plt.subplots(2, 1)
    64.     ax[0].plot(Y_history.index, Y_history.values, '.', color='red')
    65.     Y_history.min(axis=1).cummin().plot(kind='line')
    66.     plt.show()

    五、更多代码

            如果你想查看完整代码，可以访问http://t.csdn.cn/Ux4Pz，仅供学习交流使用，严禁转载商用。（一千行未维护的代码，直接用来水大作业不太行）

     

     

  • 相关阅读:
    高分辨空间代谢组学测试的样品要求以及常见问题
    SCI常用经典词和常用句型（一）
    Flask 学习-34.restful-full 请求参数自定义参数校验类型 (reqparse.RequestParser() )
    Oracle/PLSQL: Cast Function
    Java基于微信小程序的一起考研学习平台
    Python之ini配置文件详解
    TCP协议
    Azure Machine Learning - Azure AI 搜索中的矢量搜索
    pdd.order.list.get拼多多店铺订单列表查询接口（拼多多店铺订单详情接口，订单明文接口，订单解密接口，订单插旗接口，订单备注接口）代码对接教程
    Hadoop高可用环境搭建-HDFSNameNode高可用搭建、Yarn高可用搭建
  • 原文地址：https://blog.csdn.net/Zhonexixi/article/details/132857862