# 题干理解解读

(从题干出发的想法)
存在变量,keyword:历史销售,需求情况,蔬菜品种,进货价格,“成本加成定价”,商品品相,sale=f(t)sale=f(t)(需求侧出发)

que1:菜类商品不同品类或不同单品之间可能存在一定的关联关系,请分析蔬菜各品类及单品销售量的分布规律及相互关系。

que2:考虑商超以品类为单位做补货计划,请分析各蔬菜品类的销售总量与成本加成定价的关系,并给出各蔬菜品类未来一周(2023 年 7 月 1-7 日)的日补货总量和定价策略,使得商超收益最大

que3:因蔬菜类商品的销售空间有限,商超希望进一步制定单品的补货计划,要求可售单品总数控制在 27-33 个,且各单品订购量满足最小陈列量 2.5 千克的要求。根据 2023年 6 月 24-30 日的可售品种,给出 7 月 1 日的单品补货量和定价策略,在尽量满足市场对各品类蔬菜商品需求的前提下,使得商超收益最大

# que1

# 1.1自己的理解

看题目的意思应该是要去找到某个类似与矩阵的东西,里面填充各种单品之间的相关系数去得到不同单品之间的相关关系,然后写出不同单品之的分布规律

  • 分布规律:按照商品编号进行分类统计,对销售进行求和,绘制出分布型绘图描述销售量关系

  • 单品之间的相关关系:这里我考虑是用皮尔逊相关系数矩阵进行解决问题,然后需要去权衡哪些数值对他们的相关系数产生影响呢?

