SISSO构建描述符的使用

qq_46264842

已于 2024-12-02 15:21:33 修改

阅读量2.8k

点赞数 40

分类专栏： SISSO 文章标签： python 机器学习

于 2024-12-02 15:17:16 首次发布

本文链接：https://blog.csdn.net/qq_46264842/article/details/144188738

版权

SISSO 专栏收录该内容

1 篇文章

订阅专栏

SISSO构建描述符的使用

编写日期：2024年10月
编写这个的原因主要是记录SISSO的使用，由于时间长了可能忘了SISSO的操作，便编写了这个笔记。

项目开源地址：rouyang2017/SISSO

国内镜像：项目首页 - GitCode

文章链接：R. Ouyang, S. Curtarolo, E. Ahmetcik, M. Scheffler, and L. M. Ghiringhelli, Phys. Rev. Mater. 2, 083802 (2018).

学习参考：简单了解 SISSO 特征筛选方法在回归问题下应用

1. 安装

1.1 安装依赖

intel编译器
与之前编译vasp的一样，安装最新版例如：

sh l_BaseKit_p_2024.1.0.596_offline.sh -a --silent --eula accept
sh l_HPCKit_p_2024.1.0.560_offline.sh -a --silent --eula accept
sh l_fortran-compiler_p_2024.1.0.465_offline.sh -a --silent --eula accept

1.2 下载编译

例如在Ubuntu22系统下：

(base) rcz@rczhh:~/app$ git clone https://github.com/rouyang2017/SISSO.git
正克隆到 'SISSO'...
remote: Enumerating objects: 681, done.
remote: Counting objects: 100% (201/201), done.
remote: Compressing objects: 100% (159/159), done.
remote: Total 681 (delta 103), reused 87 (delta 42), pack-reused 480 (from 1)
接收对象中: 100% (681/681), 4.01 MiB | 116.00 KiB/s, 完成.
处理 delta 中: 100% (355/355), 完成.
(base) rcz@rczhh:~/app/SISSO$ ls
Changelog        LICENSE    SISSO_Guide_v3.5.pdf  utilities
input_templates  README.md  src
(base) rcz@rczhh:~/app/SISSO$ cd src/
(base) rcz@rczhh:~/app/SISSO/src$ mpiifort -fp-model precise var_global.f90 libsisso.f90 DI.f90 FC.f90 FCse.f90 SISSO.f90 -o ~/app/SISSO/bin/SISSO
ifort: remark #10448: Intel(R) Fortran Compiler Classic (ifort) is now deprecated and will be discontinued late 2024. Intel recommends that customers transition now to using the LLVM-based Intel(R) Fortran Compiler (ifx) for continued Windows* and Linux* support, new language support, new language features, and optimizations. Use '-diag-disable=10448' to disable this message.
(base) rcz@rczhh:~/app/SISSO/src$ cd ~/app/SISSO/bin
(base) rcz@rczhh:~/app/SISSO/bin$ ls
SISSO
(base) rcz@rczhh:~/app/SISSO/bin$ cd ../input_templates/
(base) rcz@rczhh:~/app/SISSO/input_templates$ mkdir test && cp train.dat_regression ./test/train.dat && cp SISSO.in ./test/
(base) rcz@rczhh:~/app/SISSO/input_templates$ cd test/ && ulimit -s unlimited
(base) rcz@rczhh:~/app/SISSO/input_templates/test$ ls
SISSO.in  train.dat
(base) rcz@rczhh:~/app/SISSO/input_templates/test$ ../../bin/SISSO > log
(base) rcz@rczhh:~/app/SISSO/input_templates/test$ ls
log  Models  SISSO.in  SISSO.out  SIS_subspaces  train.dat
(base) rcz@rczhh:~/app/SISSO/input_templates/test$ less SISSO.out
****************************************************************
  Sure Independence Screening and Sparsifying Operator (SISSO)
             Version SISSO.3.5, August, 2024.
****************************************************************

Read in data from train.dat
Standard Deviation (sqrt[sum(y-y_mean)^2/n]) of the target property:   1.07688
Read in data from SISSO.in
Property type:     1
Number of tasks:        1
Descriptor dimension:        2
Number of samples for the task(s):      5
Restarts :  0
Scheme for feature storing in memory:        1
Number of scalar features:        3
Tier of the feature space (max depth of expression tree):        2
Maximal feature complexity (number of operators in a feature):        3
Unit of input primary feature, each represented by a row vector:
    1.00    0.00
    1.00    0.00
    0.00    1.00
The feature will be discarded if the minimum of the maximal abs. value in it <    0.10000E-02
The faature will be discarded if the maximum of the maximal abs. value in it >     0.10000E+06
Size of the SIS-selected (single) subspace :  50000 50000
Operators for feature construction: (+)(-)(*)(/) (+)(-)(*)(/)
Method for sparse regression:  L0
Fitting intercept:      T
Metric for model selection: RMSE
Number of the top-ranked models to output:      100
--------------------------------------------------------------------------------

