Predicting Properties of Oxide Glasses Using Informed Neural Networks

Gregor Maier1,2,* Jan Hamaekers1,3 Dominik-Sergio Martilotti1 Benedikt Ziebarth4

1Fraunhofer Institute for Algorithms and Scientific Computing SCAI, Schloss Birlinghoven, 53757 Sankt Augustin, Germany.

2Institut für Numerische Simulation, Universität Bonn, Friedrich-Hirzebruch-Allee 7,
53115 Bonn, Germany.

3Fraunhofer Center for Machine Learning, Schloss Birlinghoven, 53757 Sankt Augustin, Germany.

4Schott AG, Hattenbergstrasse 10, 55122 Mainz, Germany.

*Corresponding author,

Many modern-day applications require the development of new materials with specific properties. In particular, the design of new glass compositions is of great industrial interest. Current machine learning methods for learning the composition-property relationship of glasses promise to save on expensive trial-and-error approaches. Even though quite large datasets on the composition of glasses and their properties already exist (i.e., with more than 350000 samples), they cover only a very small fraction of the space of all possible glass compositions. This limits the applicability of purely data-driven models for property prediction purposes and necessitates the development of models with high extrapolation power.

In this paper, we propose a neural network model which incorporates prior scientific and expert knowledge in its learning pipeline. This informed learning approach leads to an improved extrapolation power compared to blind (uninformed) neural network models. To demonstrate this, we train our models to predict three different material properties, that is, the glass transition temperature, the Young's modulus (at room temperature), and the shear modulus of binary oxide glasses which do not contain sodium. As representatives for conventional blind neural network approaches we use five different feed-forward neural networks of varying widths and depths.

For each property, we set up model ensembles of multiple trained models and show that, on average, our proposed informed model performs better in extrapolating the three properties of previously unseen sodium borate glass samples than all five conventional blind models.

§ INTRODUCTION

The development of new materials is essential for the modern-day progress in engineering applications and future-oriented technologies. Aside from ever new demands on physical and chemical materials properties, ecological issues, such as sustainability, long service life, environmental compatibility, and recyclability, are of great importance for product development in a variety of different fields. However, the common materials design process is still majorly based on the application of suitable empirical models, on past experiences and educated guesses, and on an extensive subsequent testing phase. The development of new glassy materials, in particular, would benefit to a large extent from a more resource-efficient, systematic, data-driven approach in contrast to the Edisonian trial-and-error approach which is still often used in traditional research and development <cit.>.

The space of all possible glass compositions is very large as a glass can be made from the combination of 80 chemical elements, which leads to e52 possible glass compositions <cit.>.
Moreover, since the influencing parameters are usually known only qualitatively or not at all, the optimization of glass material properties is inherently challenging. A trial-and-error approach to find a glass composition with specific properties for a certain application is time-consuming and often not feasible in practice. An expert-guided approach with integrating experiences from the past is usually not sufficient as well since there are interesting glass properties that are extremely difficult to predict. Especially when properties show nonlinearities, caused, for example, by the so-called borate anomaly in alkali borate glasses <cit.>, conventional exploration and exploitation strategies quickly reach their limits. Therefore, going beyond the area of known materials requires new approaches based on new and innovative methods. The field of machine learning (ML) provides such methods which allow to generate accurate models based on existing data in order to predict the properties of yet unseen materials.

§.§ Related Work

In recent years, ML techniques have been widely used for accelerating materials design <cit.>.
In glass science, there have been several successful attempts to use ML to predict, i.a., optical, physical, and mechanical properties of glasses <cit.>.
Most ML models perform exceptionally well in interpolating the training data.
However, given the high-dimensional search space of all possible glass compositions and its sparse coverage by experimental data, the search for new glass materials is majorly a question of designing models which possess a high extrapolation power. We refer to <cit.> and references therein for reviews of the current status of ML in glass science and future challenges.

To address the lack of extrapolation power, ordinary ML methods can be extended by integrating prior knowledge which exists independently of the learning task. This idea is termed informed machine learning and we refer to the recent survey <cit.> for a taxonomy and thorough overview of its application in current ML state-of-the-art use cases. For glass design,  this idea is utilized, e.g., in <cit.> where the empirical MYEGA formula is integrated into a neural network architecture to predict the viscosity of a glass based on its compound fractions and temperature. In <cit.>, this approach is developed further by additionally integrating prior chemical and physical knowledge of the glasses' elements into the training data. Similarly, in <cit.>, the authors use external chemical and physical knowledge to carefully design enriched descriptors of glass compositions which are used as inputs for ML models to predict properties of oxide glasses. In <cit.> and <cit.>, the authors predict the dissolution kinetics of silicate glasses in an informed manner by suitably splitting the training data and using a descriptor which encodes the glasses' network structure, respectively. They demonstrate the superior performance of the informed approach compared to the uninformed approach. This superiority is also shown in <cit.>, where the authors design a neural network model which is informed by statistical mechanics in order to predict structural properties of oxide glasses.

§.§ Contributions

In this paper, we propose a new ML model based on neural networks for the property prediction of oxide glasses which integrates prior knowledge in order to achieve a high degree of extrapolation of the training data. We modify the ideas from <cit.> in order to predict three material properties, that is, the glass transition temperature \(T_g\), the Young's modulus \(E\) (at room temperature), and the shear modulus \(G\). We focus our analysis on binary oxide glasses, that is, oxide glasses which consist of exactly two compounds. Our model is informed in the sense that we explicitly integrate prior knowledge into the design of our training data, the hypothesis set, and the final hypothesis at four major points in our learning pipeline. We place emphasis on explaining how this is done in detail in terms of the taxonomy in <cit.>. Especially the design of the network architecture to realize permutation invariance with respect to the input features seems, to the best of our knowledge, to be new in the field of glass materials modeling.