# 1.2相关论文

  • 行文特点:

    • 摘要:首先进行XX算法,算出了XX,得到了XX,然后进行XX(图形绘制数据可视化......)操作,结论
    • 正文:
      首先进行数据预处理,剔除无效数据(如退货导致的负销量数据);然后计算描述性统计量(偏度系数、峰度系数、均值、中位数、标准差等),从数字特征层面分析分布规律;接着通过绘制柱状图、折线图等可视化图表,直观展示销售数据的分布趋势;之后运用相关系数分析(如 Spearman 相关系数)探究品类及单品间的关联关系;最后采用聚类算法(如 K-means++)对单品进行分类,进一步分析不同类别单品的销售特征及关联。
  • 数据处理相关:
    考量分布规律可选取如下数据偏度系数峰度系数,最小值,最大值,均值,中位数,标准差等常见统计学指标

    同时分析规律可以从多角度出发,总销量/日销量等等

  • 算法使用:

    • Spearman相关系数

      正如我之前的想法还挺正确,然后我找到了如下对比

      特征 Pearson相关系数矩阵 Spearman相关系数矩阵
      相关性类型 线性相关 单调相关(线性或非线性)
      数据要求 要求数据正态分布 无分布要求,非参数方法
      计算基础 原始数据值 数据的秩次
      异常值敏感性 高度敏感 相对稳健
      适用场景 线性关系明显的连续变量 非线性但单调的关系或有序变量

      同时他们两个的公式:

      rXY=i=1n(XiXˉ)(YiYˉ)i=1n(XiXˉ)2i=1n(YiYˉ)2r_{XY} = \frac{\sum_{i=1}^n (X_i - \bar{X})(Y_i - \bar{Y})}{\sqrt{\sum_{i=1}^n (X_i - \bar{X})^2} \sqrt{\sum_{i=1}^n (Y_i - \bar{Y})^2}}

      ρ=16i=1ndi2n(n21)\rho = 1 - \frac{6 \sum_{i=1}^n d_i^2}{n(n^2 - 1)}

      实现方法(多组数据得到spearman矩阵):

      import pandas as pd
      from scipy.stats import spearmanr
      
      data = pd.DataFrame({
          'A': [106, 86, 100, 101, 99, 103, 97, 113, 112, 110],
          'B': [7, 0, 27, 50, 28, 29, 20, 12, 6, 17],
          'C': [10, 20, 30, 40, 50, 60, 70, 80, 90, 100]
      })
      
      corr_matrix, p_matrix = spearmanr(data)
      
      print("Spearman 相关系数矩阵:")
      print(pd.DataFrame(corr_matrix, index=data.columns, columns=data.columns))
    • Kmeans++聚类算法

      Kmeans++算法相较于传统的Kmeans聚类算法优化了聚类中心的选取

      实现方法:

      1. 确定好分成几类
      2. 确定好影响变量
      3. 丢进spss即可自动聚类生成表单
    • 自相关函数ACF(探究分布规律)

      通过计算不同滞后阶段的自相关系数,绘制 ACF 图,判定销售量时间序列的周期性。例如,发现总销量存在年周期(滞后 12 个月处相关性显著)、周周期(滞后 7 天处相关性显著)和日周期(滞后 13 小时处相关性显著)。

      自相关系数计算公式:

      ρ^(k)=1Tkt=k+1T(ytyˉ)(ytkyˉ)1Tt=1T(ytyˉ)2\hat{\rho}(k) = \frac{\frac{1}{T-k} \sum_{t=k+1}^T (y_t - \bar{y})(y_{t-k} - \bar{y})}{\frac{1}{T} \sum_{t=1}^T (y_t - \bar{y})^2}

      实现方法

      import matplotlib.pyplot as plt
      from statsmodels.graphics.tsaplots import plot_acf
      
      # 生成数据(示例:AR(1)过程)
      import numpy as np
      np.random.seed(42)
      y = np.random.normal(size=100)
      for t in range(1, 100):
          y[t] = 0.7 * y[t-1] + np.random.normal(0, 0.1)
      
      # 绘制ACF
      plot_acf(y, lags=20, alpha=0.05)  # 95%置信区间
      plt.title("样本自相关函数(ACF)")
      plt.show()
    • 时间序列分解模型(探究分布规律):

      首先明确分解是如何分解的Yt=Tt×St×Ct×EtY_t=T_t\times S_t\times C_t\times E_t,其中YtY_t为观测值TtT_t为趋势成分STS_T为季节性成分CtC_t为周期性成分EtE_t为随机噪声

      作用:乘法时间序列分解的作用主要是挖掘蔬菜销售量在时间维度上的分布规律。通过将销售量时间序列拆解为趋势、季节性、周期性和随机噪声等成分,能够清晰呈现不同时间颗粒(如年、周、日)下销量的波动特征和变化趋势

      实现方法:

      import pandas as pd
      from statsmodels.tsa.seasonal import seasonal_decompose
      import matplotlib.pyplot as plt
      
      # 准备时间序列数据(索引为datetime类型,值为销量)
      df = pd.read_csv('销量数据.csv', parse_dates=['销售日期'], index_col='销售日期')
      series = df['销量(千克)']
      
      # 乘法时间序列分解(以周周期为例,period=7)
      result = seasonal_decompose(series, model='multiplicative', period=7)
      
      # 提取各成分
      trend = result.trend  # 趋势成分
      seasonal = result.seasonal  # 季节性成分
      residual = result.resid  # 残差(随机噪声)
      
      # 可视化分解结果
      fig = result.plot()
      plt.show()

      实现效果图

    • FP-Growth关联挖掘算法

      这个算法常常用来从大量交易数据中发现频繁项集

      ex:假设有如下交易

      交易ID 购买物品
      T1 牛奶,面包,尿布
      T2 可乐,面包,尿布
      T3 牛奶,尿布
      T4 牛奶,可乐,尿布
      T5 可乐,面包

      设置最小支持度阈值为2
      FP-Growth会构建FP-tree,发现频繁项集,比如 {尿布}、{牛奶, 尿布}、{面包, 尿布} 等。

      实现方法:

      from mlxtend.preprocessing import TransactionEncoder
      from mlxtend.frequent_patterns import fpgrowth, association_rules
      
      # 交易数据
      transactions = [
          ['milk', 'bread', 'diaper'],
          ['cola', 'bread', 'diaper'],
          ['milk', 'diaper'],
          ['milk', 'cola', 'diaper'],
          ['cola', 'bread']
      ]
      
      # 1. 编码成布尔矩阵
      te = TransactionEncoder()
      te_ary = te.fit(transactions).transform(transactions)
      import pandas as pd
      df = pd.DataFrame(te_ary, columns=te.columns_)
      
      # 2. 挖掘频繁项集,设最小支持度为0.4
      freq_itemsets = fpgrowth(df, min_support=0.4, use_colnames=True)
      
      print("频繁项集:")
      print(freq_itemsets)
      
      # 3. 挖掘关联规则(最小置信度0.7)
      rules = association_rules(freq_itemsets, metric="confidence", min_threshold=0.7)
      
      print("\n关联规则:")
      print(rules[['antecedents', 'consequents', 'support', 'confidence', 'lift']])
      