Dimension:   1
-------------------
Feature Construction (FC) starts ...
Total number of features in the space Phi00:              3
Total number of features in the space Phi01:             15
Total number of features in the space Phi02:            335
Size of the SIS-selected subspace from Phi02:        141
Time (second) used for this FC:            0.00

Descriptor Identification (DI) starts ...
Total number of SIS-selected features for this DI:        141
Time (second) used for this DI:            0.00

# OK SISSO安装基本成功了

# 编译完成后 将SISSO的程序放入环境变量当中  全局使用SISSO
# 打开用户下的 .bashrc 文件  在文件后添加
export PATH=~/app/SISSO/bin:$PATH    #这里添加你编译好的SISSO文件的地址
source ~/.bashrc

遇到了问题多半是intel编译器的安装错误。

2. 使用方法

具体内容请看官方文档（SISSO_Guide_v3.5.pdf）

这里对官方文档的ai翻译，仅供参考概不负责

2.1 介绍

SISSO 是一种数据驱动的方，结合了符号回归和压缩感应技术在材料信息学中以明确的数学表达式识别模型和描述符，该方法最初在 2018 年提出 [Phys. Rev. Mater. 2, 083802 (2018)]。目前有数个可用的 SISSO 方法实现，此用户指南针对 Runhai Ouyang 在 https://github.com/rouyang2017/SISSO 创建并维护的 SISSO FORTRAN 代码。当前 SISSO 代码的功能包括：

Regression & Classification 回归与分类
[R. Ouyang, S. Curtarolo, E. Ahmetcik, M. Scheffler, L. M. Ghiringhelli, Phys. Rev. Mater. 2, 083802 (2018)]
Multi-task learning 多任务学习
[R. Ouyang, E. Ahmetcik, C. Carbogno, M. Scheffler, L. M. Ghiringhelli, J. Phys.: Mater. 2, 024002 (2019)] [J. Wang, H. Xie, Y. Wang, R. Ouyang, J. Am. Chem. Soc. 145, 11457 (2023)]
Variable selection assisted symbolic regression 变量选择辅助符号回归
[Z. Guo, S. Hu, Z.-K. Han, R. Ouyang, J. Chem. Theory Comput. 18, 4945 (2022)]

这个 SISSO 代码以 Apache 许可，版本 2.0 进行许可。想要了解有关使用、再生产和分发的条款和条件的更多详情，请参考：http://www.apache.org/licenses/

许可权限：

商业用途
修改
散布
专利使用
私人用途

限制：

× 商标使用
× 责任
× 保修

条件：

 许可证和版权声明
 状态更改

2.2 安装

要编译 SISSO 并行程序，需要一个 FORTRAN MPI 编译器。以下是使用英特尔 MPI 编译器编译该程序的两个选项（其他编译器也可能起作用，但建议使用英特尔编译器）。在 Linux 操作系统中，下载并解压软件包，然后进入 ‘src’ 文件夹，然后执行以下操作：

1)  mpiifort -fp-model precise var_global.f90 libsisso.f90 DI.f90 FC.f90 FCse.f90 SISSO.f90 -o ~/bin/SISSO
2)  mpiifort -O2 var_global.f90 libsisso.f90 DI.f90 FC.f90 FCse.f90 SISSO.f90 -o ~/bin/SISSO

选项 (1) 可实现更佳的准确性和浮点计算的单次运行到多次运行的可重复性；选项 (2) 比 (1) 快约 2 倍，但不同类型的处理器（例如 Intel 和 AMD）之间可能会出现微小的一次运行到多次运行的变化。

如果编译时出现与“mpi”相关的错误，请尝试打开文件“var_global.f90”并用“include ‘mpif.h’”替换行“use mpi”。但是，强烈建议使用“use mpi”（请参见 https://www.mpi-forum.org/docs/mpi-3.1/mpi31-report/node411.htm）。

该程序的模块：
var_global.f90 ！声明全局变量。
libsisso.f90 ！SISSO 中用于数学运算的子程序和函数。
DI.f90 ！模型稀疏化（描述符标识）。
FC.f90 ！将存储在内存中作为数值数据的特征进行特征构造。
FCse.f90 ！将存储在内存中作为表达式树的特征进行特征构建
SISSO.f90 ！主程序。

2.3 运行SISSO

要运行 SISSO 来进行回归、分类和多任务学习，两个文件，SISSO.in 和 train.dat，在工作目录中必须存在。它们的输入模板可以在文件夹 “input_templates” 中找到。简单的命令行用法：

SISSO > log

您可能需要通过运行命令“ulimit -s unlimited”来取消资源限制。
但是，建议用户在集群上运行代码，因为 SISSO 作业通常需要大量 CPU 内核和大量内存。例如，将此命令置于您的提交脚本中：

mpirun -np number_of_cores SISSO > log