To examine the extrapolation power of our models, we train and validate them on glass samples which do not contain sodium in their compositions. The trained models are then used to predict the properties of sodium borate glass compositions with varying compound fractions. For each property, we train multiple models and study the average performance of the model ensemble.
To demonstrate the superiority of the informed model ensemble compared to blind (uninformed) approaches, we perform the same experiments with five standard fully connected feed-forward neural networks of varying widths and depths without integration of any prior knowledge. We compare the results quantitatively in terms of error metrics and qualitatively in terms of a meaningful approximation of the respective composition-property curves.

The remainder of this paper is organized as follows: In Sect. <ref>, we explain our methodology. That is, in Sect. <ref>, we present our automated pipeline for collecting and preparing data for model training, validation, and testing. In Sect. <ref>, we describe the different model setups in the blind and the informed setting. In Sect. <ref>, we explain how we train and evaluate our models. We discuss the results of our experiments in Sect. <ref> and conclude our findings in Sect. <ref>.

For notational convenience, we use the letter P whenever we refer to one of the three properties \(T_g\), \(E\), or \(G\). Moreover, for all entities which exist for every property P, we use the prefix “P-” to specify the respective entity. For example, given P, we refer to the dataset that is used to train a model for predicting P by “P-training set”.

Moreover, we use the symbols \(\mathbb{N}\) and \(\mathbb{R}\) to denote the set of positive integers and the set of real numbers, respectively.

§ METHODOLOGY

The prediction quality of any data-driven machine learning algorithm in the context of supervised learning is strongly dependent on the quantity and quality of the training data. Before presenting our neural network approach to the problem of glass property prediction in detail in Sect. <ref> and Sect. <ref>, we therefore describe in the following Sect. <ref> how we collect and prepare our data.

§.§ Data Collection and Preparation

We use data from the INTERGLAD Ver. 8 database <cit.> and the SciGlass database <cit.> and merge them together into a common  database. For the identification of oxide glasses we follow the same definition as in <cit.> and only consider glasses whose mole atomic fraction of oxygen is at least 0.3 and whose compounds do not contain the chemical elements S, H, C, Pt, Au, F, Cl, N, Br, and I, which could affect the balance of oxygen. The resulting  database of oxide glasses consists of 420973 glass  samples in total. It lists the mole atomic fractions of 118 chemical elements and the mole fractions of 439 compounds, i.e., the oxides that a glass composition consists of, together with the values of 87 material properties. However, among the 118 elements, only 66 elements appear with non-vanishing fraction in at least one glass sample. Among the 439 compounds, only 183 compounds appear with non-vanishing fraction in at least one glass sample.

To obtain clean data for training, validating, and testing our models, we apply a sequence of preprocessing steps which follows in parts the procedure described in <cit.>. For each glass property P, we extract clean data from the “dirty”  database in an automated fashion in form of a preprocessing pipeline whose steps are schematically depicted in Figure <ref>. The number of samples which are dropped in each step is shown in Table <ref>.

We begin with all samples from the entire  database. As a first step, we make sure that all glass samples have numerically valid entries. That is, we first remove glass samples which have a Not-a-Number (NaN) entry for at least one compound fraction. Moreover, we drop all glasses which have NaN entries for P.[At the end of this and all the following preprocessing steps, we always drop all compounds which do not appear in any of the glass samples that are present in the dataset at the respective preprocessing stage.]

Next, we make sure that all glass samples are physically valid binary glass compositions.
For this, we first discard glasses whose compound fractions do not add up to a value in the closed range between \(0.9999\) and \(1.0001\). Then, we exclude all samples which do not consist of exactly two compounds.

To ensure physically valid property values, we fix a closed range of values between a minimum and maximum cut-off value for each property P (see Table <ref>). We determine these values by investigating the distribution of the glass samples with respect to their P-values in the datasets that result from the preprocessing pipeline up to this point. Property values outside of this range are considered non-physical but may be present in the database due to typos or other mistakes. Hence, we drop each glass sample with a P-value outside of the respective range.

As the minimum and maximum values are rather crude bounds, in a further step, we remove glasses with extreme P-values, which have a high chance of still appearing in the datasets again because of typos or other mistakes. To do so, we compute the \(0.05\)th percentile and the \(99.95\)th percentile of P-values among all remaining glass samples and subsequently discard all glasses with P-values below the lower percentile or above the upper percentile.

A lot of glass samples appear in both the INTERGLAD and the SciGlass database. Consequently, there may be many duplicates among the remaining data points at this stage of the preprocessing pipeline. We therefore apply a duplicate filter which consists of the following steps:

* We group all glasses with the same (up to the fifth decimal place) compound fractions.

* For each such group we do the following:

2.1 We drop all but the first sample which agree exactly in their values of P.

2.2 We compute the midpoint of the range of values of P among all remaining samples.

2.3 If the P-value of every sample has a distance to the midpoint smaller than a certain P-dependent threshold (see Table <ref>), then, as a representative of the group of duplicates, we select the first sample in the group, assign to it the median of the P-values of the samples in the group, and drop all other samples.
Otherwise, we discard the whole group of glass samples.

The values for the duplicate thresholds are determined by using domain know-ledge and investigating the average spread of P-values in a group of duplicates.

In the next step, we deal with compounds of low representability and iteratively drop compounds which appear in less than one percent of all remaining glass samples. This allows us to reduce the dimension of the compound space and leaves us only with glasses whose compounds are present in sufficiently many samples in order to use them for robust model training.[At the end of each iteration, we again drop those samples whose compound fractions do not add up to a value in the closed range between \(0.9999\) and \(1.0001\).]