# que2

# 2.1 自己的理解

首先题目要求求解各蔬菜品类的销售总量与成本加成定价的关系,肯定是需要用归化算法去拟合出个函数关系,其次并给出各蔬菜品类未来一周(2023 年 7 月 1-7 日)的日补货总量和定价策略,使得商超收益最大我觉得可以联系ACF值和上述求解得到的函数关系进行定义/利用归化算法解决决策问题

# 2.2 相关论文

  • 成本加成定价法

    首先要考虑成本加成定价法和各蔬菜品类的销售总量,我们得先知道什么是成本加成定价法
    成本加成定价法

  • Pearson相关系数(衡量是否线性相关)

    rXY=i=1n(XiXˉ)(YiYˉ)i=1n(XiXˉ)2i=1n(YiYˉ)2r_{XY} = \frac{\sum_{i=1}^n (X_i - \bar{X})(Y_i - \bar{Y})}{\sqrt{\sum_{i=1}^n (X_i - \bar{X})^2} \sqrt{\sum_{i=1}^n (Y_i - \bar{Y})^2}}

    • 假设有 n 个数据点(如不同时间或不同批次的销售数据),每个数据点包括:
      • XiX_i:成本加成定价
      • YiY_i:对应品类销售总量
  • Pearson相关系数检验方法(计算得到的相关系数r是否显著)=>t检验、临界值法、p值法

    | 检验方法 | t值法 | 临界值法(查表法) | p值法 |
    |----------------|-------------------------------------|---------------------------------|-------------------------------------|
    | 假设 | H0:ρ=0H_0: \rho = 0<br>H1:ρ0H_1: \rho \neq 0(双侧) | 同左 | 同左 |
    | 统计量 | t=rn21r2t = \frac{r \sqrt{n-2}}{\sqrt{1-r^2}} | 直接使用样本相关系数 rr | 同t值法,但计算 pp 值 |
    | 关键值 | 查tt分布表(df=n2df = n-2)得 tα/2t_{α/2} | 查表得 rcriticalr_{\text{critical}} | 计算 pp 值(或软件输出) |
    | 判断标准 | t>tα/2\|t\| > t_{α/2} ⇒ 拒绝H0H_0 | r>rcritical\|r\| > r_{\text{critical}} ⇒ 拒绝H0H_0 | p<αp < α ⇒ 拒绝H0H_0 |
    | 适用场景 | 需计算t统计量,适合理论推导 | 快速查表,适合手工计算 | 现代统计软件默认方法,结果更直观 |
    | 优缺点 | 严格但计算稍复杂 | 简单但需依赖临界值表 | 直接提供概率,无需查表 |

  • 线性回归拟合线性关系

    • 假定线性关系
      • cpp:成本加成定价
      • cs:蔬菜品类销售总量

    cpp=β0+β1^cscpp=\beta_{0}+\hat{\beta_{1}}cs

    • 使用最小二乘法求解系数

    β^1=i=1n(csicsˉ)(cppicppˉ)i=1n(csicsˉ)2\hat{\beta}_1 = \frac{\sum_{i=1}^n (cs_i - \bar{cs})(cpp_i - \bar{cpp})}{\sum_{i=1}^n (cs_i - \bar{cs})^2}

  • Prophet 模型

    主要用于 单变量时间序列预测(如销量、股票价格、温度等)。它 不直接用于变量关系预测(如回归分析),但可以间接变量关系
    该模型的基本形式

    y(t)=g(t)+s(t)+h(t)+ϵty(t)=g(t)+s(t)+h(t)+\epsilon_t

    g(t)g(t):趋势向(线性/非线性)

    s(t)s(t):季节性(年/月/周)

    h(t)h(t):节假日效应

    在本题中考虑到要分析各蔬菜品类的销售总量与成本加成定价的关系,于是引入外生回归项βPrice(t)\beta·Price(t),同时剔除季节性、趋势性和节假日效应,从而分离出定价对销量的净影响

    y(t)=g(t)+s(t)+h(t)+ϵt+βPrice(t)y(t)=g(t)+s(t)+h(t)+\epsilon_t+\beta·Price(t)

    (1) 数据预处理

    • 价格变量:将成本加成定价(公式 (11))作为外生变量加入模型。
    • 节假日设定:定义中国节日(含无假期的节日如元宵节)及其窗口期(如春节前后7天)。
    • 季节性配置:仅保留 周度weekly_seasonality)和 年度yearly_seasonality)周期,忽略日内波动以降低噪声。

    (2) 模型拟合
    通过 add_regressor() 引入价格变量:

    from prophet import Prophet
    
    model = Prophet(
        yearly_seasonality=True,  # 年度周期
        weekly_seasonality=True,  # 周度周期
        holidays=holidays_df      # 自定义节假日
    )
    model.add_regressor('price')  # 添加价格作为外生变量
    model.fit(df)

    (3) 成分分解
    拟合后,通过 plot_components() 可视化各成分:

    • 趋势项 g(t)g(t):反映长期需求变化(如2022年1月低谷)。
    • 季节性 s(t)s(t):周度(周末销量高)和年度(7-8月旺季)模式。
    • 节假日 h(t)h(t):国庆高峰、春节下降等。
    • 价格效应 βPrice(t)\beta \cdot \text{Price}(t):回归系数 β\beta 的符号和大小反映定价影响。

    (4) 销量预测与关系分析

    • 预测:在已知未来价格时,可预测销量(需提供未来价格数据)。
    • 关系量化:通过 β\beta 判断定价与销量的关系(如表4中花叶类 β=21.910\beta = -21.910 表示价格每升1元,销量降约21.9单位)。
  • 模拟退火算法

    1. 算法核心思想
    模拟退火算法通过模拟金属退火过程中的 温度冷却能量最小化 过程,在优化问题中实现:

    • 全局搜索能力:通过允许暂时接受劣解,避免陷入局部最优
    • 渐进收敛性:随着温度降低,逐步聚焦到高质量解区域
    • 可控随机性:通过温度参数平衡探索(exploration)与利用(exploitation)

    2. 在补货定价问题中的建模
    优化目标

    maxpiProfit=i=17[piQi(pi)(1γi)wiQi(pi)1δi]\max_{p_i} \text{Profit} = \sum_{i=1}^{7} \left[ p_i \cdot Q_i(p_i) \cdot (1-\gamma_i) - w_i \cdot \frac{Q_i(p_i)}{1-\delta_i} \right]

    其中:

    • pip_i:第ii天的定价(决策变量)
    • Qi(pi)Q_i(p_i):由Prophet模型预测的销量(价格函数)
    • wiw_i:批发价(外生输入)
    • γi\gamma_i:折扣率
    • δi\delta_i:损耗率

    3. 算法实现步骤
    步骤1:初始化

    current_p = [last_day_price] * 7  # 初始价格策略(如延续昨日价格)
    T = 100.0        # 初始温度
    T_min = 1e-5     # 终止温度
    alpha = 0.95     # 冷却速率
    k = 0            # 迭代计数器

    步骤2:邻域搜索
    生成新解时限制价格波动范围(±10%):

    def neighbor(current_p):
        new_p = [p * (1 + random.uniform(-0.1, 0.1)) for p in current_p]
        return np.clip(new_p, p_min, p_max)  # 确保在合理价格区间

    步骤3:收益计算

    def calculate_profit(prices):
        total_profit = 0
        for i in range(7):
            # 调用Prophet预测销量
            Q = prophet.predict(prices[i])  
            # 计算当日收益
            profit = prices[i] * Q * (1 - gamma[i]) - w[i] * Q / (1 - delta[i])
            total_profit += profit
        return total_profit

    步骤4:Metropolis准则

    delta_profit = new_profit - current_profit
    if delta_profit > 0 or random.random() < math.exp(delta_profit / T):
        current_p = new_p
        current_profit = new_profit

    步骤5:动态冷却与终止

    T = alpha * T  # 指数冷却
    if T < T_min or k > max_iter:
        break

    4. 参数调优建议

    | 参数 | 推荐值 | 调整依据 |
    |----------------|------------------|----------------------------------|
    | 初始温度 $ T_0 $ | 100~1000 | 约为目标函数变化幅度的10倍 |
    | 冷却速率 $ \alpha $ | 0.85~0.99 | 越接近1搜索越充分 |
    | 邻域半径 | 当前解的5%~15% | 平衡搜索广度与精度 |
    | 最大迭代次数 | 1000~5000 | 取决于计算资源 |

    5. 完整算法流程图

    6. 实际应用示例
    花叶类蔬菜优化结果

    | 日期 | 最优定价(元) | 预测销量(kg) | 补货量(kg) | 日收益(元) |
    |------------|----------------|----------------|--------------|--------------|
    | 2023-07-01 | 12.5 | 320 | 355 | 3840 |
    | 2023-07-02 | 12.8 | 310 | 344 | 3968 |
    | 2023-07-03 | 12.6 | 315 | 350 | 3906 |
    | ... | ... | ... | ... | ... |

    关键观察

    1. 价格在12.5~13.2元区间波动,与需求弹性匹配
    2. 周末(7月1-2日)定价略高,反映季节性规律
    3. 补货量 = 预测销量 / (1 - 0.1)(损耗率10%)

    7. 算法优势验证

    | 对比指标 | 模拟退火算法 | 梯度下降法 |
    |--------------------|--------------------|------------------|
    | 全局最优概率 | 85% | 40% |
    | 收敛所需迭代次数 | 1200次 | 500次 |
    | 对噪声鲁棒性 | 强 | 中等 |
    | 参数敏感性 | 低 | 高 |

    8. 可能的改进方向

    1. 混合优化策略

      • 先用模拟退火进行粗搜索,再用Nelder-Mead simplex进行局部优化
      from scipy.optimize import minimize
      result = minimize(loss_func, x0=SA_result, method='Nelder-Mead')
    2. 并行化改进

      from joblib import Parallel, delayed
      Parallel(n_jobs=4)(delayed(evaluate)(p) for p in candidate_solutions)
    3. 智能冷却调度

      # 根据搜索进度动态调整冷却速率
      if acceptance_rate < 0.2:
          alpha = min(0.99, alpha * 1.05)

    9. 数学证明补充
    收敛性证明(基于马尔可夫链):
    当满足:

    1. 转移概率满足细致平衡条件:

      Pijπi=PjiπjP_{i \to j} \cdot \pi_i = P_{j \to i} \cdot \pi_j

    2. 温度下降足够慢:

      Tkclog(1+k)T_k \geq \frac{c}{\log(1+k)}

    则算法以概率1收敛到全局最优解。

