2024辽宁省数学建模钢铁产品质量优化问题

最新推荐文章于 2024-09-11 19:01:38 发布

DeeSzeChyi

最新推荐文章于 2024-09-11 19:01:38 发布

阅读量495

点赞数 5

文章标签：数学建模

本文链接：https://blog.csdn.net/m0_64642912/article/details/140483533

版权

数据预处理：绘制折线图，缺失值无，异常值疑似5条记录

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import lightgbm as lgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
prepared=pd.read_excel(r"C:\Users\27734\Desktop\2024年辽宁省大学生数学建模竞赛赛题\2024年辽宁省数学建模竞赛赛题\2024年辽宁省大学生数学建模竞赛B题\附件1.xlsx",header=1)
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签

df = prepared  # 使用你的 DataFrame

# 生成x轴数据
x = np.arange(1, len(df) + 1)

# 列分组
column_groups = [df.columns[:4], df.columns[4:8], df.columns[8:]]

# 绘制每个列分组的折线图
for i, group in enumerate(column_groups):
    plt.figure(figsize=(14, 8))
    for column in group:
        plt.plot(x, df[column], label=column)
    # 添加图例
    plt.legend(loc='center left', bbox_to_anchor=(1, 0.5), fontsize=20)
    # 添加标题和标签
    plt.title(f"带钢参数变化折线图 - 第{i+1}组", fontsize=20)
    plt.xlabel("样本编号", fontsize=20)
    plt.ylabel("数值", fontsize=20)
    plt.xticks(fontsize=20)
    plt.yticks(fontsize=20)
    # 显示图形
    plt.show()

问题一：请帮助现场操作人员确定哪些参数对于带钢的机械性能具有重要影响？

通过计算皮尔逊相关系数，评估了各参数与性能指标之间的相关性，筛选出潜在重要的工艺参数。
然后，运用 LightGBM 模型深入分析了这些参数的影响力，发现加热炉温度，含碳
量，平整机张力，过时效炉温度，缓冷炉温度和均热炉温度对硬度的影响程度较为显
著，其他因素对硬度影响不是很显著但存在一定的影响。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import lightgbm as lgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
prepared=pd.read_excel(r"C:\Users\27734\Desktop\2024年辽宁省大学生数学建模竞赛赛题\2024年辽宁省数学建模竞赛赛题\2024年辽宁省大学生数学建模竞赛B题\附件1_new.xlsx",header=1)


plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号
plt.figure(figsize=(15,5),dpi=300)
data = prepared
# 计算皮尔逊相关系数
correlation_matrix = data.corr()

# 绘制热力图
plt.figure(figsize=(10, 8))  # 设置图形大小
sns.heatmap(correlation_matrix, annot=True, cmap='coolwarm', fmt=".2f", linewidths=0.5)
plt.title('皮尔逊相关系数热力图',fontsize=20)
x_new=np.arange(0.5,len(data.columns)+0.5)
x_label_ticks=data.columns
plt.xticks(x_new,x_label_ticks,rotation=90,fontsize=20)
plt.yticks(x_new,x_label_ticks,rotation=0,fontsize=20)
plt.show()


import pandas as pd
import lightgbm as lgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
import matplotlib.pyplot as plt
import seaborn as sns
import os

# 解决警告
os.environ["LOKY_MAX_CPU_COUNT"] = "4"  # 假设你有4个CPU核心

# 创建数据集
df=prepared

# 选取特征列和目标列
X = df.iloc[:, :-1]
y = df.iloc[:, -1]

# 获取特征名称列表
train_features = X.columns.tolist()

# 将数据分为训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# 创建LightGBM模型
clf = lgb.LGBMClassifier()

# 训练模型
clf.fit(X_train, y_train)

# 预测
y_pred = clf.predict(X_test)

# 计算准确率
accuracy = accuracy_score(y_test, y_pred)
print(f'Accuracy: {accuracy:.4f}')

# 获取特征重要性
feature_importance_values = clf.feature_importances_

# 创建一个 DataFrame 存储特征重要性
feature_importance = pd.DataFrame()
feature_importance['fea_name'] = train_features
feature_importance['fea_imp'] = feature_importance_values
feature_importance = feature_importance.sort_values('fea_imp', ascending=False)

# 绘制特征重要性的条形图
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.figure(figsize=[20, 10], dpi=100)
ax = sns.barplot(x='fea_name', y='fea_imp', data=feature_importance, palette='bright')

