主题081_机器学习预测辐射物性

kkchenjj

388人浏览 · 2026-04-29 20:46:41

kkchenjj · 2026-04-29 20:46:41 发布

主题081：机器学习预测辐射物性

摘要

传统辐射物性测量和计算方法往往耗时耗力，难以满足材料快速筛选和设计的需求。机器学习技术为辐射物性预测提供了新的解决方案。本主题系统介绍如何利用机器学习方法预测复杂材料的辐射物性，包括发射率、吸收率、反射率和散射特性等。通过构建神经网络模型，结合材料基因组学和高通量计算方法，实现辐射物性的快速准确预测，为辐射换热领域的材料设计和优化提供强大工具。

关键词

机器学习、物性预测、神经网络、材料基因组、高通量计算、发射率预测、特征工程、模型优化

在这里插入图片描述

1. 引言：机器学习在辐射物性预测中的应用

1.1 传统辐射物性预测的挑战

辐射物性是描述材料与辐射相互作用的关键参数，包括发射率、吸收率、反射率、透射率和散射特性等。准确获取这些物性对于辐射换热的工程计算和系统设计至关重要。

传统方法的局限性：

实验测量：需要昂贵的仪器设备，测量周期长，难以覆盖广阔的材料空间
理论计算：第一性原理计算精度高但计算成本极高，只能处理简单体系
经验模型：适用范围有限，难以预测新材料的物性
数据库建设：传统数据库建设缓慢，数据稀疏且不平衡

1.2 机器学习的优势

机器学习为辐射物性预测带来了革命性的变化：

数据驱动

从大量实验和计算数据中学习规律
自动发现材料特征与物性之间的复杂关系
能够处理高维、非线性的映射关系

快速预测

训练后的模型预测速度极快（毫秒级）
适合高通量筛选和大规模材料设计
可实现实时物性预测和优化

知识发现

揭示影响辐射物性的关键材料特征
发现新的物理规律和材料设计原则
指导实验设计和材料合成

1.3 本主题内容安排

本主题系统介绍机器学习预测辐射物性的方法：

第2章介绍机器学习基础知识和常用算法
第3章介绍辐射物性数据的获取与预处理
第4章介绍特征工程方法
第5章介绍神经网络模型构建与训练
第6章介绍模型评估与优化策略
第7章通过6个案例演示机器学习预测辐射物性的应用

2. 机器学习基础

2.1 机器学习概述

机器学习是人工智能的一个分支，通过算法让计算机从数据中学习规律，并用于预测和决策。

监督学习

给定输入-输出对（训练数据），学习从输入到输出的映射函数：

$f(\mathbf{x}; \theta)$

其中 $x\mathbf{x}$ 是输入特征， $y$ 是目标输出， $θ\theta$ 是模型参数。

无监督学习

从无标签数据中发现隐藏的结构和模式，如聚类、降维等。

强化学习

通过与环境交互学习最优策略，适用于控制和优化问题。

2.2 常用机器学习算法

线性回归

最简单的回归模型：

$\mathbf{w}^T \mathbf{x} + b$

通过最小化均方误差（MSE）来优化参数：

$L=1N∑i=1N(yi−y^i)2 \mathcal{L} = \frac{1}{N} \sum_{i=1}^{N} (y_i - \hat{y}_i)^2$

决策树与随机森林

决策树通过递归划分特征空间进行预测。随机森林集成多棵决策树，提高预测精度和鲁棒性。

支持向量机（SVM）

寻找最优超平面将数据分隔：

$\min_{\mathbf{w}, b} \frac{1}{2} ||\mathbf{w}||^2 + C \sum_{i=1}^{N} \max(0, 1 - y_i(\mathbf{w}^T \mathbf{x}_i + b))$

神经网络

多层非线性变换：

$\mathbf{z}^{[l]} = \mathbf{W}^{[l]} \mathbf{a}^{[l-1]} + \mathbf{b}^{[l]}$

$\mathbf{a}^{[l]} = g(\mathbf{z}^{[l]})$

其中 $g$ 是激活函数，如ReLU、sigmoid、tanh等。

2.3 深度学习框架

多层感知机（MLP）

最基本的神经网络结构，由输入层、隐藏层和输出层组成。

卷积神经网络（CNN）

擅长处理具有空间结构的数据，如材料微观结构图像。

循环神经网络（RNN）

适合处理序列数据，如时间依赖的辐射特性。

图神经网络（GNN）

直接学习分子图结构，适合预测分子和晶体材料的物性。

3. 辐射物性数据获取与预处理

3.1 数据来源

实验数据库