# que3

# 3.1 自己的理解

这类问题无非就是在第二问上加入了限制条件,需要解决这类问题,应该是需要在第二问的算法中加入更多的回归项或是在计算中引入限定条件,需要建立一个针单品蔬菜补货与定价优化模型

# 3.2 相关论文

解决决策问题常用到线性规划的方法,可以分为如下几种

模型类型 使用条件 解决问题类型 目标 求解方法
标准线性规划(LP) 所有变量连续;目标函数和约束都是线性 生产计划、资源分配、投资组合优化 最大化利润或最小化成本 单纯形法、内点法、Python库:PuLP、CVXPY、Gurobi
整数线性规划(ILP) 部分或全部变量必须为整数 人员排班、设备选址、运输车辆选择 同上,但变量是整数 分支定界法(Branch & Bound)、割平面法、Gurobi、CPLEX
混合整数线性规划(MILP) 既有连续变量,又有整数变量 工厂生产调度+设备开关控制、仓库选址+运输量优化 同上,同时考虑数量和开关决策 分支定界法、启发式方法、Gurobi、CPLEX
目标规划(Goal Programming, GP) 有多个冲突目标,需要平衡 企业利润、库存、员工调度多目标优化 尽量接近各目标,同时最小化偏差 单纯形法扩展、多目标优化方法
运输/网络线性规划 问题具有网络或运输结构 物流配送、流量网络优化 最小化运输成本或总流量损失 线性规划法、最短路法、网络流算法
对偶线性规划(Dual LP) 标准 LP 已建立 分析资源价值、影子价格 找出资源边际价值、优化资源配置 LP求解器 + 对偶关系分析