# 在每个柱子上添加文本标签
for p in ax.patches:
    ax.annotate(format(p.get_height(), '.0f'), 
                (p.get_x() + p.get_width() / 2., p.get_height()), 
                ha='center', va='center', 
                xytext=(0, 10), 
                textcoords='offset points',
                fontsize=20)

ax.set_xticklabels(labels=feature_importance['fea_name'], rotation=45, fontsize=20)
plt.yticks(fontsize=20)
plt.xlabel('特征', fontsize=20)
plt.ylabel('特征重要性', fontsize=20)
plt.title('特征重要性可视化', fontsize=20)
plt.tight_layout()
plt.show()

问题二：请帮助现场操作人员建立一个数据驱动的带钢产品质量在线检测模型，
并分析该模型的性能；

本文采用随机森林模型，表现效果中规中矩

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import lightgbm as lgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
prepared=pd.read_excel(r"C:\Users\27734\Desktop\2024年辽宁省大学生数学建模竞赛赛题\2024年辽宁省数学建模竞赛赛题\2024年辽宁省大学生数学建模竞赛B题\附件1_new.xlsx",header=1)


import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.ensemble import RandomForestRegressor
from sklearn import metrics
import shap

# 假设你已经准备好了数据集 prepared
dataset = prepared

# 准备训练数据
X = dataset.iloc[:,:-1].values
y = dataset.iloc[:,-1].values

# 将数据分为训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=0)

# 特征缩放
sc1 = StandardScaler()
X_train = sc1.fit_transform(X_train)
X_test = sc1.transform(X_test)

# 训练随机森林回归模型，扩大树的数量和最大深度，并行处理
regressor1 = RandomForestRegressor(n_estimators=200, random_state=0, n_jobs=-1)
regressor1.fit(X_train, y_train)
y_pred = regressor1.predict(X_test)

# 评估回归性能
print('Mean Absolute Error:', metrics.mean_absolute_error(y_test, y_pred))
print('Mean Squared Error:', metrics.mean_squared_error(y_test, y_pred))
print('Root Mean Squared Error:', np.sqrt(metrics.mean_squared_error(y_test, y_pred)))
print('R^2:', metrics.r2_score(y_test, y_pred))

问题三：现场操作人员通常是根据个人经验对带钢产品的工艺参数进行设定，你
能否帮助他们建立一个带钢工艺参数优化的解决方案？

利用粒子群优化算法搜索工艺过程的最优参数，同时利用质量检测模型检测在搜索过程中的每组参数预计会得到钢材质量，当钢材质量最优时，则其参数为工艺过程的最佳参数。
选用的钢材的硬度指标的高低来判断钢材的质量好坏，因此设定优化的目标是最大化钢
材硬度。由于并不知道钢材硬度与 12 个工艺指标之间的具体关系，于是可以使用问题二中
的检测模型来构造一个目标函数（使用问题二中的随机森林模型来评估每个粒子的当前适应度值。），同时各个工艺指标的取值范围为数据集中对应的最小值到最大值的区间。

其中 H(x) 为钢材的硬度，x1 至 x12 为各个工艺过程的代号

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import lightgbm as lgb
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号
prepared=pd.read_excel(r"C:\Users\27734\Desktop\2024年辽宁省大学生数学建模竞赛赛题\2024年辽宁省数学建模竞赛赛题\2024年辽宁省大学生数学建模竞赛B题\附件1_new.xlsx",header=1)
data=prepared


import pandas as pd

# 获取每列的最大值和最小值
min_values = data.min().tolist()
max_values = data.max().tolist()

# 输出结果
print("最小值列表:", min_values)
print("最大值列表:", max_values)


import pandas as pd
import numpy as np
from sklearn.ensemble import RandomForestRegressor
import pyswarms as ps
import matplotlib.pyplot as plt

# 准备数据
features = data[['带钢厚度', '带钢宽度', '碳含量', '硅含量', '带钢速度', '加热炉温度', '均热炉温度', '缓冷炉温度','过时效炉温度', '快冷炉温度', '淬火炉温度', '平整机张力']]
target = data['硬度']

# 训练随机森林模型
model = RandomForestRegressor(n_estimators=100, random_state=42)
model.fit(features, target)