As a final step, we apply an outlier detection based on a one-class support vector machine (SVM) followed by an outlier detection based on Gaussian process (GP) regression.[We fit a Gaussian process to the data and drop samples with a too large deviation in their P-value from the respective mean curve.]

The resulting cleaned datasets encompass all problem-specific information which is available for each glass property P. Given P, we denote the elements and compounds which are present (with non-vanishing value in at least one glass sample) in the corresponding cleaned dataset as P-elements and P-compounds, respectively. The cleaned datasets are subsequently split into training, validation, and test sets as described in Sect. <ref>.

Steps in the preprocessing pipeline as described in Sect. <ref>.

Data reduction in each step of the preprocessing pipeline. The first and last row show the number of samples (#S) and number of compounds (#C) which are present in the dirty and cleaned dataset, respectively. The rows in between show the number of samples and compounds which are dropped in the respective preprocessing steps. We remark that in the cleaned datasets, the number of elements appearing with non-vanishing fraction in at least one glass sample is 32 for \(T_g\), 23 for \(E\), and \(24\) for \(G\) (and therefore coincides with #C).

P     2r\(T_g\)     2r\(E\)     2r\(G\)

1r#S     1r#C     1r#S     1r#C     1r#S     1r#C

Dirty dataset     420973     439     420973     439     420973     439

Drop NaNs     344247     283     396460     329     410589     356

Check compound fractions     17483     1     5482     0     1867     0

Filter binary glasses     50681     83     16651     57     6386     38

Min-max filter     16     0     21     0     15     0

Drop extreme values     10     0     4     0     4     0

Duplicate filter     5902     0     1577     0     1253     2

Drop rare compounds     229     40     79     30     67     19

One-class SVM     22     0     5     0     10     0

GP regression     205     0     62     0     63     0

Cleaned dataset     2178     32     632     23     719     24

Minimum and maximum cut-off values and duplicate thresholds used in the preprocessing pipeline.

P          Min. cut-off        Max. cut-off        Duplicate threshold

\(T_g\)     (°C)     50     1.8e3     5.0

\(E\)     (GPa)     5.0     2.0e2     1.5

\(G\)     (GPa)     0.10     2.0e2     0.75

§.§ Model Setups

We use neural networks for the approximation of the composition-property relationship of binary oxide glasses. The target quantity is given by one of the respective properties P. The composition of a glass can be represented in various ways. Designing a representation in form of a feature vector that encodes a given glass composition in a way that is suitable as input for a neural network is an essential part of the modeling process and is one of the key differences between the blind (uninformed) and our informed learning approach. A second difference lies in the design of the model architecture where the black-box modeling approach of standard blind feed-forward neural networks can be leveraged in the informed setting by integrating prior scientific knowledge.

In the following Sect. <ref> and Sect. <ref>, we describe in detail the choice of the feature vectors and the network architectures for the blind and informed models, respectively, and highlight their differences.

§.§.§ Blind Models

In the blind approach, only the available problem-specific data is used to design a suitable ML model for the composition-property relationship of oxide glasses. This approach is blind or uninformed in the sense that no prior knowledge that exists independently of the learning task is integrated into the model setup.

Feature Vectors

Each glass composition is, by definition, uniquely determined by its compound fractions. It is therefore natural to use the compound fractions, grouped together in a feature vector for a given glass composition, as input for a neural network model to predict one of the glass's properties.

Network Architectures

If no further information is available, the standard architectural design of a neural network is given by a (fully connected) feed-forward neural network (FFNN) <cit.>. This class of models satisfies the universal approximation theorem, that is, for any continuous function on a compact domain there exists a FFNN which approximates the function within a given arbitrary tolerance <cit.>.
This result justifies the usage of the set of FFNNs as hypothesis space. In the context of glass materials research, this approach is followed for example in <cit.> to model several different properties of oxide glasses.

A fully connected FFNN is characterized by (i) its input and output dimensions, i.e., the number of units in its input and output layer, respectively, (ii) its depth, i.e., the number of layers (without counting the input layer), and (iii) the width, i.e., the number of units, of each hidden layer. An example architecture of a fully connected FFNN is shown in Fig. <ref>.

node=[very thick,circle,draw=myblue,minimum size=22,inner sep=0.5,outer sep=0.6]
connect=[->,thick,mydarkblue,shorten >=1]

node 1/.style=node,mydarkgreen,draw=mygreen,fill=mygreen!25,
node 2/.style=node,mydarkblue,draw=myblue,fill=myblue!20,
node 3/.style=node,mydarkred,draw=myred,fill=myred!20,

ı[evaluate==̧int(ı==); =/2-ı;
=; =;
] in 1,...,

[node ] (N-ı) at (,) ;

ȷin 1,...,[]
[white,line width=1.2,shorten >=1] (N-ȷ) – (N-ı);
[connect] (N-ȷ) – (N-ı);

=[connect] (N-ı) –++ (0.5,0);

[connect] (0.5,) – (N-ı);

[above=3,align=center,mydarkgreen] at (N1-1.90) Input
[above=2,align=center,mydarkblue] at (N2-1.90) Hidden
layer 1;
[above=2,align=center,mydarkblue] at (N3-1.90) Hidden
layer 2;
[above=2,align=center,mydarkblue] at (N4-1.90) Hidden
layer 3;
[above=3,align=center,mydarkred] at (N-1.90) Output

Schematic architecture of a fully connected FFNN with input dimension \(3\), output dimension \(2\), depth \(4\), and constant width \(4\). Each layer represents an affine linear function. The additive bias nodes are not shown. In case of the blind models, the input nodes store the compound fractions of a given glass sample and the output node provides the predicted value for P. Image adapted from <cit.>.

We use a variety of different FFNNs as benchmark models which we compare our informed model to. To capture the main architectural trends of designing a FFNN and their effects on the prediction quality, we consider five different FFNNs with depths \(L = 2, 4, 8, 16, 32\) for each property P. The input dimensions are determined by the number of respective P-compounds and the output dimension is always one as we predict scalar-valued properties.
For each P-model, we choose the width to be constant for all hidden layers such that the total number of trainable parameters is roughly the same among all P-models including the informed model which we describe in Sect. <ref>. Each hidden layer is a linear layer with an additive bias term. In accordance to the universal approximation theorem, the output layer is a linear layer with no additive bias term. The exact dimensions of all models are summarized in Table <ref>. In all models, we use the rectified linear unit (ReLU) as activation function.

Model hyperparameters

2*\(T_g\)     5c2*Blind     3cInformed

5c     Former     Non-former     Down

Input dim.     32     32     32     32     32     28     28     32

Length \(L\)     2     4     8     16     32     4     4     4

Width \(W\)     304     63     38     26     18     32     32     32

Output dim.     1     1     1     1     1     14     18     1

#Parameters     10336     10206     10184     10712     10872     3c10336

2*\(E\)     5c2*Blind     3cInformed

5c     Former     Non-former     Down

Input dim.     23     23     23     23     23     26     26     23

Length \(L\)     2     4     8     16     32     4     4     4

Width \(W\)     385     63     38     25     17     32     32     32

Output dim.     1     1     1     1     1     10     13     1

#Parameters     9625     9639     9842     9725     9605     3c9623

2*\(G\)     5c2*Blind     3cInformed

5c     Former     Non-former     Down

Input dim.     24     24     24     24     24     26     26     24

Length \(L\)     2     4     8     16     32     4     4     4

Width \(W\)     372     63     38     25     17     32     32     32

Output dim.     1     1     1     1     1     10     14     1

#Parameters     9672     9702     9880     9750     9622     3c9688

§.§.§ Informed Model

Rather than just using the compound fractions as input features for a neural network, we can increase the informational capacity of a glass sample's representation by utilizing characteristic chemical and physical quantities of each element which is present in the given glass sample and provide them as additional inputs to a neural network. Features which are carefully engineered in such an informed manner can lead to an improved prediction quality of the model, given that the model's expressive power is large enough. The latter issue is a question of the model's architecture. If there are too few parameters, the model will underfit the training data no matter how carefully we designed the input features. If there are too many parameters, however, the model might overfit the training data and pick up on spurious patterns and noise in the input features. In general, by building as much prior information as possible into the model's architecture we expect to obtain a more robust inference behavior, especially in the extrapolation regime.

Feature Vectors

Compared to the uninformed approach, we change our viewpoint and identify a given glass composition not by the fractions of its compounds but by the mole atomic fractions of its elements. For each element, there is additional extensive scientific knowledge about its chemical and physical properties, which exists independently of our learning problem. According to the taxonomy in <cit.> this knowledge is represented as a weighted graph. Its nodes are given by elements and properties and each element node is connected with a property node via an edge which is weighted by the element's respective property value. We integrate this knowledge into our training data by designing feature vectors in a hybrid fashion. We partly follow the approach used in <cit.>, where, for each element, the authors extract physical and chemical properties, i.a., from the Python library  <cit.>. They design and select feature vectors for a neural network model in a way that allows them to complement information from the space of chemical compositions with information from the space of chemical and physical properties.

In our case, we first extract for each glass property P and each P-element a list of characteristic chemical and physical properties from . We drop properties with non-numeric values and only keep those which are available for all P-elements. We also drop properties which we consider to be unrelated to the elements' influence on the glass material properties, namely, the elements' abundances in the earth crust and the elements' dipole polarizability uncertainties.
We refer to <cit.> and references therein for a detailed explanation of all the available properties in the  library.

Among the remaining properties, we drop those which are highly correlated. More specifically, we compute the standard pairwise Pearson correlation coefficients and iteratively, for each pair of properties whose coefficient is larger than \(0.95\), we only keep one of the two properties. The resulting list of properties is given in Table <ref>.

For each glass property P and each P-element, we group the resulting element properties into vectors, which we call P-element-vectors. These element-vectors represent all prior knowledge about the elements which exists independently of the learning problem.

In the next step, we complement this information with our given problem-specific data. For each glass property P and a given glass sample, we consider the collection of P-element-vectors corresponding to the sample's elements and extend each of them by one more component which lists the mole atomic fraction of the respective element in the glass sample. Eventually, we obtain for each glass sample a collection of feature vectors \((v_1,\dots,v_M)\in\mathbb{R}^{d\times M}\), where we arrange the tuple in lexicographical order of the elements' symbols. Here, \(M = M(\Vec{P})\) is the number of P-elements and \(d = d(\Vec{P})\) is the number of properties resulting from the property extraction process described above (including the entry with the mole atomic fraction). Each glass sample can thus be represented as a point in a subspace \(\Omega = \Omega(\Vec{P})\subset \mathbb{R}^{d\times M}\).

Network Architecture

For each glass property P, we want to design a neural network which approximates the functional relationship \(f: \Omega \to \mathbb{R}, \Vec{V}\mapsto \Vec{P}(\Vec{V})\), where \(\Vec{P}(\Vec{V})\) is the value of property P for the glass sample with representation \(\Vec{V} = (v_1,\dots,v_M)\).
We design the architecture of the network by two leading principles in the spirit of informed learning.

First, we observe that the order in which the feature vectors are passed to the function \(f\) actually does not matter. That is, the function \(f\) is permutation invariant with respect to the order of the input vectors. More specifically, \(f\) is a function on sets of the form \(\{v_1,\dots, v_M\}\). In terms of the taxonomy in <cit.> we use this scientific knowledge, which is represented as a spatial invariance, and directly integrate it into the architecture of the network which we use to approximate \(f\). It is shown in <cit.> that such a function \(f\) on sets can be written as

\begin{equation}\label{eq: func on sets}
f(\{v_1, \ldots, v_M\}) =  \psi\left(\sum_{i=1}^M \phi(v_i)\right) \, ,
\end{equation}

where \(\phi : \mathbb{R}^d \to \mathbb{R}^N\) denotes an inner embedding function with \(N\in\mathbb{N}\) being an appropriately chosen embedding dimension, and \(\psi : \mathbb{R}^N \to \mathbb{R}\) denotes an outer (downstream) function. Here, \(\phi\) and \(\psi\) can be approximated by neural networks. Using the universal approximation theorem of neural networks, the right-hand-side of (<ref>) yields architectures of neural networks which, in principal, can approximate \(f\) arbitrarily well.

In our specific use case, we can refine the network's architecture even further by integrating prior chemical knowledge. Glass oxides can be categorized in three groups <cit.>. Glass formers are oxides that can readily form a glassy material and build the backbone of a glass's network structure. Glass modifiers are oxides that cannot form a glassy material by themselves but influence its material properties when mixed with a glass former. Glass intermediates are oxides which can act both as a glass former as well as a glass modifier depending on the respective cation's oxidation number. For our purposes, we only differentiate between oxides which are glass formers and oxides which are not glass formers. We refer to the latter group as glass non-formers.
We use the classification proposed in <cit.> to determine for every element whether its oxide is a glass former or a glass non-former. The classification is shown in Table <ref>.

Classification of elements based on the glass-forming and glass-non-forming properties of their oxides. The last row shows all elements whose oxides are glass formers according to the classification in <cit.>. The second column lists, for given P, the respective P-elements whose oxides are glass formers and which are a subset of the elements in the last row. The fourth column lists all respective P-elements whose oxides are glass non-formers. We classify oxygen as a glass non-former.

P     Formers     #Formers     Non-formers     #Non-formers

\(T_g\)     As, B, Bi, Ge, Mo, P, Pb, Sb, Si, Sn, Te, Tl, V, W     14     Ag, Al, Ba, Ca, Cs, Cu, Fe, Ga, K, La, Li, Mg, Na, O, Rb, Sr, Ti, Zn     18

\(E\)     B, Bi, Ge, Mo, Nb, P, Pb, Si, Te, V     10     Al, Ba, Ca, Co, Cs, K, Li, Mg, Na, O, Sr, Ti, Zn     13

\(G\)     B, Bi, Ge, Mo, Nb, P, Pb, Si, Te, V     10     Al, Ba, Ca, Co, Cs, K, Li, Mg, Na, O, Rb, Sr, Ti, Zn     14

All     As, B, Bi, Ge, Mo, Nb, P, Pb, Sb, Se, Si, Sn, Ta, Te, Tl, V, W      17

The scientific knowledge whether an element's oxide has glass-forming or glass-non-forming ability is naturally represented as a simple knowledge graph, where each element is represented by a node. There is also a glass former node and a glass non-former node. Each element node is connected via an edge with the glass former node or the glass non-former node depending on whether the element's oxide is a glass former or a glass non-former. Due to the largely different influence on a glass's properties, we integrate this prior knowledge additionally into our hypothesis set by using two functions to treat glass formers and non-formers separately. The glass former network receives as input only feature vectors of elements whose oxides are glass formers. The glass non-former network receives as input all other elements whose oxides, by definition, are glass non-formers. The outputs of the glass former network are added together, as are the outputs of the glass non-former network. The results are concatenated and then used as input for the downstream network which yields the final prediction for the respective property P.

More specifically, let \(\Omega = \Omega_f \cup \Omega_{nf}\) be the decomposition of \(\Omega\) into the space \(\Omega_f\) of representations of glass formers and the space \(\Omega_{nf}\) of representations of glass non-formers. We then replace the inner function \(\phi\) in (<ref>) by two separate functions, \(\phi_{f}: \Omega_f \to \mathbb{R}^{N_f}\), \(N_f\in\mathbb{N}\), for the glass former network and \(\phi_{nf}:\Omega_{nf}\to\mathbb{R}^{N_{nf}}\), \(N_{nf}\in\mathbb{N}\), for the glass non-former network. Permutation invariance then holds only within the feature vectors \(v_1,\dots,v_{M_f}\) corresponding to the \(M_f\) glass formers and within the \(M_{nf} := M - M_f\) feature vectors \(v_{M_f + 1}, \dots, v_M\) corresponding to the glass non-formers. The resulting representation of \(f\) then has the following form,

\begin{equation}\label{eq: informed func on sets}
f(\{v_1, \ldots, v_M\}) = \psi\left(\sum_{i=1}^{M_f} \phi_f(v_i), \sum_{i=M_f + 1}^M \phi_{nf}(v_i)\right) \, ,
\end{equation}

where, under slight abuse of notation, we used the same notation \(\psi\) for the downstream function as in (<ref>). As glass former network, glass non-former network, and downstream network we use three separate ReLU-FFNNs to approximate the functions \(\phi_f, \phi_{nf}\), and \(\psi\) in (<ref>), respectively. Their widths and depths are listed in Table <ref>. The overall network architecture of our informed model is illustrated in Fig. <ref>.

The embedding dimensions \(N_f\) and \(N_{nf}\) are hyperparameters of the glass former and non-former network, respectively. It is shown in <cit.> that in the scalar case, \(d=1\), the choice \(N = M\) in (<ref>) is a sufficient and necessary condition in order to approximate the function \(f\) arbitrarily well by a neural network whose architecture is given by the right-hand side in (<ref>). In the vector-valued case, \(d > 1\), to the best of our knowledge, no non-trivial necessary condition on the embedding dimension is known so far. In <cit.>, the authors prove a sufficient condition in form of an upper bound on the embedding dimension, which, however, is very pessimistic. Based on the results in the one-dimensional case, we choose \(N_f = M_f\) and \(N_{nf} = M_{nf}\) in (<ref>).

Architecture of the informed model. The feature vectors of the elements are split according to whether the elements' oxides are glass formers or glass non-formers and input to separate neural networks. The results of the latter are first summed individually and then concatenated to a vector which is used as input for the final downstream neural network that predicts a value for property P.

§.§ Model Training and Evaluation

Recall that we consider three different glass material properties: glass transition temperature \(T_g\), Young's modulus \(E\) at room temperature, and shear modulus \(G\). For each of these properties P, we split the cleaned datasets from Sect. <ref> further into datasets for training, validation, and testing. We then apply the blind models and the informed model discussed in Sect. <ref> and Sect. <ref>.

For data management and visualization, we use the Python libraries  <cit.>,  <cit.>, and  <cit.>, respectively.
The neural network models are built using the  <cit.> module.

We describe the data splitting in a bit more detail. For each property P, we apply the following steps. First, we apply the preprocessing pipeline described in Sect. <ref>. Then, we apply the feature design and selection processes for the blind and the informed models as described in Sect. <ref> and Sect. <ref>, respectively. Next, we split up the resulting cleaned dataset into those glass samples which contain sodium and those which do not contain sodium. From the samples which contain sodium we extract only those binary oxides which consist of B2O3 as glass former and Na2O as glass non-former. The resulting dataset is our P-test set. The other glass samples, which do not contain sodium, are randomly split for each model into a P-training set and a P-validation set using a \(80\%/20\%\) ratio.
The dimensions of the resulting datasets are shown in Table <ref>.
We emphasize that sodium as an element is totally absent in the training and validation sets and only present in the test sets. Examining the performance of the trained models on the test sets therefore allows us to properly evaluate their extrapolation power.

Dimensions of the training, validation, and test sets.

P     #Training samples         #Validation samples        #Test samples

\(T_g\)     1385     347     125

\(E\)     415     104     42

\(G\)     477     120     73

We use the bagging method from the field of ensemble learning <cit.>. For each model setup and each property P, we train 50 models. Their architectures and weight initializations are the same, but each model is trained and validated on a different random \(80\%/20\%\)-split into training set and validation set. We therefore end up with a model ensemble of 50 different models.

For training, we use the ADAM optimizer with default settings <cit.> and weight decay of e-5 and train for a maximum of 1000 epochs with a batch size of 8. We start with a learning rate of 0.001 and multiply it by a factor of 0.5 if the model's performance on the validation set in terms of the mean squared error (MSE) does not improve over the course of 50 epochs. Moreover, to avoid overfitting, we use the early stopping criterion and stop training if the model's MSE on the validation set does not improve over the course of 100 epochs.

After training, we apply a post-processing step. For each property P, we discard those models whose predictions for P on the whole test set can be considered to be constant. More specifically, we first compute for each sample in the test set the mean value of the models' predictions for P. Then, we drop those models where the deviation of the predicted property values for all samples in the test set from the respective mean value is less than or equal to the P-duplicate threshold from Table <ref>. This is in alignment with informed learning since we know a priori that for each property P, not all Na2O-B2O3 glass samples have the same P-value. In terms of the taxonomy in <cit.>, we thus use this expert knowledge, which is represented as algebraic equations, i.e., being of non-constant value, and integrate it into our final hypothesis.

Among the remaining models we compute the mean and the \(95\%\)-confidence interval of the predictions. This yields the final prediction of the model ensemble and quantifies its uncertainty.
We compare the ensembles' performances quantitatively in terms of their root mean squared errors (RMSE), mean absolute errors (MAE), and maximum errors (MAX) on the respective P-test sets, which are summarized in Table <ref>. All blind 32-layer networks yield constant predictions for all three properties and are therefore discarded as non-physical in the post-processing step. Nevertheless, we still record their respective ensemble errors in Table <ref> to get a more conclusive picture. However, when talking about the best and worst error values, we only consider the ensembles of blind models with depths \(L=2,4,8,16\) and neglect the values of the models with 32 layers.

To also get a qualitative picture of the ensembles' extrapolation performances, we plot the composition-property curves of the ensembles' averaged predictions on the P-test sets in Figs. <ref>–<ref>. Recall that the test sets consist of the Na2O-B2O3 glass samples together with their respective P-values. We also plot the predictions of the best and worst performing model of each P-ensemble in terms of the RMSE on the test set. Moreover, we plot the property values of all available alkali borate glasses, that is, binary glasses which consist of B2O3 as glass former and Na2O, Li2O, and Rb2O as glass non-former, respectively. It is known that these glass compositions have similar material properties <cit.>.

Concerning Figs. <ref>–<ref>, a few remarks are in order. First, there are no Rb2O-B2O3 glass samples available for \(E\). Next, since we only consider binary oxide glasses, knowing the compound fraction of B2O3 completely determines the compound fraction of the respective alkali oxide as glass non-former as well. Finally, since the blind \(32\)-layer networks are discarded, their predictions are not shown.

Results and errors of the model ensembles' averaged predictions on the test sets. Among the blind ensembles with depths \(L=2,4,8,16\), bold blue numbers denote the lowest error values, bold red numbers the highest ones. Bold black numbers denote the lowest error values among all model ensembles. The last column shows the relative improvement in the error when comparing the error value of the informed ensemble to the lowest error value (blue) among the blind ensembles. The blind models with 32 layers are not considered in the error analysis as these models only yield constant predictions and are therefore discarded as non-physical.

2*\(T_g\) (°C)     6cBlind

Depth     2     4     8     16     32     -2*  Informed       -2*Rel. improv.

RMSE          86.2     111      112     109     86.5     62.5     27%

MAE          63.0     77.2     81.7      84.4     73.6     44.4     30%

MAX          265      341     338     311     205     186     30%

#Non-const. predictions          50     50     50     48     0     24

2*\(E\) (GPa)     6cBlind
Depth     2     4     8     16     32     -2*  Informed       -2*Rel. improv.

RMSE          9.88     10.7      12.0     11.6     15.6     4.67     53%

MAE          8.58     9.44      10.7     10.3     12.5     3.50     59%

MAX          16.8     18.2     19.5      19.8     33.7     12.5     26%

#Non-const. predictions         50     50     50     43     0     48

2*\(G\) (GPa)     6cBlind
Depth     2     4     8     16     32     -2*  Informed       -2*Rel. improv.

RMSE          2.92     3.65      4.02     3.81     5.94     1.39     52%

MAE          2.43     3.09      3.39     3.13     5.08     1.12     54%

MAX          6.10     7.24     7.89      8.37     12.2     2.95     52%

#Non-const. predictions         50     50     50     39     0     48

§ RESULTS AND DISCUSSION

We first compare the models' performances quantitatively in terms of their average errors in Table <ref>. Among the blind models, we note that the ensembles of shallow two-layer networks perform best for all three properties in terms of RMSE, MAE, and MAX. Considering the worst performing ensembles, we note that the deeper networks with depths of 8 and 16 layers perform worst, on average, in terms of almost all three error metrics for all three properties. Only for \(T_g\) in the case of MAX, the ensemble of models with only four layers performs worst. We note that in this case, the ensemble of models with a depth of 32 layers actually performs best in terms of MAX. In all other cases, however, the 32-layer network ensembles perform worst for all three properties when compared to the other blind models. We conclude that, in general, increasing network complexity in terms of increasing depth tends to lead to worse performing models.

The number of models with non-constant predictions clearly decays with increasing network depth for all three properties. Whereas the networks with depths of \(2, 4\), and \(8\) layers lead to no constant predictions, the \(16\)-layer networks lead to some constant predictions. There is a steep decay when increasing the number of layers from \(16\) to \(32\), where all models for all properties lead to only constant predictions. A possible explanation for this phenomenon could be that the models' loss landscapes become more and more rugged with increasing network depth yielding constant predictions to be local minima which are hard to escape during the optimization routine. This matches the observation from above that increasing network depth generally tends to lead to worse performing models.

Invoking now the errors of the informed models, we see that they perform best, on average, in terms of all three error metrics for all three properties. They lead to a relative improvement in the errors between 26% up to 59%. For \(E\) and \(G\), only two models yield constant predictions, whereas for \(T_g\) more than half of all models do. Again, this could indicate that the loss landscape of the informed networks is much more rugged for \(T_g\) than for the other two properties.

To get a more conclusive qualitative picture of the extrapolation behavior of all models, we take a closer look at Figs. <ref>–<ref>. We observe that the blind networks in terms of the ensembles' means as well as the best and worst performing models are not able to qualitatively capture the trend of the Na2O-B2O3 curves correctly and instead generally deviate from the test points to a large extent. However, the ensembles' predictions for all blind networks seem to be quite close to each other for all three properties. This is reflected by the small width of the confidence band around the mean curves as well as the similar shape of the mean curves and the curves of the best and worst performing models. This indicates that the blind models are robust with respect to training.

The mean curves of the informed model ensembles qualitatively capture the trend of the Na2O-B2O3 curves to a more acceptable degree. This is most noticeable in the cases of \(T_g\) and \(G\) where the mean curves are able to capture the nonlinearity of the respective Na2O-B2O3 curve, which the blind networks are not capable of. For \(E\), the informed model ensemble yields more accurate trajectories than the blind ensembles in the linear regime with B2O3-fractions between 0.7 and 1.0, but the kink in the Na2O-B2O3 curve at a B2O3-fraction of around 0.7 is not captured. This could, in parts, be due to the small training and validation sets which are available for \(E\) (see Table <ref>) and, in particular, to the lack of Rb2O-B2O3 glass samples which the models could base their predictions on. We explain the latter point in more detail below.
Nevertheless, whereas the mean curve for \(E\) shows at least a physically reasonable trajectory in the region of low B2O3-fractions, where there are no data points of alkali borate glasses available, the mean curves for \(T_g\) and \(G\) show a non-physical incline for glasses of B2O3-fractions of less than 0.2 and 0.4, respectively.
As a further observation, we note that, in general, for all three properties, the uncertainty of the model ensembles' predictions in terms of the width of the confidence band around the mean curve is much higher than in the blind settings, especially in the regions of low B2O3-fractions where there are only few or no alkali borate glass samples available. This indicates less robustness with respect to training the models and is also most noticeably reflected by the large deviation of the worst performing model's curve from the mean curve for all three properties.

As most probable explanation, we suspect these observations to be caused by the choice of our training and test sets. As already indicated in Sect. <ref>, we note that the curves of all alkali borate glasses show a similar trajectory for all three properties since these glasses have similar material properties. We also note that only Na2O-B2O3 glasses are not present in the training and validation sets. In regions where the other alkali borate glasses are available in the training and validation sets, the models are thus, in principal, able to learn the properties of the Na2O-B2O3 glasses based on the other alkali borate glass samples. In regions where there are many of these samples available and where their property curves are very close to the Na2O-B2O3 curve, the informed models' predictions thus tend to be quite accurate.
In regions where only few or no alkali borate glasses are available in the training and validation sets, the models are prone to base their predictions on other spurious or noisy features. This leads to non-physical predictions with a high uncertainty. This phenomenon is amplified by a large feature set and is therefore pronounced to a much higher degree in the informed setting than in the blind one.

In summary, the informed model shows, on average in the ensemble setting, clear superior performance than all considered blind (uninformed) models in extrapolating the property curves of Na2O-B2O3 binary glasses for all three properties \(T_g\), \(E\), and \(G\). This is in terms of quantitative error measurements on the test sets as well as in the qualitative approximation of the property curves.

Finally, we emphasize the importance of the ensemble setting. Whereas single models might yield bad predictions, averaging multiple trained models, as we observe in our specific use case, often yields good approximations of the target quantities <cit.>.

§ CONCLUSION AND OUTLOOK

In this paper, we presented an informed neural network approach for the prediction of three material properties of binary oxide glasses, that is, glass transition temperature \(T_g\), Young's modulus \(E\) (at room temperature), and shear modulus \(G\). We compared this approach to five different blind (uninformed) models for all three properties and demonstrated its superior average extrapolation power when applied in an ensemble setting to alkali borate glass samples which contain sodium as previously unseen element.

In terms of the taxonomy of informed machine learning introduced in <cit.>, we integrated prior knowledge into our learning pipeline at four major points. We integrated scientific knowledge, represented as a weighted graph, knowledge graph, and spatial invariance in the training data and in the hypothesis set, respectively. Moreover, we integrated expert knowledge, represented as algebraic equations, into the final hypothesis.

Our informed neural network model could be improved in various ways. First, the list of chemical and physical element features could be extended. Second, instead of classifying glass oxides into formers and non-formers, we could follow the refined classification into formers, modifiers, and intermediates and treat these three classes by three separate neural networks. Third, in this paper, we did not tune any of the models' hyperparameters. A thorough hyperparameter study probably leads to improved model performance. Finally, by relying on further expert knowledge, we could potentially filter out even more predicted property curves in the post-processing step than just constant predictions. This might improve the final predictions even further.

Our results show that our informed neural network model is capable of meaningfully extrapolating various properties of binary glass samples with previously unseen compounds. As a next step, we plan to scale up our approach in order to make it applicable to oxide glass samples with three or more compounds. We also plan to make it more universal, such that it can accurately predict more material properties.

This work was supported in part by the BMBF-project 05M2AAA MaGriDo (Mathematics for Machine Learning Methods for Graph-Based Data with Integrated Domain Knowledge), by the Fraunhofer Cluster of Excellence Cognitive Internet Technologies, and by the Deutsche Forschungsgemeinschaft (DFG, German Research
Foundation) via project 390685813 - GZ 2047/1 - Hausdorff Center for Mathematics (HCM).

\(T_g\)-values of binary alkali borate glasses. Scattered points represent \(T_g\)-values given in the cleaned \(T_g\)-dataset. Solid lines show the predictions for the \(T_g\)-value of Na2O-B2O3 glass samples of the blind and informed models, respectively. The model ensembles' mean curves are shown as blue solid lines with the shaded blue area depicting the 95% confidence band. The predictions of the best and worst performing models in the ensembles are shown as green and red solid lines, respectively.

\(E\)-values of binary alkali borate glasses. Scattered points represent \(E\)-values given in the cleaned \(E\)-dataset. Solid lines show the predictions for the \(E\)-value of Na2O-B2O3 glass samples of the blind and informed models, respectively. The model ensembles' mean curves are shown as blue solid lines with the shaded blue area depicting the 95% confidence band. The predictions of the best and worst performing models in the ensembles are shown as green and red solid lines, respectively.

\(G\)-values of binary alkali borate glasses. Scattered points represent \(G\)-values given in the cleaned \(G\)-dataset. Solid lines show the predictions for the \(G\)-value of Na2O-B2O3 glass samples of the blind and informed models, respectively. The model ensembles' mean curves are shown as blue solid lines with the shaded blue area depicting the 95% confidence band. The predictions of the best and worst performing models in the ensembles are shown as green and red solid lines, respectively.

Chemical and physical properties extracted from the  library which are used for the correlation study described in Sect. <ref>. Properties that are marked with  are dropped before the correlation study as they are not available for all respective P-elements. Properties that are marked with  are not highly correlated among each other and are used as final features. All unmarked properties are dropped due to too high correlation with other properties. We refer to <cit.> and the references therein for detailed explanations of the properties.

Properties from      \(T_g\)     \(E\)     \(G\)

Atomic number

Atomic radius

Atomic radius by Rahm et al.

Atomic volume

Atomic weight

Boiling temperature

\(C_6\) dispersion coefficient by Gould and Bučko

Covalent radius by Cordero et al.

Single bond covalent radius by Pyykko et al.

Double bond covalent radius by Pyykko et al.

Dipole polarizability

Electron affinity

Electron affinity in the Allen scale

Electron affinity in the Ghosh scale

Electron affinity in the Pauling scale

Glawe’s number

Group in periodic table

Heat of formation

First ionization energy

Lattice constant

Maximum coordination number

Maximum oxidation state

Melting temperature

Mendeleev's number

Minimum coordination number

Minimum oxidation state

Period in periodic table

Pettifor scale

Index to chemical series

Number of valence electrons

Van der Waals radius

Van der Waals radius according to Alvarez

Van der Waals radius according to Batsanov

Van der Waals radius from the MM3 FF

Van der Waals radius from the UFF

Number \(d\) of all features  (including mole atomic fractions)     28     26     26