VS-SISSO 算法支持对大量输入特征执行符号回归，程序 VarSelect_SISSO.py（位于“utilities”文件夹中）被创建来完成此工作。在 VarSelect_SISSO.py 中有由用户定义的计算参数。在工作目录中有已准备的 SISSO.in、train.dat 和 VarSelect_SISSO.py 的情况下，用户可以通过提交脚本中的命令来运行 VS-SISSO：

Python3 VarSelect_SISSO.py

2.4 输入与输出文件

2.4.1 输入文件

SISSO.in

文件 SISSO.in 中包含了用于定义 SISSO 任务的大多数控制参数，关于机器学习（回归、分类和多任务学习）类型、模型复杂度的信息，数据和输入特征的数量，稀疏化算符等等。该代码逐行读取文件 SISSO.in 以检测关键字和相应的值。以感叹号 (!) 开头的空白行、空格和注释将被忽略。用户可以在文件夹 “input_templates” 中找到SISSO.in 的模板。关于如何设置控制参数的说明在第 5 节中进行了描述。

train.dat

在“input_templates”文件夹中包含了几个 train.dat 文件的模板，可根据特定的机器学习任务使用，包括 train.dat_regression、train.dat_classification、train.dat_regression_multitask 和train.dat_classification_multitask。也就是说，用户可根据 SISSO.in 中定义的 SISSO 任务使用相应的 train.dat 模板来准备他们的训练数据，然后在运行 SISSO 之前将该文件重命名为 train.dat。

在 train.dat_regression 中，第一行为解释/命名每一列的字符串。从第二行开始，每一行描述一个数据，第一个条目为解释该数据的字符串，例如化合物名称，第二个条目为该化合物的目标属性值，所有其余条目均为输入特征的值。每个特征名称都是一个最大长度为 30 个字符的字符串。不允许空格、逗号、制表符和数学运算符等特殊字符出现在特征名称中。

在 train.dat_classification 中，格式与 train.dat_regression 中的格式相同，只不过去掉了目标属性的第二列。分类不需要属性值。SISSO 只需要数据集的类别信息，这将在 SISSO.in 中提供。例如，假设有三个类别 A、B、C，它们分别有 50、60、70 个数据，则需要在 SISSO.in 中设置“nsample = (50, 60, 70)”。相应地，需要在 train.data 中按顺序排列这三组数据，即 A-50 数据从第 2 行到第 51 行，B-60 数据从第 52 行到第 111 行，C-70 数据从第 112 行到第 181 行。

在train.dat_regression_multitask和train.dat_classification_multitask中，这些格式分别与train.dat_regression和train.dat_classification的格式相同，除了数据以任务为单位按顺序排列之外。例如，假设在回归中有三个任务T1、T2、T3，它们分别有100、200、300个数据，那么您将在SISSO.in中设置“nsample = 100, 200, 300”。相应地，您需要将T1-100数据从第2行安排到第101行，将T2-200数据从第102行安排到第301行，将T3-300数据从第302行安排到第601行。在多任务分类中，假设有T1和T2两个任务，其中T1的数据为(50, 60, 70)，T2的数据为(80, 90, 70)，那么您将在SISSO.in中设置“nsample = (50, 60, 70), (80, 90, 70)”，并将T1-(50, 60, 70)数据从第2行安排到第181行，并将T2-(80, 90, 70)数据从第182行安排到第421行。

2.4.2 输出文件

SISSO.out

控制参数、特征空间以及已识别最佳模型的信息输出到文件 SISSO.out。对于 SISSO-nD 计算（在 SISSO.in 中将模型维度设置为 n），将生成从 1 到 n 维度的模型。然而，仅生成 mD 模型（1≤m\u003cn），因为它们对于生成 nD 模型是必需的。这些来自 SISSO-nD 计算的 mD 模型往往不是最理想的，这与直接针对 mD 进行计算时不同。原因在于，对于 SISSO-mD 作业，SIS 子空间的计算可承受大小可以远大于 SISSO-nD 作业，因此获得更好模型的可能性也更高。因此，如果您想要获得 xD 模型，且 SIS 子空间大小适当，建议执行 SISSO-xD 计算。有可能当 SIS 子空间大小超过某个阈值时，模型和/或精度不再随大小而变化，从而达到收敛。

为避免混淆，从 3.5 版开始，在 nD 计算中，SISSO.out 中仅输出 nD 模型。

在回归中，nD模型的功能形式为 y=∑ 1n 𝛽i𝑑i+𝛽0，其中 {𝑑i} 称为nD描述符，{𝛽i} 是相应的系数，𝛽0 是截距。在SISSO.out中，每个𝑑i之后的 features_ID 是 SIS_subspaces 文件夹中的 Uspace.expressions 文件中此特征的行号。在多任务学习回归中，所有任务共享相同的nD描述符，但每个任务都有一组不同的系数。