NIST数据库：包含大量材料的光学常数
商业数据库：如Sopra、Palik等
文献数据：从已发表论文中提取

计算数据

第一性原理计算：DFT计算介电函数
分子动力学模拟：计算近场辐射特性
蒙特卡洛模拟：求解辐射传递问题

高通量计算

材料基因组计划（MGI）
AFLOW、Materials Project等数据库
自动化计算流程

3.2 数据预处理

数据清洗

处理缺失值：删除、插值或预测
异常值检测：使用统计方法或孤立森林
数据去重：去除重复样本

数据归一化

将特征缩放到相同范围：

Min-Max归一化：

$\frac{x - x_{min}}{x_{max} - x_{min}}$

Z-score标准化：

$\frac{x - \mu}{\sigma}$

数据划分

训练集：用于训练模型（通常60-80%）
验证集：用于调参和早停（通常10-20%）
测试集：用于最终评估（通常10-20%）

3.3 数据增强

物理约束增强

利用物理规律生成额外数据：

基尔霍夫定律： $ε=α\varepsilon = \alpha$ （热平衡时）
能量守恒： $R + A + T = 1$
互易关系：光谱对称性

合成数据生成

高斯噪声添加
插值和外推
生成对抗网络（GAN）

4. 特征工程

4.1 材料特征描述

组成特征

元素组成：原子序数、电负性、原子半径
化学计量比：各元素的比例
分子量：摩尔质量

结构特征

晶体结构：晶格常数、空间群
配位数：最近邻原子数
键长键角：原子间距离和角度

电子特征

能带结构：带隙、能带宽度
态密度：电子态分布
介电函数：实部和虚部

4.2 特征选择

相关性分析

计算特征与目标变量的相关系数：

$\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}}$

特征重要性

基于树模型的特征重要性
置换重要性（Permutation Importance）
SHAP值分析

降维方法

主成分分析（PCA）
t-SNE可视化
自编码器

4.3 案例：发射率预测的特征工程

问题描述

预测金属材料在红外波段的发射率。

特征构建

原子特征：原子序数Z、电负性χ、原子半径r
电子特征：价电子数、电离能、电子亲和能
晶体特征：晶格类型、配位数、密度
光谱特征：等离子体频率、阻尼系数

特征选择结果

最重要的特征：等离子体频率、阻尼系数、电负性
次要特征：原子序数、晶格类型
冗余特征：原子半径（与电负性相关）

5. 神经网络模型构建

5.1 网络架构设计

输入层

节点数等于特征维度
接收归一化后的特征向量

隐藏层

层数和节点数：通常2-5层隐藏层，每层50-500节点
激活函数：ReLU（隐藏层），Sigmoid/Tanh（输出层）

输出层

节点数等于预测目标数
发射率预测：1个节点（0-1范围）

5.2 损失函数与优化

损失函数

均方误差（MSE）：

$LMSE=1N∑i=1N(yi−y^i)2 \mathcal{L}_{MSE} = \frac{1}{N} \sum_{i=1}^{N} (y_i - \hat{y}_i)^2$

平均绝对误差（MAE）：

$LMAE=1N∑i=1N∣yi−y^i∣ \mathcal{L}_{MAE} = \frac{1}{N} \sum_{i=1}^{N} |y_i - \hat{y}_i|$

Huber损失：

结合MSE和MAE的优点：

