Machine Learning Crystal Structure Prediction

Machine Learning Crystal Structure Prediction

Hardness, which is defined as the quantitative value of resistance to permanent or plastic deformation, plays a dominant role in materials design for many industrial applications, for example, drilling, boring, abrasives, medical implants, and protective coatings.

Boosting in manufacturing demand spurs the requirements of new materials of diverse hardness.

Which makes it necessary to understand the physical origin of these properties. Hardness testing such as the Mohs scratch test is non-destructive and one material is scratched with another material with a specified hardness number between 1 and 10. Easily scratched materials are given a low Mohs number, whereas materials that possess high resistance to plastic deformation and are difficult to scratch are given a high Mohs number. For example, diamonds have a Mohs number of 10. (Garnett, 2019) 

Many factors determine the Mohs hardness of a material. The hardness of minerals, which are single-crystal brittle materials, depends on chemical bonding. While only van der Waals bonding between molecular chains for sheets, the materials typically have low Mohs hardness. Materials with larger Mohs hardness usually possess ionic or electronic bonging. For materials at the top of the Mohs scale, they have large covalent bonds which restrict the start of dislocations under stress that produce resistance to plastic deformation.

Besides bonding, hardness is also relevant to the correlation between composition and bond strength. Light elements normally have short and strong bonds. Higher Mohs hardness is related to the high average bond length, the high number of bonds per unit volume, and a higher average number of valence electrons per atom. (Garnett, 2019) 

Currently, computational methodologies for hardness deformation and compositional factors include molecular dynamics, density functional theory, and machine learning. However, all methods, especially the frequently used molecular dynamics method, encounter a major challenge of the large computational cost when calculating various materials with varied chemical and structures at the same time. Moreover, especially when considering the cost of  optimization and energy determination of various deformation across a single material class.  (Garnett, 2019)

Therefore, if machine learning approach implements on the naturally occurring ceramic minerals to create a model based on the compositional-based features to predict Mohs hardness across a broad range of structural and chemical space of artificial crystal with similar atomic and structural characteristics, it will make hardness prediction easily and drastically reduced the computational force than the molecular dynamics’ method. Three data analysis methods, principal component analysis (PCA), random forest (RF), and artificial neural networks (ANN) used to analyze 622 naturally occurring ceramic minerals, and models created. The efficiency and results of those models in predicting the Mohs hardness of 52 synthetic crystals with similar atomic and structural characteristics discussed later.

Methods 

Datasets 

The dataset used for training the models originated from experimental Mohs hardness  data, with their crystal classes, and chemical compositions of naturally occurring ceramic  minerals reported in the Physical and Optical Properties of Minerals CRC Handbook of  Chemistry and Physics and the American Mineralogist Crystal Structure Database. (Garnett,  2019) The testing dataset contains the composition, crystal structure, and Mohs hardness values  of 52 synthetic single crystals reported in the literature. Both datasets were analyzed by Python  scripts. The compositional features used to describe the characteristic of naturally occurring  materials, which were obtained from the Physical and Optical Properties of Minerals CRC  Handbook, were used to train models to predict hardness across a wide variety of mineral  compositions and crystal classes and were listed in Table 1 (Garnett, 2019) 

Table 1 Primary Features 

Principal Component Analysis (PCA) 

As listed in Table 1, multiple features involved in predicting the Mohs hardness.  PCA is an appropriate technique to reduce the dimensions of the data because it can simplify  the complexity in high-dimensional data while retaining trends and patterns. It achieves this by  transforming the data into a few dimensions, which perform as summaries of features. As an  unsupervised learning method, PCA is similar to clustering, and it can find patterns without  previous knowledge about the origin treatment group of the samples or whether it has  phenotypic differences.

The principle of PCA is to reduce data by finding the directions of  maximum variance in high-dimensional data and geometrically projecting them onto lower  dimensions named principal components (PCs), and a limited number of PCs are chosen to  achieve the best summary of the data. All PCs are geometrically orthogonal to make sure of the  uncorrelation between PCs. (Lever, Krzywinski, & Altman, 2017) 