在分类中，SISSO.out 中只有描述符，没有系数。nD 描述符表示 n 维映射。目前，最多已针对分类实施 3D，但 3D 的凸包构建不稳定，并且仍存在一些错误需要修复。如第 6.2 节所述，分类错误定义为所有重叠区域中的数据量。重叠大小用作模型选择中的第二个度量。如果没有所有域之间的任何重叠，则重叠大小将以负数给出，其大小表示域之间的最小分离距离。使用已识别的 nD 描述符，可以直接应用支持向量分类来查找分隔超平面，或决策树来查找分隔规则。

Models

此文件夹包含 n 个顶级模型/描述符的列表，其中 n 可在 SISSO.in 中设置。每个模型中的特征由其 \spaces 文件夹中 Uspace.expressions 文件中的行号。对于回归，还以 (c_i,i=0,n)_j,j=1,ntask 的格式在 topn_Dm_coeff 文件中提供系数，其中 c_0 是截距，c_1 到 c_n 是线性模型 y=∑1n 𝑐i𝑑i+𝑐0 中与 d_1 到 d_n 对应的系数。

在“模型”内，还包含了文件夹“data_top1”（v3.5 之前的版本中称为“desc_dat”）。它包含了列表中第一个模型的数据文件（这也是 SISSO.out 中显示的模型）。给出了 y_true（来自训练数据）、y_pred（来自模型预测）的值，以及每个训练样本的描述符。文件“desc_Dn_tm.dat”表示维度 n 和任务 m 处的数据。

SIS_subspaces

在 SISSO-nD 计算中，共有 n 个 SIS 子空间。所有这些子空间的并集，包括特征表达式和特征数据，都保存在 SIS_subspaces 文件夹中。Uspace.expressions 文件包含特征表达式，并且这些表达式的相应数据提供了 Uspace_tn.dat 文件中，其中“tn”表示任务 n。例如，Uspace.expressions 中从第 1 行到第 N 行的特征在 Uspace_t001.dat 中从第 2 列到第 N+1 列（对于回归，目标属性占据第一列）有其相应的对应数据，对于分类，则从 1 到 N。Uspace_t001.dat 中的“t001”表明这些数据对应于任务 1。Uspace.expressions 中的特征根据 SIS 评分进行排序，即每个特征与目标属性之间的相关性。

convexnD_hull

在分类中，每个类别域都由数据组的凸包逼近，convex𝑛D_hull文件包含构建数据组在已识别的𝑛D描述符构成的空间中的𝑛D体积信息。

log

日志文件（SISSO > 日志）具有 SISSO 作业的进度信息。可以使用它来监控进度并估算作业的剩余时间

VS_results

执行 VS-SISSO 的结果会输出到文件 VS_results，其中包含有关模型误差变更、选定变量和其他详细信息的信息以及迭代。

2.5 控制参数

2.5.1 SISSO.in

在 SISSO.in 中的关键词可以被用于回归或分类，或两者。那些只适用于回归任务的被标示为 ®，而那些只适用于回归任务的被标记为 ©。(R&C) 表示关键词可适用于回归和分类任务。

ptype

ptype = 1 或 2。目标属性的类型。1：对于连续属性进行回归，2：对于分类属性进行分类。

ntask

(R&C)。ntask=整数型。任务数。ntask = 1：普通的单任务机器学习；ntask>1：多任务学习 (MTL)。

task_seighting

®。task_weighting = 1（默认值）或 2。MTL 回归的权重。1：无权重（所有任务都同等对待，而不考虑其相应数据集的大小）；2：每个任务的权重是各个任务中 $N^t/\sum_t N^t$ 。即此项中的数据量除以全部数据量。

scmt

®。scmt = 真或假。如果为真，则调用符号受限 MTL 回归。

desc_dim

（R 和 C）。desc_dim = 整数。描述符/模型的维度。

nsample

® nsample = 整数，整数，…

train.dat 中示例的数量。

对于单任务回归，输入只是一个整数。对于 MTL 回归，输入是逗号分隔的多个整数，用于表示各个任务的样本数量。
对于单任务分类，会仅有一个圆括号，括号内是逗号分隔的整数，用于表示相应类别的数量。对于 MTL 分类，会有多个以逗号分隔的圆括号。每个括号定义了一个任务，并且每个括号中的整数定义了类别的数量。

restart

(R&C)。restart = 0 或 1。0：从头开始作业；1：继续作业（未完成作业的进度信息存储在文件 CONTINUE 中）

fstore

（R&C）。fstore = 1 或 2。1：特征按数据（一个特征，一列数据）存储在内存中；2：特征以 S 表达式树存储在内存中。只要需要，便可以通过对表达树调用一个评估函数来动态计算特征值。因此，fstore=1 速度快、内存需求高，而 fstore=2 内存需求低，但速度可能比 fstore=1 慢几倍。下图显示了不同特征存储内存方案之间的效率比较。如果因为数据集太大而出现内存问题，则设置 fstore=2；否则，为了获得较高的计算速度，建议使用 fstore=1。