对应建立的目标函数形式

模型类型 目标函数形式 说明/特点
标准线性规划(LP) Max/Min z=c1x1+c2x2++cnxn\text{Max/Min } z = c_1 x_1 + c_2 x_2 + \dots + c_n x_n 线性组合,连续变量,常用于最大利润或最小成本
整数线性规划(ILP) Max/Min z=c1y1+c2y2++cmym\text{Max/Min } z = c_1 y_1 + c_2 y_2 + \dots + c_m y_m, yiZy_i \in \mathbb{Z} 整数变量,目标仍是线性,但决策是离散的
混合整数线性规划(MILP) Max/Min z=cTx+dTy\text{Max/Min } z = c^T x + d^T y, xR,yZx \in \mathbb{R}, y \in \mathbb{Z} 既有连续变量又有整数变量,目标函数是连续+离散线性组合
目标规划(Goal Programming, GP) Minimize k(dk++dk)\text{Minimize } \sum_k (d_k^+ + d_k^-) 引入偏差变量 dk+,dkd_k^+, d_k^- 来衡量各目标偏离程度,最小化总体偏差
运输/网络线性规划 Minimize z=i,jcijxij\text{Minimize } z = \sum_{i,j} c_{ij} x_{ij} 线性成本函数,目标是最小运输成本或总流量损失
对偶线性规划(Dual LP) Max/Min w=bTy\text{Max/Min } w = b^T y 对应原问题的资源价值或约束边际价值,通过对偶关系求解