Random Forest (RF) 

Since RF is a supervised machine learning algorithm used widely in solving  classification and regression problems, it is chosen to conduct regression on all those  compositional features to predict the Mohs hardness. It is powerful because it is the  combination of many decision trees on the different samples with random bagging and feature,  and the outcome of it is by taking the majority vote of all decision trees for classification and  average in case of regression. (“Random Forest | Introduction to Random Forest Algorithm,”  n.d.)

Artificial Neural Networks (ANNs) 

ANNs was chosen because it is a subset of machine learning which is the heart of deep  learning algorithms. Inspired by the human brain, ANNs mimic the way that biological neurons  signal to one another. It is composed of a node layer, containing an input layer, one or more  hidden layers, and an output layer. Each node, which is also called an artificial neuron, connects  with the associated weight and threshold. The node is activated once the output of any  individual node exceeds the specified threshold value and then sends data to the network’s next  layer. Otherwise, no data is passed to the next layer of the network. (“What are Neural  Networks? | IBM,” n.d.) 

Result 

PCA Results on Dimension Reduction  

Standardization of the dataset was performed before finding the PCA components to  avoid the influence of the difference in scales of features on the result. Then explained variance  and explained variance ratio of principal components, which can evaluate the significance of  principal components, were calculated and the results are shown in Table 2.  

Table 2 Explained Variance and Explained Variance Ratio of PCs 

Explained variance indicates the amount of variance explained by each of the selected components and explained variance ratio is the percentage of variance explained by each of the selected components.

As shown in Table 2, 64.35% of the variance is contained in the first two  principal components and 79.68% in the first three principal components. The dominance of  the first three components in the variance is also revealed by the cumulative explained variance and the explained variance ratio plots in Figure 1 and Figure 2, where cumulative explained  variance line increase dramatically in the first three components, and explained variance and  explained variance ratio drop sharply in first three components. 

The reduced dimensional representation of data, which is also called score plot, using  only PC1 and PC2, and PC2 and PC3 are presented in Figure 3 and Figure 4 respectively.  Accompany with the loading plots of PC2 and PC1 in Figure 5 which represents the location  of all eigenvectors on PC1 and PC2 and shows the contribution that each variable makes to  each eigenvector, it is reasonable to conclude that a few minerals have a high atomic average  first ionization energy, atomic average Pauling electronegativity of the most common oxidation  state, and atomic average number of valence electrons because a few data point are cluster on  the bottom which are the positive directions of all those three features.

What’s more, since the first three PCs contribute a high variance of 79.68%, 3D loading plots of PC1, PC2, and PC3  were plotted and can also be found in Figure 6.  

Biplot combines both the principal component scores and the loading vectors in a single plot and it is a powerful way for visualization of results from PCA. 2D biplot of PC1 and PC2, and 3D biplot of PC1, PC2, and PC3 are shown in Figure 7 and Figure 8. Biplot is more  powerful in interpreting the result because the orientation, length in space, and angles between  vectors all possess meanings. The more parallel the vector to a principal component axis, the  more contribution it makes to that PC. (“Interpretation and Visualization • SOGA • Department of  Earth Sciences,” n.d.)

Therefore, from Figure 7, since vectors of atomic average elemental density, atomic average covalent atomic radius, the atomic average number of electrons, and  atomic average atomic weight are close to PC1 axis orientation, it is reasonable to conclude  those features contribute more to PC1, and this is verified in Figure 9, which shows the effect  of variables on each component, where they all show correlation value greater than 0.4.

The feature that contributes more to PC2 is the atomic average Pauling electronegativity of the most common oxidation state, which has an absolute correlation of 0.51 reported in Figure 9. The  length of vectors in the space indicates the variability of this variable represented by the two  displayed principal components. Besides that, the longer the vector, the more variability. 

Shorter vectors mean those variables are better represented in another dimension. 

(“Interpretation and Visualization • SOGA • Department of Earth Sciences,” n.d.) Therefore,  it makes sense that total elemental density and the total number of electrons have short vectors  because they have more effect on PC3, where their correlation exceeds 0.6 as indicated in  Figure 9.