(R&C).nsf=整数。train.dat文件中提供的标量特征数。“标量特征”指每种特征在数据集里是一列。

ops

(R&C) 中，ops = ‘(operator)(operator)(operator)…’。用于特征构建的数学运算符。用户可以从列表中自定义运算符：{+, -, *, /, exp, exp-, ^-1, ^2, ^3, sqrt, cbrt, log, |-|, scd, ^6, sin, cos}，“exp-”、“cbrt”、“|-|”、“scd”分别表示 exp(-x)、三次根号、|x1-x2|和标准柯西分布。

在特征构造期间，存在不同的阶数（表达式树的最大深度），以生成特征空间的不同层级。rung=0（树深度 0）对应于 fcomplexity = 0（输入特征的原始空间）；rung=1（树深度 1）对应于 fcomplexity = 1（最大特征复杂度，即特征中的运算符数量，为 1）；rung=2（树深度 2）对应于 fcomplexity = 2 或 3；rung=3（树深度 3）对应于 fcomplexity = 4 或 5 或 6 或 7。较高的阶数需要巨大的内存，可以使用 fstore=2。

在特征构造中，每个递归中运算符可以被设置为相同或不同。例如，如果 rung =3（fcomplexity 设置为 4 或 5 或 6 或 7），则 ops = ‘(+)(-)(*)(/)’表示 +, -, *, / 运算符将在所有递归的特征构造中都使用；ops = ‘(+)(-)(*)(/)’，‘(+)(exp)’，‘(/)(^2)’ 表示 +, -, *, / 运算符将用于第一个递归，+ 和 exp 将用于第二个递归，/ 和 ^2 将用于第三个递归。

fcomplexity

(R&C)。fcomplexity=整数。特征复杂性定义为特征中的操作符数量。通过输入变量之间数学组合生成的所有输入变量和表达式都被称为特征。因此，fcomplexity=0意味着特征空间中只有输入变量，而fcomplexity=3意味着特征空间中所有特征的复杂性都不大于3。fcomplexity\u003e7对应于树深度\u003e3，需要大量内存。在这种情况下，fstore=2可用于减少内存需求。

funit

(R&C). funit=(整数：整数)(整数：整数) … 每个括号 ( ) 表示一种类型的单位，括号中的整数表示该文件中的 n1 到 n2 之间特征具有相同的单位。例如，如果在 train.dat 中第一个至第五个特征具有单位为电子伏特 (eV) 的能量维度，第六个至第九个特征具有单位为埃的长度维度，第十个特征无量纲，第 11 个特征具有其他不同的单位，则可以在 SISSO.in 中设置 funit = (1:5)(6:9)(10:10)(11:11)。未包含在 funit 中的特征被视为无量纲。如果所有输入特征都是无量纲的，只需注释掉该行或删除关键词 funit。该关键词的作用是使 SISSO 能够执行量纲分析，并避免“能量 + 长度”和“能量 - 长度”之类的无意义组合。

在 SISSO 中，每个特征的单位由一个 jD 向量表示，其中 j 是单位类型的数量（funit 中输入括号的数量）。例如，如果 train.dat 中有 6 个特征，其单位分别为 eV、eV、eV、Å、Å 和 Å，这对应于 SISSO.in 中的 funit = (1:3)(4:6)，则特征单位表示为：

feature1: 1 0
feature2: 1 0
feature3: 1 0
feature4: 0 1
feature5: 0 1
feature6: 0 1

在特征构造过程中，仅允许具有相同单位向量的两个特征之间的线性组合。通过对现有特征进行数学运算而生成的特征的单位向量 (v) 定义为：

feature1 + feature2 和 feature1 - feature2 -> new_feature: v = v1（或 v2）
feature1*feature2 -> new_feature: v = v1 + v2
feature1/feature2 -> new_feature: v = v1 - v2

log(特征1), exp(特征1), sin(特征1), cos(特征1) -> 新特征: v = 0
(特征1)n -> v = v1*n

SISSO 还接受来自名为 “feature_units” 的文件中输入的单位。例如，如果 6 个输入特征的单位为 eV、eV / 、eV、Å、Å、Å，那么用户可以在工作目录中准备包含以下内容的文件：
1.0 0.0
0.5 0.0
1.0 0.0
0.0 2.0
0.0 1.0
0.0 1.0

请注意，当使用“feature_units”文件时，用户必须在SISSO.in中定义funit参数，以便括号参数的数量等于“feature_units”文件中的列数。

fmax_min

(R&C). fmax_min = 浮点值。这个阈值表示：如果一个特征数据中最大绝对值的绝对值小于 fmax_min，则意味着特征内所有数字的更小，将视为零特征并被去除。

fmax_max

（R&C）。 fmax_max = 浮点。如果特征中数据的最大绝对值大于 fmax_max，则表示该特征中某些数字的幅度非常大，该特征将被视为无穷大特征并被丢弃。

nf_sis