# 定义约束范围
lb = [8020.0, 180.0, 42.0, 0.0, 451.0, 673.0, 619.5, 573.5, 312.0, 40.0, 28.0, 1962.35958]  # 下界
ub = [9570.0, 234.0, 615.0, 22.0, 679.0, 812.6, 725.0, 720.5, 416.5, 72.0, 60.0, 3064.760135] # 上界
bounds = (lb, ub)

# 定义PSO的目标函数
def objective(params):
    hardness = []
    for param_set in params:
        param_set_df = pd.DataFrame([param_set], columns=features.columns)
        hardness.append(model.predict(param_set_df)[0])
    return -np.array(hardness)  # 取负值是因为PSO是最小化目标函数

# 初始化PSO
options = {'c1': 0.5, 'c2': 0.3, 'w': 0.9}
optimizer = ps.single.GlobalBestPSO(n_particles=100, dimensions=12, options=options, bounds=bounds)

# 运行PSO，迭代200次
cost, pos = optimizer.optimize(objective, iters=300)

# 输出结果
print(f"最佳参数集：{pos}")
print(f"最佳硬度值：{-cost}")

# 可视化收敛过程
plt.plot(optimizer.cost_history)
plt.xlabel('迭代次数')
plt.ylabel('硬度指标的相反数')
plt.title('粒子群优化收敛过程')
plt.show()

也可以使用SA-PSO组合模型：

import pandas as pd
import numpy as np
from sklearn.ensemble import RandomForestRegressor
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号

# 假设 data 已经定义并包含数据
# 准备数据
features = data[['带钢厚度', '带钢宽度', '碳含量', '硅含量', '带钢速度', '加热炉温度', '均热炉温度', '缓冷炉温度', '过时效炉温度', '快冷炉温度', '淬火炉温度', '平整机张力']]
target = data['硬度']

# 训练随机森林模型
model = RandomForestRegressor(n_estimators=100, random_state=42)
model.fit(features, target)

# 定义约束范围
lb = [8020.0, 180.0, 42.0, 0.0, 451.0, 673.0, 619.5, 573.5, 312.0, 40.0, 28.0, 1962.35958]  # 下界
ub = [9570.0, 234.0, 615.0, 22.0, 679.0, 812.6, 725.0, 720.5, 416.5, 72.0, 60.0, 3064.760135] # 上界
bounds = (lb, ub)

# 定义PSO的目标函数
def objective(params):
    hardness = []
    for param_set in params:
        param_set_df = pd.DataFrame([param_set], columns=features.columns)
        hardness.append(model.predict(param_set_df)[0])
    return -np.array(hardness)  # 取负值是因为PSO是最小化目标函数

# 定义并初始化PSO
def pso(iter_max):
    m = 30  # 粒子数量
    d = 12  # 粒子维度
    ac1 = 2.0
    ac2 = 0.5
    w_max = 0.9
    w_min = 0.4

    # 初始化粒子位置和速度
    x = np.random.rand(m, d) * (np.array(ub) - np.array(lb)) + np.array(lb)
    v = np.random.rand(m, d) * (np.array(ub) - np.array(lb)) * 0.1

    # 初始化个体最优位置和最优值
    p = x.copy()
    pbest = objective(x)

    # 初始化全局最优位置和最优值
    g = x[np.argmin(pbest)].copy()
    gbest = np.min(pbest)

    # 记录收敛过程
    gb_history = []

    for i in range(iter_max):
        for j in range(m):
            # 更新个体最优值
            fitness = objective([x[j]])
            if fitness < pbest[j]:
                p[j] = x[j].copy()
                pbest[j] = fitness.item()  # 提取单个元素

            # 更新全局最优值
            if pbest[j] < gbest:
                g = p[j].copy()
                gbest = pbest[j]

            # 更新粒子速度和位置
            w = w_max - (w_max - w_min) * i / iter_max
            c1 = ac1 * (1 - i / iter_max)
            c2 = ac2 * (1 - i / iter_max)
            v[j] = w * v[j] + c1 * np.random.rand() * (p[j] - x[j]) + c2 * np.random.rand() * (g - x[j])
            x[j] = x[j] + v[j]

            # 边界处理
            x[j] = np.clip(x[j], lb, ub)
            v[j] = np.clip(v[j], -abs(np.array(ub) - np.array(lb)), abs(np.array(ub) - np.array(lb)))

        # 记录当前全局最优值
        gb_history.append(gbest)
        
    return gb_history