Furthermore, the angles between vectors of different variables show their correlation  in this space. Small angles represent high positive correlation, right angles represent lack of  correlation, and opposite angles represent the high negative correlation. (“Interpretation and  Visualization • SOGA • Department of Earth Sciences,” n.d.) Thus, the small angle between  atomic average Pauling electronegativity of the most common oxidation state and atomic  average first ionization energy indicates they are highly correlated. Atomic average atomic  number to mass number ratio and atomic average covalent atomic radius are uncorrelated  because they are orthogonal to each other and this makes physical sense.  

Model Performance and Feature Reduction of RF 

Taking the data of 622 natural occurrences of ceramic minerals as training data and data  of 52 synthetics crystals as the testing model and passing them to random forest regression. The  result shows that the root-means-square error of test data is 0.93176, the R-squared value of  training data is 0.92966, and the R-squared value of testing data is 0.61858, which indicates  that the accuracy of that random forest regression model is not high enough to use as a valid  model. This result can be visualized in Figure 11. Though all data points are clustered around  the red line, which represents the perfect prediction of testing data, they all have some distance  to the line, which also proves that the random forest model is not accurate enough.  

Singletree and random forests with multiple numbers of trees were also investigated to  see whether random forests help boost the accuracy of prediction and root mean square error of testing data, R-squared value of testing and training data of each model are listed in Table  3. 

Table 3 R-squared results for different numbers of trees

Consistent with the trend shown in Figure 10, which plots the R-squared value versus  the number of trees, the accuracy of random forest does not have obvious improvement after the  number of trees reached 13, and the R-squared value fluctuates around 0.6 when the number  of trees exceeds 20. This result indicates that random forest performs high accuracy than a  single tree, but though adding a few numbers of trees contributes to a drastic boost inaccuracy of the model, adding too many trees would not help increase accuracy continuously, and  contrarily R-squared value would reach a ceiling. 

Feature importance can describe which features are relevant and how much each feature decreases the variance in the tree.

Though one of the two correlated features shows up high while the other is low. This effect can be reduced by random selection of features during each  node creation. Moreover, the importance of features is ranked in Table 4. The visualization of the trend  in feature importance can be found in Figure 12 and Figure 13. Which are the line plot of feature  importance versus variables and bar plot of feature importance of each variable respectively. 

Table 4 Feature Importance 

The result shows that none of the features have low importance, which is normally defined as lower than 0.005. For this dataset, 0.05 is set as the lowest bar of importance. Moreover, further analysis was performed after removing the data of atomicweight_Average, and the dendrogram is presented in Figure 14. 

Performance of ANNs 

Rectified linear unit, which returns 0 if it receives a negative input and returns the same positive value if it receives a positive value, was used as the activation function because it is a commonly used activation function in deep learning models to account for non-linearities and  specific interaction effects. (“Rectified Linear Units Definition | DeepAI,” n.d.) Mean squared  error (MSE), which is the sum of the squared differences between the prediction and true values, was chosen as the cost function. Regressor ANN was performed on the training and testing  dataset.

The prediction values of hardness based on synthetic crystal’s compositional features were compared to the testing values and the result is illustrated in Figure 15. Datapoint is clustered but distant from the red line, which represents a perfect prediction of values.  Therefore, it is reasonable to conclude that the ANN regressor cannot predict the hardness  values of artificial crystal accurately though both the test and training loss manifested a  significant decrease over epochs, which can be found in Figure 16.

Discussion 

Himanen’s article provides an overview of historical development and the current state  of data-driven materials science and reviews the importance of successes and challenges so far.  In this article, the FAIR Data Principles were mentioned, which were established to provide  unifying guidelines for the widely varying groups interested in Open Data and to aid in data  management development. (Himanen, Geurts, Foster, & Rinke, 2019) 