(R&C). nf_sis = 整型（或整型_1，整型_2，…，整型_desc_dim）。SIS 子空间的大小。对于 SISSO-nD 计算，将有 n 个 SIS 子空间。如果输入只是一个整数，则 n 个子空间的大小相同。否则，用户可以输入 n 个整数，用逗号分隔，以指定不同子空间的大小。

method_so

（R&C 为“L0”，仅 R 为“L1L0”）。method_so 为“L0”或“L1L0”。“L0”为 l-范数最小化稀疏化方法。“L1L0”是 Ghiringhelli 等人提出的方法 [PRL 114, 105503 (2015)]，即首先通过 LASSO 筛选确定数量的特征（如 30），然后使用 l 找出最佳 1D、2D、3D、… 模型/描述符。“L0”因其在低维模型上的精度而优于“L1L0”。

n|1|0

®. nl1l0= integer. LASSO 之后选定的特性的数量ℓ0 。 SISSO.in 模板中未明确包含用于 LASSO 的其他控制参数。这些参数的默认值是：

L1_max_iter =1e6 LASSO 中最大迭代次数。

L1_dens=120
惩罚参数λ的采样密度。𝜆的取值通过以下表达式确定：
log10(𝜆) = log10(𝜆_max) – (i-1)*[log10(𝜆_max) - log10(𝜆_min)]/L1_dens。这里，𝜆_min是所有系数为零的最小𝜆 [J. Friedman, T. Hastie, R. Tibshirani, J. Stat. Softw. 33, (2010)]，而𝜆_max定义为0.001∗𝜆_min。

L1_tole =1e-6
给定𝜆时，LASSO优化的容差。

L1_minrmse=1e-3
LASSO在新的𝜆试验中停止的条件是，若模型LASSO_RMSE < L1_minrmse。

L1_warm_start = .true.
使用上一个𝜆的解来启动下一个新𝜆的优化。

fit_intercept

（R）。fit_intercept = .true. 或 .false. 为线性模型拟合非零/零截距。

metric

（R）。metric = ‘RMSE’ 或 ‘MaxAE’。回归模型选择所用的指标：RMSE（均方根误差）或 MaxAE（最大绝对误差）。

nmodels

(R&C)。nmodels = 整数。输出排名前几的模型数量（参见“模型”文件夹）。

isconvex

bwidth

（C）. bwidth = 浮点数。每个域的边界容差，用于包含那些位于域外但非常接近域内的数据。

2.5.2 VarSelect_SISSO.py

如果使用 VarSelect_SISSO.py 来完成 VS-SISSO，则可以在 VarSelect_SISSO.py 中定义附加控制参数。有关 VS-SISSO 算法的详细信息，请参阅第 6.3 节。

n_init = integer # Initial size for the S.
n_RS = integer # The size of the S .
n_max = integer # Maximal size of the S (S=S +S ).
nstep_max = integer #maximal iterations of the VS-SISSO
nstep_converge = integer #stop if the best model unchanged for some num. of steps.
restart = 0 or 1 #0 : start from scratch, 1: continue the job
runSISSO = ‘command-to-run-SISSO’ # E.g. ‘mpirun -np 64 SISSO > log’.

2.6 算法 ---- 不翻译

2.6.1 回归

2.6.2 分类

2.6.3 变量选择辅助符号回归

2.6.4 多任务学习

2.6.5 有约束多任务学习

2.7 实用程序

有关 SISSO 数据处理和结果分析的有用工具位于“utilities”文件夹中。可在自述文件中找到有关如何使用这些工具的说明。

文件夹中的Readme：

注意：这些工具可能无法及时随 SISSO 代码更新，因此可能会发现由于 SISSO 输入和输出文件中的格式和参数演进，其中一些工具无法正常工作。

工具的使用（另请参见 SISSO_guide.pdf）：

- SISSO_predict.f90: 从SISSO.out中读取模型，并对未知/测试材料进行预测。
1. Installation: ifort SISSO_predict.f90 -o ~/bin SISSO_predict
2. Use: 在你的工作目录下运行命令 'SISSO_predict'
3. Input files: predict.dat（包含测试材料和特征，格式与train.dat相同）; 
                输入参数可以通过交互式提供，也可以从一个名为SISSO_predict_para的文件中获取。
                (see the provided template).
4. Output files: predict_X.out（predict.dat中材料的描述符坐标）; 
                 predict_Y.out (预测材料在predict.dat中的属性).
5. Note: 避免使用运算符作为变量名；否则，从字符串到公式的转换将会失败。


- SVC.py: 找到基于SISSO描述符的SVC分类超平面在包含文件train.dat的工作目录中运行'python SVC.py >out'。请仔细检查并根据您的应用需求对SVC.py进行必要的修改。