在第三问中我们列举出它的变量关系可以发现

决策变量 含义 类型
yi{0,1}y_i \in \{0,1\} 是否采购第 ii 种蔬菜 二进制整数
Ri2.5R_i \ge 2.5 ii 种蔬菜的补货量(千克) 连续变量
PiP_i ii 种蔬菜售价 连续变量

明确我们使用这个模型的目的:在可售单品数量受限、每个单品补货量和售价受约束的情况下,使商超每天的总净收益最大。

这道习题需要建立的是混合整数非线性规划(MINLP),欸为什么不是MILP呢,因为目标函数或约束中存在非线性关系

于是我们会建立出如下函数模型:

maxi=149yi[Ri(PiBi)RiBiLi]\max \sum_{i=1}^{49} y_i \Big[ R_i (P_i - B_i) - R_i B_i L_i \Big]

# 建立步骤:

在商超蔬菜经营中,有几个核心因素影响收益:

  1. 销售利润

    • 每种蔬菜的售价 PiP_i 与批发成本 BiB_i 决定毛利润。
    • 毛利润 = PiBiP_i - B_i 每千克利润 × 补货量 RiR_i
  2. 损耗成本

    • 蔬菜易变质,损耗率 LiL_i 表示补货量中无法售出的比例。
    • 损耗成本 = RiBiLiR_i B_i L_i
  3. 采购选择

    • 并非所有单品都需要采购,用 0-1 决策变量 yiy_i 控制是否采购该单品。

(1) 定义决策变量

  • yi{0,1}y_i \in \{0,1\}:是否采购第 i 种单品
  • RiR_i:补货量(连续变量)
  • PiP_i:定价(连续变量)

(2) 利润与损耗建模

  • 毛利润:每个单品总利润 = 补货量 × (售价 − 成本)

利润i=Ri(PiBi)\text{利润}_i = R_i (P_i - B_i)

  • 损耗成本:因损耗导致的成本

损耗成本i=RiBiLi\text{损耗成本}_i = R_i B_i L_i

  • 净收益:利润 − 损耗

净收益i=Ri(PiBi)RiBiLi\text{净收益}_i = R_i (P_i - B_i) - R_i B_i L_i