FAIR in the context of the database means findability, accessibility, interoperability,  and reuse of digital assets. Which is an idea that was first proposed in an articled named ‘FAIR  Guiding Principles for scientific data management and stewardship’ published in Scientific  Data in 2016 (“FAIR Principles – GO FAIR,” n.d.). This principle could be simply interpreted  as the convenience of machine actionability on scientific data, which can minimize human  intervention.  

Furthermore, Himanen also emphasized the importance of the Open Science movement.

Which is the foundation of data-driven materials science. Grateful for the Open Data, it is convenient to find a published paper in materials sciethathich offers the data behind it currently. In addition, as an illustration, while confining the search period into the last five years, after reviewing six papers, the article named ‘Interaction of a tripeptide with titania surfaces: RGD adsorption rutile TiO2(110) and model dental implant surfaces’ showed up with VAMAS data file of raw data for the SLA and SLActive smiles following edition of RGD at different pH levels provided  on Mendeley Data. 

Table 1 in Himanen’s article listed the major materials data infrastructures at the time of the paper, divided into non-commercial and commercial. The situation concerning data does not undergo a drastic change from 2019 when the paper was published, but with the help of the rapid development of technology, searching and sharing data becomes more convenient. 

The Open Quantum Materials Database (OQMD) is one of the services provided by  various materials data resources centers listed in Table 2. The OQMD is a database of DFT  calculated thermodynamic and structural properties of 1,022,603 materials. There are three  data analysis tools provided on this database. The first one is Grand Canonical Linear  Programming (GCLP), which is a method developed by Wolverton et. al. for analyzing the  complex ground state thermodynamics of metal hydrides. The method works by mapping a  free energy minimization problem to a linear algebra problem. Using linear programming  allows this problem to be solved nearly instantaneously even with a very large number of  phases in very high-dimensional phase space. 

Additionally, it can use the qhull algorithm to create phase diagrams for arbitrary regions of phase space.

Depending on the dimensionality of the requested phase diagram, different plots are made, where 2D phase diagrams are the convex hull, 3D phase diagrams are the Gibbs’ triangle, and the 4D phase diagram is an  isometric shot of the Gibbs tetrahedron. The formation energies supplied to this tool were  calculated using elemental reference states which are a combination of DFT total energies of  elemental phases, fit to experimental formation energies, corrected to remove the effects of  DFT+U. The last tool is to visualize a crystal structure, which is accomplished by inputting the  structure data in a format readable by JSMol. (“OQMD | The Open Quantum Materials  Database,” n.d.) 

An overview of the main types of machine learning and the illustration of how they can add value to material infrastructures are mentioned in section & of Himanen’s article. It emphasizes that supervised learning is the most mature and powerful among the subfields of machine learning and is widely used in the majority of machine learning studies in the physical  sciences, this technique is typically applied on predictions of physical properties, which is  exactly the method like RF and ANNs used by this research to predict hardness of artificial crystals. 

Additionally, Himanen’s paper discussed the function of unsupervised learning, which often can be utilized to preprocess input data such as PCA, one of the simplest forms of feature learning, which is also used in this study to reduce dimension. Since the goal of PCA is to find  the minimum number of linear combinations of the variables that explain a large fraction of  the data, the first three components for this dataset, which contribute a total of 79.68% of the  variance, are an appropriate approximation of the variance present in an original 11D dataset. 

Furthermore, PCA result reveals that atomic average elemental density, atomic average  covalent atomic radius, the atomic average number of electrons, and atomic average atomic  weight contribute more to PC1, atomic average Pauling electronegativity of the most common  oxidation state contributes more to PC2, and total elemental density and the total number of  electrons have more effect on PC3. 

Both RF and ANN analysis failed to provide a mature model from naturally occurring crystals to predict the hardness of artificial crystals. RF analysis manifests that as the number of trees grows to around 20, the R-squared value fluctuates around 0.6. The possible cause of reaching the ceiling of R-squared which is not high enough and the inaccurate ANN  model could be the integration of different crystal structures in a single training data without  identifying the difference crystal structure by class.