- k-fold-cv & leave-percent-out-cv: 为交叉验证创建子数据集，使用SISSO方法
1. 在文件 k-fold-cv.f90 中定义输入参数 (leave-percent-out-cv.f90)
2. 编译程序: ifort k-fold-cv -o kfoldcv (ifort leave-percent-out-cv.f90 -o lpocv)
3. 在工作目录中，使用全部数据的 train.dat 和相应的 SISSO.in，运行 kfoldcv（或 lpocv）来生成 k 个文件夹。
4. 仔细检查并提交这k个SISSO任务
5. 收集并分析来自k个文件夹的cv结果

- Ionic_Radii
The table of extended Shannon radii (R. Ouyang, Chem. Mater. 32, 595 (2020).

- af2traindat.f90: 自动根据提供的化合物的原子特征创建train.dat文件
1. 准备两个输入文件：atom_features（原子特征列表）和samplelist（训练样本列表）。
   查看提供的模板
2. Compilation: ifort af2traindat.f90 -o af2traindat
3. Use: ./af2traindat

- VarSelect.py: 变量选择辅助的SISSO。
请查看VarSelect.py程序中的使用说明。

3. 使用案例

3.1 单原子催化剂 N原子的吸附能

找出线性模型。

3.1.1 数据

单原子	N吸附能	Nd（d电子数）	EM(原子质量)	EX(原子电负性)
Sc	-4.14856	1	44.96	1.36
Ti	-6.28459	2	47.87	1.54
V	-8.06858	3	50.94	1.63
Cr	-8.46626	5	52	1.66
Mn	-8.50702	5	54.94	1.55
Fe	-7.54049	6	55.85	1.83
Co	-6.14985	7	58.93	1.88
Ni	-5.20835	8	58.69	1.91

3.1.2 SISSO的输入文件

SISSO.in ---- 构建回归的模型（简单的案例）
注意当train.dat 改变时，SISSO的参数也要检查一下对应的数据进行修改。

ptype=1                !1：对于连续属性进行回归，2：对于分类属性进行分类。
ntask=1                !(R&C) ntask = 1：普通的单任务机器学习；ntask>1：多任务学习 (MTL)。
scmt=.false.           !(R) Sign-Constrained MTL is invoked if .true.
desc_dim=1            !(R&C) Dimension of the descriptor/model.
nsample=8              
!nsample=(n1,n2,...)   
restart=0              

fstore=1              
nsf= 3                 ! 3个特征
ops='(+)(-)(*)(/)(exp)(^-1)(^2)(^3)(sqrt)(log)(|-|)'    
fcomplexity=3         
!funit=(1:2)(5:5)    ! 注释掉   
fmax_min=1e-3        
fmax_max=1e5        
nf_sis=20         !(R&C) Number of features in each of the SIS-selected subspace. 

nl1l0= 1
method_so= 'L0'        !(R&C) 'L0' or 'L1L0'(LASSO+L0). 'L0' is always recommended.
fit_intercept=.true.  
metric= 'RMSE'         
nmodel=100            
!isconvex=(1,1,...)    
bwidth=0.001           !(C) Boundary tolerance for classification

train.dat

Atom	Enads	Nd	EM	EX
Sc	-4.14855663	1	44.96	1.36
Ti	-6.28458782	2	47.87	1.54
V	-8.06857765	3	50.94	1.63
Cr	-8.46626119	5	52	1.66
Mn	-8.50701896	5	54.94	1.55
Fe	-7.54048986	6	55.85	1.83
Co	-6.14984802	7	58.93	1.88
Ni	-5.20835431	8	58.69	1.91

执行命令：

~/app/SISSO/bin/SISSO > log

3.1.3 SISSO输出

SISSO.out

****************************************************************
  Sure Independence Screening and Sparsifying Operator (SISSO)  
             Version SISSO.3.5, August, 2024.                
****************************************************************

Read in data from train.dat
Standard Deviation (sqrt[sum(y-y_mean)^2/n]) of the target property:   1.50415
Read in data from SISSO.in
Property type:     1
Number of tasks:        1
Descriptor dimension:        1
Number of samples for the task(s):      8
Restarts :  0
Scheme for feature storing in memory:        1
Number of scalar features:        3
Tier of the feature space (max depth of expression tree):        2
Maximal feature complexity (number of operators in a feature):        3
Unit of input primary feature, each represented by a row vector: 
    0.00
    0.00
    0.00
The feature will be discarded if the minimum of the maximal abs. value in it <    0.10000E-02
The faature will be discarded if the maximum of the maximal abs. value in it >     0.10000E+06
Size of the SIS-selected (single) subspace :     20
Operators for feature construction: (+)(-)(*)(/)(exp)(^-1)(^2)(^3)(sqrt)(log)(|-|) (+)(-)(*)(/)(exp)(^-1)(^2)(^3)(sqrt)(log)(|-|) 
Method for sparse regression:  L0        
Fitting intercept:      T
Metric for model selection: RMSE
Number of the top-ranked models to output:      100
--------------------------------------------------------------------------------

Dimension:   1
-------------------
Feature Construction (FC) starts ...
Total number of features in the space Phi00:              3
Total number of features in the space Phi01:             40
Total number of features in the space Phi02:           4782
Size of the SIS-selected subspace from Phi02:         20
Time (second) used for this FC:            0.01

Descriptor Identification (DI) starts ...
Total number of SIS-selected features for this DI:         20
 
  1D descriptor/model(y=sum(ci*di)+c0) : 
================================================================================
     d001 = abs((Nd)^3-(Nd*EM))    feature_ID:000001
 
        coeff.(ci):    -0.3750038240E-01
                c0:    -0.3171784216E+01
    RMSE and MaxAE:     0.3518172639E+00    0.6717443960E+00
================================================================================
Time (second) used for this DI:            0.00
Total time (second):            0.01
Have a nice day !

SIS_subspaces/Uspace.expressions文件

abs((Nd)^3-(Nd*EM))     SIS_score =      0.9723
abs(log(Nd)-sqrt(EX))     SIS_score =      0.9229
abs((Nd-EX)-(EX)^2)     SIS_score =      0.8949
(sqrt(Nd)-log(Nd))     SIS_score =      0.8911
((Nd)^-1+EX)     SIS_score =      0.8868
abs((EX/Nd)-log(EX))     SIS_score =      0.8751
abs((Nd*EX)-sqrt(EM))     SIS_score =      0.8720
(log(EX)/(EX)^2)     SIS_score =      0.8635
((EX/Nd)+sqrt(EM))     SIS_score =      0.8560
((EX)^2/log(EX))     SIS_score =      0.8547
abs(log(EM)-Nd)     SIS_score =      0.8511
abs((Nd)^2-(EM/Nd))     SIS_score =      0.8485
((Nd+EM)+(EM/Nd))     SIS_score =      0.8373
((Nd)^3-(Nd*EM))     SIS_score =      0.8351
((EX/Nd)+EX)     SIS_score =      0.8163
abs((Nd)^2-(EM/EX))     SIS_score =      0.8158
(log(Nd)/(Nd*EX))     SIS_score =      0.7933
abs((Nd-EX)-(Nd/EX))     SIS_score =      0.7914
(log(Nd)/(Nd+EX))     SIS_score =      0.7855
((EM/Nd)+(Nd*EX))     SIS_score =      0.7789

3.1.4 SISSO得到的描述符

直接从 SISSO.out 中得到x是什么与c1、c0的数值。
$y = - 0.03750038240 * x - 3.171784216$
$x=|(Nd)^3 - Nd *EM|$

3.1.5 展示效果

使用 chatgpt 构建python代码：
在jupyter notebook执行下方代码

# 库文件没有的需要自己pip下载

import numpy as np
import matplotlib.pyplot as plt
from sklearn.metrics import r2_score

# 数据
atoms = ["Sc", "Ti", "V", "Cr", "Mn", "Fe", "Co", "Ni"]
Enads = [-4.14855663, -6.28458782, -8.06857765, -8.46626119, -8.50701896, -7.54048986, -6.14984802, -5.20835431]
Nd = [1, 2, 3, 5, 5, 6, 7, 8]
EM = [44.96, 47.87, 50.94, 52, 54.94, 55.85, 58.93, 58.69]
EX = [1.36, 1.54, 1.63, 1.66, 1.55, 1.83, 1.88, 1.91] 

X = np.abs(np.array(Nd)**3 - np.array(Nd) * np.array(EM))
Y = np.array(Enads)
Y_pred = -0.03750038240 * X - 3.171784216

r2 = r2_score(Y, Y_pred)

plt.scatter(X, Y, label="Data points")
plt.plot(X, Y_pred, color='red', label="y = -0.03750038240 * x - 3.171784216")
plt.text(0.5, -5, f"R² = {r2:.4f}", fontsize=12, color='blue')
# 散点标签
for i, atom in enumerate(atoms):
    plt.annotate(atom, (X[i], Y[i]), textcoords="offset points", xytext=(0,5), ha='center', fontsize=9)

plt.xlabel("X = abs((Nd)^3 - (Nd * EM))")
plt.ylabel("Enads")
plt.legend()
plt.show()

很多文献中给出的单原子催化剂吸附N原子的吸附能是火山图的关系，且描述符大多数是Nd*EX
测试一下:

确实有火山图的趋势。

3.1.6 其他注意事项

高维的只需要更改参数的desc_dim：

desc_dim=1            !(R&C) Dimension of the descriptor/model.

特征属性中两两无意义的运算，使用funit区分开：

funit=(1:2)(5:5)      !(R&C) Feature unit: (n1:n2), features from n1 to n2 in train.dat have the same unit

样本数量：

nsample=8              !(R) Number of samples in train.dat. Set nsample=N1,N2,... for MTL.

特征数量：

nsf= 3                 !(R&C) Number of scalar features provided in the file train.dat

需要值得关注的就这几个参数，其他功能需要自己去测试使用。