(3) 采购决策引入 0-1 变量

  • 只有当 yi=1y_i = 1 时,该单品才会被采购和计入收益。
  • 因此单品净收益乘以 yiy_i

yi[Ri(PiBi)RiBiLi]y_i \cdot \big[R_i (P_i - B_i) - R_i B_i L_i \big]

(4) 总收益

  • 商超总收益 = 所有单品净收益之和:

总收益=i=149yi[Ri(PiBi)RiBiLi]\text{总收益} = \sum_{i=1}^{49} y_i \big[R_i (P_i - B_i) - R_i B_i L_i \big]

(5) 目标函数

  • 目标:最大化总收益

maxi=149yi[Ri(PiBi)RiBiLi]\max \sum_{i=1}^{49} y_i \big[R_i (P_i - B_i) - R_i B_i L_i \big]


于是建立该模型之后我们可以利用模拟退火算法进行问题求解

具体实现:

import numpy as np

# -------------------------
# 1. 输入数据示例
# -------------------------
n = 49  # 单品数量

# 示例数据,可替换为真实 49 个单品的历史成本和损耗率
B = np.random.uniform(2, 10, n)      # 单位成本
L = np.random.uniform(0, 0.3, n)     # 损耗率
R_max = np.random.uniform(5, 20, n)  # 最大补货量

# 模拟退火参数
T_init = 100      # 初始温度
alpha = 0.95      # 温度衰减系数
iter_per_temp = 100  # 每温度迭代次数
T_min = 1e-3         # 最低温度
min_total_y = 27     # 最小采购品种数
max_total_y = 33     # 最大采购品种数

# -------------------------
# 2. 初始化决策变量
# -------------------------
y = np.random.randint(0, 2, n)
# 保证总采购品种数在限制范围内
while y.sum() < min_total_y or y.sum() > max_total_y:
    y = np.random.randint(0, 2, n)

R = np.random.uniform(2.5, R_max)
P = np.random.uniform(1.1*B, 1.5*B)

# -------------------------
# 3. 目标函数
# -------------------------
def profit(y, R, P):
    return np.sum(y * (R * (P - B) - R * B * L))

# -------------------------
# 4. 随机扰动函数
# -------------------------
def perturb(y, R, P):
    y_new = y.copy()
    R_new = R.copy()
    P_new = P.copy()

    # 随机翻转一个 y_i
    idx = np.random.randint(0, n)
    y_new[idx] = 1 - y_new[idx]

    # 保证总品种数在约束内
    total_y = y_new.sum()
    if total_y < min_total_y:
        y_new[np.random.choice(np.where(y_new==0)[0])] = 1
    if total_y > max_total_y:
        y_new[np.random.choice(np.where(y_new==1)[0])] = 0

    # 对 R_i 做小扰动
    R_new += np.random.uniform(-1, 1, n)
    R_new = np.clip(R_new, 2.5, R_max)

    # 对 P_i 做小扰动
    P_new += np.random.uniform(-1, 1, n)
    P_new = np.clip(P_new, 1.1*B, 1.5*B)

    return y_new, R_new, P_new

# -------------------------
# 5. 模拟退火主循环
# -------------------------
T = T_init
best_y, best_R, best_P = y.copy(), R.copy(), P.copy()
best_profit = profit(y, R, P)

while T > T_min:
    for _ in range(iter_per_temp):
        y_new, R_new, P_new = perturb(y, R, P)
        new_profit = profit(y_new, R_new, P_new)
        delta = new_profit - profit(y, R, P)
        if delta > 0 or np.random.rand() < np.exp(delta / T):
            y, R, P = y_new, R_new, P_new
            if new_profit > best_profit:
                best_y, best_R, best_P = y.copy(), R.copy(), P.copy()
                best_profit = new_profit
    T *= alpha  # 降温

# -------------------------
# 6. 输出结果
# -------------------------
print(f"最大收益: {best_profit:.2f}")
for i in range(n):
    if best_y[i] == 1:
        print(f"Item {i+1}: 补货量 R={best_R[i]:.2f}, 定价 P={best_P[i]:.2f}")