Furthermore, Himanen mentioned in his paper that the building machine learning model can be broken down into four key steps: data acquisition, feature engineering, model building,  and analysis. Compared to the analysis in this study, the feature engineering step is missing. 

Feature engineering refers to the act of introducing domain expertise to the learning model by affecting which features are used.

Feature extraction, which is problem-specific feature transformations, is also beneficial. It can make it easier for the model to learn a unique mapping  by exploiting known symmetries in the input exhibiting symmetries concerning translation, rotation, and index permutation in the Cartesian coordinates representing a geometry, physical properties transform can make the input invariant concerning these symmetries and will help the  learning process by creating a unique mapping from an atomic geometry to its properties. Such  structural descriptors have been successfully applied in the prediction of molecular and  crystalline properties. (Himanen et al., 2019) Since this study aims to use the compositional  feature to predict the hardness of crystals, applying to transform on structure features may boost  the accuracy of the machine learning result. 

What’s more, Himanen also emphasized that the four key steps are independent, and successful machine learning systems require multiple iterations of each step. Therefore,  repetition of each step in RF and ANN may also contribute to the increased inaccuracy of the  current result. 

Conclusions 

PCA manifests successful reduction of the 11D features into 3D, which contributes 79.68% of the variance in the dataset. Also, the result shows that atomic average elemental density, atomic average covalent atomic radius, the atomic average number of electrons, and atomic average atomic weight contribute more to PC1, atomic average Pauling electronegativity of the most common oxidation state contributes more to PC2, and total elemental density and the total number of electrons have more effect on PC3. 

RF model illustrates the boosting in the accuracy of mode until the number of trees grows beyond 20, where the R-squared value starts to fluctuate around 0.6. ANNs results also show a poor fit of predicted values to actual values, which indicates this model has even worse accuracy than RF. 

For future work, information about crystal structure should be acquired from Garnett’s group. Moreover, to perform further classification analyses like classification ANN and RF. Also, as mentioned in Himanen’s paper, the four key steps, data acquisition. In addition, feature engineering, model building, and analysis. Are so important to the success of the machine learning model and need  to be repeated individually. Laslty, therefore, feature engineering and feature extraction processes should be included to boost the accuracy of the model.

Back To News

Machine Learning Crystal Structure Prediction References 

FAIR Principles – GO FAIR. (n.d.). Retrieved March 5, 2022, from https://www.go fair.org/fair-principles/ 

Garnett, J. C. (2019). Prediction of Mohs Hardness with Machine Learning Methods Using  Compositional Features [Chapter]. ACS Symposium Series, 1326, 23–48.  https://doi.org/10.1021/bk-2019-1326.ch002 

Himanen, L., Geurts, A., Foster, A. S., & Rinke, P. (2019). Data-Driven Materials Science:  Status, Challenges, and Perspectives. Advanced Science, 6(21).  https://doi.org/10.1002/advs.201900808 

Interpretation and Visualization • SOGA • Department of Earth Sciences. (n.d.). Retrieved  March 4, 2022, from https://www.geo.fu-berlin.de/en/v/soga/Geodata-analysis/Principal Component-Analysis/principal-components-basics/Interpretation-and 

visualization/index.html 

Lever, J., Krzywinski, M., & Altman, N. (2017). Points of Significance: Principal component  analysis. Nature Methods, 14(7), 641–642. https://doi.org/10.1038/nmeth.4346 

OQMD | The Open Quantum Materials Database. (n.d.). Retrieved March 5, 2022, from  https://oqmd.org/analysis/visualize 

Random Forest | Introduction to Random Forest Algorithm. (n.d.). Retrieved March 3, 2022,  from https://www.analyticsvidhya.com/blog/2021/06/understanding-random-forest/ 

Rectified Linear Units Definition | DeepAI. (n.d.). Retrieved March 4, 2022, from  https://deepai.org/machine-learning-glossary-and-terms/rectified-linear-units 

What are Neural Networks? | IBM. (n.d.). Retrieved March 3, 2022, from  https://www.ibm.com/cloud/learn/neural-networks

Machine Learning Crystal Structure Prediction