# 定义并初始化SA
def sa(iter_max):
    t_star = 100
    t_end = 0.01
    sa_k = 0.98  # 衰减因子
    t_cur = t_star
    
    x = np.random.rand(1, d) * (np.array(ub) - np.array(lb)) + np.array(lb)
    v = np.random.rand(1, d) * (np.array(ub) - np.array(lb)) * 0.1

    pbest = objective(x)[0]

    # 记录收敛过程
    gb_history = []

    for i in range(iter_max):
        for j in range(1):
            # 更新个体最优值
            fitness = objective([x[j]])
            deta_p = fitness - pbest
            if deta_p < 0 or np.exp(-deta_p / t_cur) > np.random.rand():
                pbest = fitness.item()  # 提取单个元素
                x[j] = x[j] + v[j]

            # 更新速度和位置
            w = w_max - (w_max - w_min) * i / iter_max
            c1 = ac1 * (1 - i / iter_max)
            c2 = ac2 * (1 - i / iter_max)
            v[j] = w * v[j] + c1 * np.random.rand() * (p[j] - x[j]) + c2 * np.random.rand() * (g - x[j])
            x[j] = x[j] + v[j]

            # 边界处理
            x[j] = np.clip(x[j], lb, ub)
            v[j] = np.clip(v[j], -abs(np.array(ub) - np.array(lb)), abs(np.array(ub) - np.array(lb)))

        # 记录当前全局最优值
        gb_history.append(pbest)
        t_cur *= sa_k
    
    return gb_history

# 初始化SA-PSO
def sa_pso(iter_max):
    m = 30  # 粒子数量
    d = 12  # 粒子维度
    ac1 = 2.0
    ac2 = 0.5
    sa_k = 0.98  # 衰减因子
    w_max = 0.9
    w_min = 0.4
    t_star = 100
    t_end = 0.01

    # 初始化粒子位置和速度
    x = np.random.rand(m, d) * (np.array(ub) - np.array(lb)) + np.array(lb)
    v = np.random.rand(m, d) * (np.array(ub) - np.array(lb)) * 0.1

    # 初始化个体最优位置和最优值
    p = x.copy()
    pbest = objective(x)

    # 初始化全局最优位置和最优值
    g = x[np.argmin(pbest)].copy()
    gbest = np.min(pbest)

    # 记录收敛过程
    gb_history = []

    # 模拟退火和粒子群优化结合
    t_cur = t_star
    for i in range(iter_max):
        for j in range(m):
            # 更新个体最优值
            fitness = objective([x[j]])
            if fitness < pbest[j]:
                p[j] = x[j].copy()
                pbest[j] = fitness.item()  # 提取单个元素
            else:
                deta_p = fitness - pbest[j]
                if deta_p < 0 or np.exp(-deta_p / t_cur) > np.random.rand():
                    p[j] = x[j].copy()
                    pbest[j] = fitness.item()  # 提取单个元素

            # 更新全局最优值
            if pbest[j] < gbest:
                g = p[j].copy()
                gbest = pbest[j]

            # 更新粒子速度和位置
            w = w_max - (w_max - w_min) * i / iter_max
            c1 = ac1 * (1 - i / iter_max)
            c2 = ac2 * (1 - i / iter_max)
            v[j] = w * v[j] + c1 * np.random.rand() * (p[j] - x[j]) + c2 * np.random.rand() * (g - x[j])
            x[j] = x[j] + v[j]

            # 边界处理
            x[j] = np.clip(x[j], lb, ub)
            v[j] = np.clip(v[j], -abs(np.array(ub) - np.array(lb)), abs(np.array(ub) - np.array(lb)))

        # 记录当前全局最优值
        gb_history.append(gbest)
        t_cur *= sa_k
        
    return gb_history

# 设置迭代次数
iter_max = 100

# 获取收敛过程
pso_history = pso(iter_max)
sa_history = sa(iter_max)
sa_pso_history = sa_pso(iter_max)

# 绘制收敛过程比较图
plt.plot(pso_history, label='粒子群优化 (PSO)')
plt.plot(sa_history, label='模拟退火 (SA)')
plt.plot(sa_pso_history, label='模拟退火-粒子群优化 (SA-PSO)')
plt.xlabel('迭代次数')
plt.ylabel('硬度指标的相反数')
plt.title('算法收敛过程比较')
plt.legend()
plt.show()