$Lδ={12(y−y^)2if ∣y−y^∣≤δδ(∣y−y^∣−12δ)otherwise \mathcal{L}_{\delta} = \begin{cases} \frac{1}{2}(y - \hat{y})^2 & \text{if } |y - \hat{y}| \leq \delta \\ \delta(|y - \hat{y}| - \frac{1}{2}\delta) & \text{otherwise} \end{cases}$

优化算法

SGD：随机梯度下降，简单但收敛慢
Adam：自适应学习率，收敛快且稳定
L-BFGS：适合小规模问题

5.3 正则化技术

L2正则化（权重衰减）

$\mathcal{L} = \mathcal{L}_{data} + \lambda \sum_{j} w_j^2$

Dropout

训练时随机丢弃部分神经元，防止过拟合。

早停（Early Stopping）

验证集损失不再下降时停止训练。

批量归一化（Batch Normalization）

对每层输入进行归一化，加速训练。

6. 模型评估与优化

6.1 评估指标

回归问题指标

决定系数（R²）：

$R2=1−∑i=1N(yi−y^i)2∑i=1N(yi−yˉ)2 R^2 = 1 - \frac{\sum_{i=1}^{N} (y_i - \hat{y}_i)^2}{\sum_{i=1}^{N} (y_i - \bar{y})^2}$

均方根误差（RMSE）：

$RMSE=1N∑i=1N(yi−y^i)2 RMSE = \sqrt{\frac{1}{N} \sum_{i=1}^{N} (y_i - \hat{y}_i)^2}$

平均绝对百分比误差（MAPE）：

$MAPE=100%N∑i=1N∣yi−y^iyi∣ MAPE = \frac{100\%}{N} \sum_{i=1}^{N} \left|\frac{y_i - \hat{y}_i}{y_i}\right|$

6.2 交叉验证

K折交叉验证

将数据分为K份，轮流作为验证集：

$\frac{1}{K} \sum_{k=1}^{K} \mathcal{L}_k$

留一法（LOO）

每次留一个样本作为验证集，适合小数据集。

6.3 超参数优化

网格搜索

遍历所有超参数组合，计算成本高但全面。

随机搜索

随机采样超参数组合，效率更高。

贝叶斯优化

利用先验知识指导搜索，高效找到最优解。

常用超参数

学习率：0.001 - 0.1
批量大小：16 - 256
隐藏层节点数：50 - 500
正则化系数：1e-6 - 1e-2

7. Python仿真案例

案例1：发射率数据生成与预处理

生成模拟的发射率数据集，展示数据清洗、归一化和划分。

案例2：神经网络模型训练

构建并训练神经网络模型预测材料发射率。

案例3：模型性能评估与比较

比较不同机器学习算法（线性回归、随机森林、神经网络）的预测性能。

案例4：特征重要性分析

使用SHAP值分析各特征对发射率预测的贡献。

案例5：高通量筛选应用

利用训练好的模型对大量候选材料进行发射率预测和筛选。

案例6：机器学习预测辐射物性的动态演示

创建动画展示机器学习预测辐射物性的完整流程。

附录：Python代码说明

本主题配套的Python程序run_simulation.py实现了6个机器学习预测辐射物性的仿真案例：

依赖库

numpy：数值计算
matplotlib：数据可视化
scikit-learn：机器学习算法

运行方法

python run_simulation.py

输出文件

case1_data_preprocessing.png：数据预处理结果
case2_neural_network_training.png：神经网络训练过程
case3_model_comparison.png：不同模型性能比较
case4_feature_importance.png：特征重要性分析
case5_high_throughput_screening.png：高通量筛选结果
ml_emissivity_animation.gif：机器学习预测动态演示

#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""
主题081：机器学习预测辐射物性 - Python仿真程序
包含6个案例：
1. 发射率数据生成与预处理
2. 神经网络模型训练
3. 模型性能评估与比较
4. 特征重要性分析
5. 高通量筛选应用
6. 机器学习预测辐射物性的动态演示
"""

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle, Circle, FancyArrowPatch
import matplotlib.animation as animation
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.preprocessing import StandardScaler, MinMaxScaler
from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor
from sklearn.neural_network import MLPRegressor
from sklearn.metrics import mean_squared_error, r2_score, mean_absolute_error
import warnings
warnings.filterwarnings('ignore')

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans']
plt.rcParams['axes.unicode_minus'] = False


def generate_emissivity_data(n_samples=1000, noise_level=0.05, random_state=42):
    """
    生成模拟的发射率数据集
    基于Drude模型简化：发射率与等离子体频率和阻尼系数相关
    """
    np.random.seed(random_state)
    
    # 特征1: 等离子体频率 (eV)
    omega_p = np.random.uniform(5, 15, n_samples)
    
    # 特征2: 阻尼系数 (eV)
    gamma = np.random.uniform(0.05, 0.5, n_samples)
    
    # 特征3: 电负性
    electronegativity = np.random.uniform(0.8, 4.0, n_samples)
    
    # 特征4: 原子序数 (归一化)
    atomic_number = np.random.uniform(10, 80, n_samples)
    
    # 特征5: 晶格常数 (Å)
    lattice_constant = np.random.uniform(3, 6, n_samples)
    
    # 基于物理模型计算发射率（简化公式）
    # 发射率与阻尼系数正相关，与等离子体频率负相关
    emissivity = (
        0.3 + 
        0.4 * (gamma / 0.5) + 
        0.2 * (electronegativity / 4.0) -
        0.1 * ((omega_p - 10) / 5) +
        0.1 * ((atomic_number - 45) / 35) +
        0.05 * ((lattice_constant - 4.5) / 1.5)
    )
    
    # 添加噪声
    emissivity += np.random.normal(0, noise_level, n_samples)
    
    # 限制在0-1范围内
    emissivity = np.clip(emissivity, 0.01, 0.99)
    
    # 组合特征矩阵
    X = np.column_stack([omega_p, gamma, electronegativity, atomic_number, lattice_constant])
    
    feature_names = ['等离子体频率(eV)', '阻尼系数(eV)', '电负性', '原子序数', '晶格常数(Å)']
    
    return X, emissivity, feature_names


def case1_data_preprocessing():
    """案例1：发射率数据生成与预处理"""
    print("\n" + "=" * 70)
    print("案例1：发射率数据生成与预处理")
    print("=" * 70)
    
    # 生成数据
    X, y, feature_names = generate_emissivity_data(n_samples=1000, noise_level=0.05)
    
    print(f"\n数据集信息:")
    print(f"  样本数: {X.shape[0]}")
    print(f"  特征数: {X.shape[1]}")
    print(f"  特征列表: {', '.join(feature_names)}")
    
    # 数据划分
    X_train, X_test, y_train, y_test = train_test_split(
        X, y, test_size=0.2, random_state=42
    )
    
    print(f"\n数据划分:")
    print(f"  训练集: {X_train.shape[0]} 样本")
    print(f"  测试集: {X_test.shape[0]} 样本")
    
    # 数据归一化
    scaler = StandardScaler()
    X_train_scaled = scaler.fit_transform(X_train)
    X_test_scaled = scaler.transform(X_test)
    
    # 可视化
    fig, axes = plt.subplots(2, 3, figsize=(15, 10))
    
    # 子图1: 发射率分布
    ax1 = axes[0, 0]
    ax1.hist(y, bins=30, color='skyblue', edgecolor='black', alpha=0.7)
    ax1.set_xlabel('发射率', fontsize=12)
    ax1.set_ylabel('频数', fontsize=12)
    ax1.set_title('发射率分布', fontsize=14)
    ax1.axvline(x=np.mean(y), color='r', linestyle='--', label=f'均值: {np.mean(y):.3f}')
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 子图2-6: 各特征与发射率的关系
    for i, (ax, name) in enumerate(zip(axes.flat[1:], feature_names)):
        ax.scatter(X[:, i], y, alpha=0.5, s=20)
        ax.set_xlabel(name, fontsize=12)
        ax.set_ylabel('发射率', fontsize=12)
        ax.set_title(f'{name} vs 发射率', fontsize=12)
        ax.grid(True, alpha=0.3)
        
        # 计算相关系数
        corr = np.corrcoef(X[:, i], y)[0, 1]
        ax.text(0.05, 0.95, f'r = {corr:.3f}', transform=ax.transAxes, 
                fontsize=10, verticalalignment='top',
                bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
    
    plt.tight_layout()
    plt.savefig('case1_data_preprocessing.png', dpi=150, bbox_inches='tight')
    print("\n结果已保存至: case1_data_preprocessing.png")
    plt.close()
    
    # 打印统计信息
    print("\n特征统计信息:")
    print("-" * 70)
    print(f"{'特征':<20} {'均值':<12} {'标准差':<12} {'最小值':<12} {'最大值':<12}")
    print("-" * 70)
    for i, name in enumerate(feature_names):
        print(f"{name:<20} {np.mean(X[:, i]):<12.3f} {np.std(X[:, i]):<12.3f} "
              f"{np.min(X[:, i]):<12.3f} {np.max(X[:, i]):<12.3f}")
    
    return X_train_scaled, X_test_scaled, y_train, y_test, feature_names, scaler


def case2_neural_network_training(X_train, X_test, y_train, y_test):
    """案例2：神经网络模型训练"""
    print("\n" + "=" * 70)
    print("案例2：神经网络模型训练")
    print("=" * 70)
    
    # 创建神经网络模型
    mlp = MLPRegressor(
        hidden_layer_sizes=(100, 50, 25),  # 三层隐藏层
        activation='relu',
        solver='adam',
        alpha=0.001,  # L2正则化
        learning_rate_init=0.001,
        max_iter=500,
        early_stopping=True,
        validation_fraction=0.1,
        random_state=42,
        verbose=False
    )
    
    print("\n神经网络架构:")
    print(f"  输入层: {X_train.shape[1]} 个节点")
    print(f"  隐藏层1: 100 个节点")
    print(f"  隐藏层2: 50 个节点")
    print(f"  隐藏层3: 25 个节点")
    print(f"  输出层: 1 个节点")
    print(f"  激活函数: ReLU")
    print(f"  优化器: Adam")
    print(f"  学习率: 0.001")
    
    # 训练模型
    print("\n开始训练...")
    mlp.fit(X_train, y_train)
    
    # 预测
    y_train_pred = mlp.predict(X_train)
    y_test_pred = mlp.predict(X_test)
    
    # 计算指标
    train_r2 = r2_score(y_train, y_train_pred)
    test_r2 = r2_score(y_test, y_test_pred)
    train_rmse = np.sqrt(mean_squared_error(y_train, y_train_pred))
    test_rmse = np.sqrt(mean_squared_error(y_test, y_test_pred))
    
    print(f"\n训练完成!")
    print(f"  训练集 R²: {train_r2:.4f}")
    print(f"  测试集 R²: {test_r2:.4f}")
    print(f"  训练集 RMSE: {train_rmse:.4f}")
    print(f"  测试集 RMSE: {test_rmse:.4f}")
    print(f"  迭代次数: {mlp.n_iter_}")
    
    # 可视化
    fig, axes = plt.subplots(1, 2, figsize=(14, 6))
    
    # 子图1: 训练损失曲线
    ax1 = axes[0]
    ax1.plot(mlp.loss_curve_, 'b-', linewidth=2)
    ax1.set_xlabel('迭代次数', fontsize=12)
    ax1.set_ylabel('损失函数值', fontsize=12)
    ax1.set_title('神经网络训练损失曲线', fontsize=14)
    ax1.grid(True, alpha=0.3)
    
    # 子图2: 预测 vs 真实值
    ax2 = axes[1]
    ax2.scatter(y_test, y_test_pred, alpha=0.6, s=30, label='测试集')
    ax2.plot([0, 1], [0, 1], 'r--', linewidth=2, label='理想预测')
    ax2.set_xlabel('真实发射率', fontsize=12)
    ax2.set_ylabel('预测发射率', fontsize=12)
    ax2.set_title(f'预测性能 (R² = {test_r2:.4f})', fontsize=14)
    ax2.legend()
    ax2.grid(True, alpha=0.3)
    ax2.set_xlim(0, 1)
    ax2.set_ylim(0, 1)
    
    plt.tight_layout()
    plt.savefig('case2_neural_network_training.png', dpi=150, bbox_inches='tight')
    print("\n结果已保存至: case2_neural_network_training.png")
    plt.close()
    
    return mlp


def case3_model_comparison(X_train, X_test, y_train, y_test):
    """案例3：模型性能评估与比较"""
    print("\n" + "=" * 70)
    print("案例3：模型性能评估与比较")
    print("=" * 70)
    
    # 定义模型
    models = {
        '线性回归': LinearRegression(),
        '随机森林': RandomForestRegressor(n_estimators=100, random_state=42),
        '神经网络': MLPRegressor(hidden_layer_sizes=(100, 50), alpha=0.001, 
                               max_iter=500, random_state=42)
    }
    
    results = {}
    
    print("\n模型训练与评估:")
    print("-" * 70)
    print(f"{'模型':<15} {'训练R²':<12} {'测试R²':<12} {'训练RMSE':<12} {'测试RMSE':<12}")
    print("-" * 70)
    
    for name, model in models.items():
        # 训练
        model.fit(X_train, y_train)
        
        # 预测
        y_train_pred = model.predict(X_train)
        y_test_pred = model.predict(X_test)
        
        # 计算指标
        train_r2 = r2_score(y_train, y_train_pred)
        test_r2 = r2_score(y_test, y_test_pred)
        train_rmse = np.sqrt(mean_squared_error(y_train, y_train_pred))
        test_rmse = np.sqrt(mean_squared_error(y_test, y_test_pred))
        
        results[name] = {
            'train_r2': train_r2,
            'test_r2': test_r2,
            'train_rmse': train_rmse,
            'test_rmse': test_rmse,
            'y_test_pred': y_test_pred
        }
        
        print(f"{name:<15} {train_r2:<12.4f} {test_r2:<12.4f} "
              f"{train_rmse:<12.4f} {test_rmse:<12.4f}")
    
    # 可视化
    fig, axes = plt.subplots(1, 3, figsize=(15, 5))
    
    for idx, (name, result) in enumerate(results.items()):
        ax = axes[idx]
        ax.scatter(y_test, result['y_test_pred'], alpha=0.6, s=30)
        ax.plot([0, 1], [0, 1], 'r--', linewidth=2)
        ax.set_xlabel('真实发射率', fontsize=12)
        ax.set_ylabel('预测发射率', fontsize=12)
        ax.set_title(f'{name}\n(R² = {result["test_r2"]:.4f})', fontsize=12)
        ax.grid(True, alpha=0.3)
        ax.set_xlim(0, 1)
        ax.set_ylim(0, 1)
    
    plt.tight_layout()
    plt.savefig('case3_model_comparison.png', dpi=150, bbox_inches='tight')
    print("\n结果已保存至: case3_model_comparison.png")
    plt.close()
    
    return results


def case4_feature_importance(X_train, y_train, feature_names):
    """案例4：特征重要性分析"""
    print("\n" + "=" * 70)
    print("案例4：特征重要性分析")
    print("=" * 70)
    
    # 使用随机森林计算特征重要性
    rf = RandomForestRegressor(n_estimators=100, random_state=42)
    rf.fit(X_train, y_train)
    
    importances = rf.feature_importances_
    indices = np.argsort(importances)[::-1]
    
    print("\n特征重要性排名:")
    print("-" * 50)
    for i, idx in enumerate(indices):
        print(f"{i+1}. {feature_names[idx]:<20} {importances[idx]:.4f}")
    
    # 计算相关系数
    correlations = []
    for i in range(X_train.shape[1]):
        corr = np.corrcoef(X_train[:, i], y_train)[0, 1]
        correlations.append(abs(corr))
    
    # 可视化
    fig, axes = plt.subplots(1, 2, figsize=(14, 6))
    
    # 子图1: 随机森林特征重要性
    ax1 = axes[0]
    sorted_names = [feature_names[i] for i in indices]
    sorted_importances = importances[indices]
    
    bars = ax1.barh(range(len(sorted_names)), sorted_importances, color='steelblue')
    ax1.set_yticks(range(len(sorted_names)))
    ax1.set_yticklabels(sorted_names)
    ax1.set_xlabel('重要性', fontsize=12)
    ax1.set_title('随机森林特征重要性', fontsize=14)
    ax1.invert_yaxis()
    ax1.grid(True, alpha=0.3, axis='x')
    
    # 添加数值标签
    for i, (bar, val) in enumerate(zip(bars, sorted_importances)):
        ax1.text(val + 0.01, i, f'{val:.3f}', va='center', fontsize=10)
    
    # 子图2: 相关系数
    ax2 = axes[1]
    corr_indices = np.argsort(correlations)[::-1]
    sorted_corr_names = [feature_names[i] for i in corr_indices]
    sorted_correlations = [correlations[i] for i in corr_indices]
    
    bars2 = ax2.barh(range(len(sorted_corr_names)), sorted_correlations, color='coral')
    ax2.set_yticks(range(len(sorted_corr_names)))
    ax2.set_yticklabels(sorted_corr_names)
    ax2.set_xlabel('|相关系数|', fontsize=12)
    ax2.set_title('特征与发射率的相关性', fontsize=14)
    ax2.invert_yaxis()
    ax2.grid(True, alpha=0.3, axis='x')
    
    # 添加数值标签
    for i, (bar, val) in enumerate(zip(bars2, sorted_correlations)):
        ax2.text(val + 0.01, i, f'{val:.3f}', va='center', fontsize=10)
    
    plt.tight_layout()
    plt.savefig('case4_feature_importance.png', dpi=150, bbox_inches='tight')
    print("\n结果已保存至: case4_feature_importance.png")
    plt.close()


def case5_high_throughput_screening(model, scaler, feature_names):
    """案例5：高通量筛选应用"""
    print("\n" + "=" * 70)
    print("案例5：高通量筛选应用")
    print("=" * 70)
    
    # 生成大量候选材料
    n_candidates = 10000
    np.random.seed(123)
    
    # 随机生成候选材料
    omega_p_candidates = np.random.uniform(5, 15, n_candidates)
    gamma_candidates = np.random.uniform(0.05, 0.5, n_candidates)
    electronegativity_candidates = np.random.uniform(0.8, 4.0, n_candidates)
    atomic_number_candidates = np.random.uniform(10, 80, n_candidates)
    lattice_constant_candidates = np.random.uniform(3, 6, n_candidates)
    
    X_candidates = np.column_stack([
        omega_p_candidates, gamma_candidates, electronegativity_candidates,
        atomic_number_candidates, lattice_constant_candidates
    ])
    
    # 标准化
    X_candidates_scaled = scaler.transform(X_candidates)
    
    # 预测发射率
    emissivity_pred = model.predict(X_candidates_scaled)
    
    # 筛选高发射率材料（发射率 > 0.7）
    high_emissivity_mask = emissivity_pred > 0.7
    high_emissivity_indices = np.where(high_emissivity_mask)[0]
    
    print(f"\n高通量筛选结果:")
    print(f"  候选材料总数: {n_candidates}")
    print(f"  高发射率材料数: {len(high_emissivity_indices)}")
    print(f"  筛选比例: {len(high_emissivity_indices)/n_candidates*100:.2f}%")
    
    # 可视化
    fig, axes = plt.subplots(1, 2, figsize=(14, 6))
    
    # 子图1: 发射率分布
    ax1 = axes[0]
    ax1.hist(emissivity_pred, bins=50, color='lightblue', edgecolor='black', alpha=0.7)
    ax1.axvline(x=0.7, color='r', linestyle='--', linewidth=2, label='筛选阈值 (0.7)')
    ax1.set_xlabel('预测发射率', fontsize=12)
    ax1.set_ylabel('材料数量', fontsize=12)
    ax1.set_title('候选材料发射率分布', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 子图2: 高发射率材料的特征分布
    ax2 = axes[1]
    
    if len(high_emissivity_indices) > 0:
        # 绘制高发射率材料的特征范围
        high_emissivity_features = X_candidates[high_emissivity_indices]
        
        for i, name in enumerate(feature_names):
            ax2.scatter(
                [i] * len(high_emissivity_indices),
                high_emissivity_features[:, i],
                alpha=0.3, s=10
            )
        
        ax2.set_xticks(range(len(feature_names)))
        ax2.set_xticklabels([name.split('(')[0] for name in feature_names], rotation=45, ha='right')
        ax2.set_ylabel('特征值', fontsize=12)
        ax2.set_title('高发射率材料的特征分布', fontsize=14)
        ax2.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('case5_high_throughput_screening.png', dpi=150, bbox_inches='tight')
    print("\n结果已保存至: case5_high_throughput_screening.png")
    plt.close()
    
    # 打印前10个最佳候选材料
    if len(high_emissivity_indices) > 0:
        top_indices = high_emissivity_indices[np.argsort(emissivity_pred[high_emissivity_indices])[-10:]]
        
        print("\n前10个高发射率候选材料:")
        print("-" * 80)
        print(f"{'排名':<6} {'发射率':<10} {'ωp(eV)':<10} {'γ(eV)':<10} {'电负性':<10} {'原子序数':<10} {'晶格常数':<10}")
        print("-" * 80)
        for rank, idx in enumerate(reversed(top_indices), 1):
            print(f"{rank:<6} {emissivity_pred[idx]:<10.4f} "
                  f"{X_candidates[idx, 0]:<10.2f} {X_candidates[idx, 1]:<10.3f} "
                  f"{X_candidates[idx, 2]:<10.2f} {X_candidates[idx, 3]:<10.1f} "
                  f"{X_candidates[idx, 4]:<10.2f}")


def case6_ml_animation():
    """案例6：机器学习预测辐射物性的动态演示"""
    print("\n" + "=" * 70)
    print("案例6：机器学习预测辐射物性的动态演示")
    print("=" * 70)
    print("\n正在生成动画...")
    
    fig = plt.figure(figsize=(14, 10))
    gs = fig.add_gridspec(2, 2)
    
    ax1 = fig.add_subplot(gs[0, 0])  # 数据收集
    ax2 = fig.add_subplot(gs[0, 1])  # 特征工程
    ax3 = fig.add_subplot(gs[1, 0])  # 模型训练
    ax4 = fig.add_subplot(gs[1, 1])  # 预测应用
    
    n_frames = 50
    
    def animate(frame):
        # 清除之前的图形
        ax1.clear()
        ax2.clear()
        ax3.clear()
        ax4.clear()
        
        t = frame / n_frames
        
        # 子图1: 数据收集过程
        ax1.set_xlim(0, 10)
        ax1.set_ylim(0, 10)
        
        # 模拟数据点逐渐增加
        n_points = int(50 * t)
        np.random.seed(42)
        x_data = np.random.random(n_points) * 8 + 1
        y_data = np.random.random(n_points) * 8 + 1
        
        ax1.scatter(x_data, y_data, c='blue', s=50, alpha=0.6)
        ax1.set_xlabel('特征1')
        ax1.set_ylabel('特征2')
        ax1.set_title(f'数据收集 ({n_points} 个样本)')
        ax1.grid(True, alpha=0.3)
        
        # 子图2: 特征重要性
        features = ['等离子体\n频率', '阻尼\n系数', '电负性', '原子\n序数', '晶格\n常数']
        importances = [0.35, 0.28, 0.18, 0.12, 0.07]
        
        # 动态显示重要性
        current_importances = [imp * min(1, t * 2) for imp in importances]
        bars = ax2.bar(features, current_importances, color='steelblue')
        ax2.set_ylabel('重要性', fontsize=12)
        ax2.set_title('特征重要性分析', fontsize=14)
        ax2.set_ylim(0, 0.5)
        ax2.grid(True, alpha=0.3, axis='y')
        
        # 添加数值标签
        for bar, val in zip(bars, current_importances):
            if val > 0.05:
                ax2.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,
                        f'{val:.2f}', ha='center', fontsize=9)
        
        # 子图3: 训练损失曲线
        epochs = np.linspace(0, 500, 100)
        loss = np.exp(-epochs / 100) * (1 + 0.1 * np.sin(epochs / 10))
        current_epoch = int(500 * t)
        
        ax3.plot(epochs[:current_epoch+1], loss[:current_epoch+1], 'b-', linewidth=2)
        ax3.set_xlabel('迭代次数', fontsize=12)
        ax3.set_ylabel('损失函数值', fontsize=12)
        ax3.set_title(f'模型训练 (迭代 {current_epoch})', fontsize=14)
        ax3.grid(True, alpha=0.3)
        ax3.set_xlim(0, 500)
        ax3.set_ylim(0, 1.2)
        
        # 子图4: 预测应用
        ax4.set_xlim(0, 10)
        ax4.set_ylim(0, 10)
        
        # 绘制预测结果
        n_predictions = int(20 * t)
        true_values = np.random.random(n_predictions)
        pred_values = true_values + np.random.normal(0, 0.05, n_predictions)
        pred_values = np.clip(pred_values, 0, 1)
        
        ax4.scatter(true_values * 8 + 1, pred_values * 8 + 1, c='green', s=80, alpha=0.6)
        ax4.plot([1, 9], [1, 9], 'r--', linewidth=2, label='理想预测')
        ax4.set_xlabel('真实发射率', fontsize=12)
        ax4.set_ylabel('预测发射率', fontsize=12)
        ax4.set_title(f'预测应用 ({n_predictions} 个预测)', fontsize=14)
        ax4.legend()
        ax4.grid(True, alpha=0.3)
        
        plt.tight_layout()
        
        return ax1, ax2, ax3, ax4
    
    # 创建动画
    anim = animation.FuncAnimation(fig, animate, frames=n_frames, interval=200, blit=False)
    anim.save('ml_emissivity_animation.gif', writer='pillow', fps=5, dpi=100)
    plt.close()
    print("动画已保存至: ml_emissivity_animation.gif")


def main():
    """主函数 - 运行所有案例"""
    print("\n" + "=" * 70)
    print("主题081：机器学习预测辐射物性 - Python仿真")
    print("=" * 70)
    
    # 案例1: 数据生成与预处理
    X_train, X_test, y_train, y_test, feature_names, scaler = case1_data_preprocessing()
    
    # 案例2: 神经网络训练
    mlp_model = case2_neural_network_training(X_train, X_test, y_train, y_test)
    
    # 案例3: 模型比较
    results = case3_model_comparison(X_train, X_test, y_train, y_test)
    
    # 案例4: 特征重要性
    case4_feature_importance(X_train, y_train, feature_names)
    
    # 案例5: 高通量筛选
    case5_high_throughput_screening(mlp_model, scaler, feature_names)
    
    # 案例6: 动画演示
    case6_ml_animation()
    
    print("\n" + "=" * 70)
    print("所有仿真案例已完成！")
    print("=" * 70)
    print("\n生成的文件：")
    print("  1. case1_data_preprocessing.png")
    print("  2. case2_neural_network_training.png")
    print("  3. case3_model_comparison.png")
    print("  4. case4_feature_importance.png")
    print("  5. case5_high_throughput_screening.png")
    print("  6. ml_emissivity_animation.gif")


if __name__ == "__main__":
    main()

AtomGit开源社区

AtomGit 是由开放原子开源基金会联合 CSDN 等生态伙伴共同推出的新一代开源与人工智能协作平台。平台坚持“开放、中立、公益”的理念，把代码托管、模型共享、数据集托管、智能体开发体验和算力服务整合在一起，为开发者提供从开发、训练到部署的一站式体验。

更多推荐

【Qt教程29】Qt5和Qt6版本对比

原创作者：郑同学的笔记原文链接：https://zhengjunxue.blog.csdn.net/article/details/159691226Qt6 与之前的版本相比，无论是在底层架构还是上层 API 上，都带来了根本性的变革。为了帮你快速了解，这里从整体和内部两个层面，梳理了 Qt6 各版本间的核心差异。